w1n5t0n/memlnaut-nisps
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

feat(nisps/ml): MLP library with fixed-architecture template + spread-aware RL (meml-wmh)

Browse source 
Stream 2 of the clean-slate rewrite: nisps/ml/ replaces src/memlp/ with a
header-only, heap-free MLP that satisfies nisps::core::MLEngine.

Files (nisps/ml/):
- activations.hpp — ReLU (leaky 0.01 for parity), sigmoid, tanh
- loss.hpp — MSE per-sample (fixes meml-ues double-scaling: returns the
  sample's MSE without an extra 1/N multiplication; the training loop
  averages explicitly)
- init.hpp — uniform/Xavier/spread-aware weight init
- training.hpp — gradient clip helper (±10.0 matches legacy)
- rl.hpp — move_weights with per-layer Xavier scaling, weight decay
  (10% * spread), gaussian noise via the deterministic Rng (matches the
  legacy JS sum-of-three-uniforms shape); draw_weights also spread-aware
- stats.hpp — per-layer mean/max/dead/saturating diagnostics
- mlp.hpp — 4-layer (3 hidden + sigmoid output) MLP class with
  std::array-backed weights, biases, gradient accumulators, dataset
  ring buffer (default 128 examples), loss history (default 4096 iters).
  Bias is a separate per-layer parameter — no input-vector mutation.
  Flat get_weights/set_weights layout: weights all layers (row-major,
  layer order), then biases all layers.

Tests (tests/cpp/, all 50 passing under -Wall -Wextra -Werror -Wpedantic):
- test_mlp_init.cpp — deterministic seeding, spread regimes,
  static_assert MLEngine concept satisfied
- test_mlp_inference.cpp — golden hand-computed forward pass match,
  sigmoid output range, set_input bounds
- test_mlp_training.cpp — XOR convergence (loss < 0.01 in <2k iters),
  ring-buffer eviction
- test_mlp_loss.cpp — meml-ues regression test: reported loss equals
  hand-computed average MSE without extra 1/N scaling; sample weights
  honoured
- test_mlp_rl.cpp — move_weights respects output_pin_mask (final-layer
  rows + biases preserved); spread regimes; grad clear after draw_weights
- test_mlp_serialize.cpp — get_weights/set_weights round-trip preserves
  inference exactly; eval_loss is non-mutating; infer_batch matches
  individual inference

Verification:
- Clean build, no warnings
- 50 tests pass (22 prior + 28 new)
- No std::vector / new / malloc in nisps/ml/
- All float literals .f-suffixed in code (comments excepted)
This commit is contained in:
w1n5t0n 2026-04-29 15:55:43 +03:00
parent 55e7bc9654
commit 825ed6ad33
14 changed files with 1734 additions and 0 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

6
nisps/CMakeLists.txt
View file

@ -51,6 +51,12 @@ if(NOT EMSCRIPTEN)
${NISPS_TEST_DIR}/test_ring_buffer.cpp
${NISPS_TEST_DIR}/test_rng.cpp
${NISPS_TEST_DIR}/test_math.cpp
${NISPS_TEST_DIR}/test_mlp_init.cpp
${NISPS_TEST_DIR}/test_mlp_inference.cpp
${NISPS_TEST_DIR}/test_mlp_training.cpp
${NISPS_TEST_DIR}/test_mlp_loss.cpp
${NISPS_TEST_DIR}/test_mlp_rl.cpp
${NISPS_TEST_DIR}/test_mlp_serialize.cpp
)
target_link_libraries(nisps_core_tests PRIVATE nisps_core)

91
nisps/ml/activations.hpp Normal file
View file

@ -0,0 +1,91 @@
// nisps/ml/activations.hpp — activation functions used by the MLP.
//
// Three activations are needed for the firmware default: ReLU (hidden layers),
// sigmoid (output layer), and tanh (alternative). Each entry exposes the
// activation and its derivative-given-pre-activation; we keep both around
// because the backprop path multiplies by the derivative of the
// pre-activation, which is cheaper to compute on the inner product than on
// the post-activation in some cases.
//
// Notes on the leak in ReLU:
// The legacy firmware uses leaky ReLU with slope 0.01 everywhere it claims
// to use "ReLU". We preserve that behavior for parity — a true zero-slope
// ReLU would silently change training dynamics in ways that aren't
// documented. The leak slope is a `static const float` so it lives in SRAM
// per Chris's rule.
//
// Sigmoid:
// We use the exact std::exp version here (NOT fast_sigmoid). The MLP
// training loop multiplies by the sigmoid derivative, and the ~1.2% bias
// in fast_sigmoid would compound during training. fast_sigmoid is fine for
// inference-only paths but not for the gradient path. The ML output
// activation is the place where we want monotonic-AND-accurate.
//
// All literals carry the .f suffix; constants used inside loops live in
// `static const float` so they are hoisted to SRAM on RP2350.
#pragma once
#include <cmath>
#include "../core/perf.hpp"
namespace nisps::ml {
// Match the legacy firmware leaky-ReLU slope. 0.01 matches PyTorch default
// and what `src/memlp/Utils.h::kReLUSlope` had before.
inline constexpr float kReluLeakSlope = 0.01f;
NISPS_FORCE_INLINE float relu(float x) noexcept {
return (x > 0.f) ? x : kReluLeakSlope * x;
}
NISPS_FORCE_INLINE float relu_deriv_pre(float pre_activation) noexcept {
return (pre_activation > 0.f) ? 1.f : kReluLeakSlope;
}
NISPS_FORCE_INLINE float sigmoid(float x) noexcept {
// Saturate inputs to avoid expf overflow / underflow noise. Mirrors
// nisps::exact_sigmoid clamps from core/math.hpp.
if (x > 40.f) return 1.f;
if (x < -40.f) return 0.f;
return 1.f / (1.f + std::exp(-x));
}
NISPS_FORCE_INLINE float sigmoid_deriv_pre(float pre_activation) noexcept {
const float s = sigmoid(pre_activation);
return s * (1.f - s);
}
NISPS_FORCE_INLINE float tanh_act(float x) noexcept {
return std::tanh(x);
}
NISPS_FORCE_INLINE float tanh_deriv_pre(float pre_activation) noexcept {
const float t = std::tanh(pre_activation);
return 1.f - t * t;
}
// Activation kind, dispatched at compile time per layer (see Layer template
// in mlp.hpp). We don't use a runtime tag because activation choice is
// architectural, not per-call.
enum class Activation : int {
ReLU = 0,
Sigmoid = 1,
Tanh = 2,
};
template <Activation A>
NISPS_FORCE_INLINE float activate(float x) noexcept {
if constexpr (A == Activation::ReLU) return relu(x);
if constexpr (A == Activation::Sigmoid) return sigmoid(x);
if constexpr (A == Activation::Tanh) return tanh_act(x);
return x; // unreachable
}
template <Activation A>
NISPS_FORCE_INLINE float activate_deriv_pre(float pre) noexcept {
if constexpr (A == Activation::ReLU) return relu_deriv_pre(pre);
if constexpr (A == Activation::Sigmoid) return sigmoid_deriv_pre(pre);
if constexpr (A == Activation::Tanh) return tanh_deriv_pre(pre);
return 1.f; // unreachable
}
} // namespace nisps::ml

72
nisps/ml/init.hpp Normal file
View file

@ -0,0 +1,72 @@
// nisps/ml/init.hpp — weight initialization strategies.
//
// Three strategies are exposed:
// - `uniform_init`: w ~ U[-1, 1]. Matches the legacy
// `utils::gen_rand<T>()` baseline; produces high-magnitude pre-activations
// that drive sigmoid output saturation.
// - `xavier_init`: w ~ U[-1, 1] / sqrt(fan_in). Centered pre-activations,
// better for sigmoid output layers.
// - `spread_init`: linearly interpolates per-layer scale between the two.
// scale = (1 - spread) + spread / sqrt(fan_in). spread=0 ⇒ uniform,
// spread=1 ⇒ Xavier. This is the playground knob.
//
// All three operate on a flat row-major weight buffer of size fan_in*fan_out.
// Biases are initialized separately and ALWAYS to zero — the legacy code
// initialized biases to zero, and the playground spread parameter never
// touches biases. We keep that.
#pragma once
#include <cmath>
#include <cstddef>
#include <span>
#include "../core/rng.hpp"
namespace nisps::ml {
// Compute the spread-aware weight scale for one layer.
// spread=0 → 1.0 (uniform [-1, 1])
// spread=1 → 1/sqrt(fan_in) (Xavier)
inline float spread_scale(float spread, std::size_t fan_in) noexcept {
if (fan_in == 0u) return 1.f;
const float inv_sqrt = 1.f / std::sqrt(static_cast<float>(fan_in));
return (1.f - spread) + spread * inv_sqrt;
}
// Initialize one layer's weights (flat row-major) and biases.
// weights: span of size fan_in * fan_out
// biases: span of size fan_out (zeroed)
// spread: see spread_scale
// rng: state advanced; caller owns it
inline void spread_init(std::span<float> weights,
std::span<float> biases,
std::size_t fan_in,
float spread,
Rng& rng) noexcept {
const float scale = spread_scale(spread, fan_in);
for (std::size_t i = 0; i < weights.size(); ++i) {
weights[i] = rng.next_float_signed() * scale;
}
for (std::size_t i = 0; i < biases.size(); ++i) {
biases[i] = 0.f;
}
}
// Convenience aliases for the endpoints. Cheap to build on top of
// spread_init; useful for tests that want to assert a specific regime.
inline void uniform_init(std::span<float> weights,
std::span<float> biases,
std::size_t fan_in,
Rng& rng) noexcept {
spread_init(weights, biases, fan_in, 0.f, rng);
}
inline void xavier_init(std::span<float> weights,
std::span<float> biases,
std::size_t fan_in,
Rng& rng) noexcept {
spread_init(weights, biases, fan_in, 1.f, rng);
}
} // namespace nisps::ml

62
nisps/ml/loss.hpp Normal file
View file

@ -0,0 +1,62 @@
// nisps/ml/loss.hpp — MSE loss with a deliberate fix for the meml-ues
// double-scaling bug.
//
// THE BUG (in `src/memlp/Loss.h::MSE` and `src/memlp/MLP.cpp::Train`):
// The legacy MSE computed loss as
// sum((diff^2) / NOut) * sampleSizeReciprocal in Loss.h
// and then `Train()` did
// current_iteration_cost_function *= sampleSizeReciprocal; // ← second
// So the unweighted average loss was scaled by 1/N twice — once in MSE()
// and once in the outer training loop. Loss values for the same training
// set get reported 1/N too small, which is misleading for diagnostics.
//
// THE FIX:
// Loss is computed PER-SAMPLE here. The training loop is responsible for
// averaging across samples once. There is no implicit per-sample weight in
// the MSE function — that is the caller's policy decision.
//
// For a single sample: mse = (1/NOut) * sum_j (label_j - pred_j)^2
// Derivative wrt pred_j: -(2/NOut) * (label_j - pred_j)
//
// When sample weights are supplied to `train()` they are applied at the
// sample level (loss_total = sum_i weight_i * mse_i), and the gradient
// per-sample is scaled by weight_i. Sample weights must sum to 1.0 — the
// caller normalizes (matches legacy contract).
#pragma once
#include <cstddef>
#include <span>
#include "../core/perf.hpp"
namespace nisps::ml {
// Compute MSE loss for one sample and write the per-output derivative into
// `loss_deriv`. Returns the scalar loss for this sample.
//
// loss = (1/NOut) * sum_j (label_j - pred_j)^2
// d(loss)/d(pred_j) = -(2/NOut) * (label_j - pred_j)
//
// `loss_deriv` and `pred` must both have size NOut. `label` likewise.
//
// NB: We deliberately do NOT multiply by any sampleSizeReciprocal here. The
// caller is responsible for averaging across samples (and scaling by sample
// weights, if any) — see training.hpp.
NISPS_FORCE_INLINE float mse_per_sample(std::span<const float> label,
std::span<const float> pred,
std::span<float> loss_deriv) noexcept {
const std::size_t n = pred.size();
if (n == 0u) return 0.f;
const float inv_n = 1.f / static_cast<float>(n);
float accum = 0.f;
for (std::size_t j = 0; j < n; ++j) {
const float diff = label[j] - pred[j]; // sign chosen so deriv is wrt pred
accum += diff * diff * inv_n;
loss_deriv[j] = -2.f * inv_n * diff;
}
return accum;
}
} // namespace nisps::ml

570
nisps/ml/mlp.hpp Normal file
View file

@ -0,0 +1,570 @@
// nisps/ml/mlp.hpp — fixed-architecture MLP, four layers (three hidden +
// output). All buffers are template-sized std::array; zero heap allocation
// in inference, training, and the dataset path.
//
// ARCHITECTURE
// MLP<NIn, NHidden1, NHidden2, NHidden3, NOut, NMaxExamples = 128>
// ┌──────┐ Linear+Bias ┌────────┐ ReLU ┌────────┐ ReLU ┌────────┐ Sigmoid
// │ NIn │ ─────────────▶ │ NH1 │ ──────▶ │ NH2 │ ──────▶ │ NH3 │ ──────▶ NOut
// └──────┘ └────────┘ └────────┘ └────────┘
// Layer 0 (NIn → NH1) ReLU
// Layer 1 (NH1 → NH2) ReLU
// Layer 2 (NH2 → NH3) ReLU
// Layer 3 (NH3 → NOut) Sigmoid
//
// We support exactly three hidden layers. The legacy firmware default is
// [10, 10, 14], so the MVP signature directly matches `MLP<NIn, 10, 10, 14,
// NOut>`. Variable layer count is deferred — see architecture.md.
//
// MEMORY MODEL
// Per layer L_k with fan_in = N_in[k], fan_out = N_out[k]:
// std::array<float, fan_in*fan_out> weights // row-major
// std::array<float, fan_out> biases
// std::array<float, fan_out> pre_activation // cached for backprop
// std::array<float, fan_out> activation // cached for backprop
// std::array<float, fan_in*fan_out> grad_w_accum // for backprop
// std::array<float, fan_out> grad_b_accum
//
// Per MLP:
// std::array<float, NIn> input_buffer (current set_input values)
// std::array<float, NOut> output (post-final-activation; outputs())
// std::array<float, NMaxExamples * NIn> dataset_features
// std::array<float, NMaxExamples * NOut> dataset_labels
// std::size_t dataset_count, dataset_head (FIFO ring buffer)
// std::array<float, NMaxIter> loss_history (max iters from train())
// Rng rng_
// std::array<float, NOut> bp_err_buf, bp_delta_buf (backprop scratch)
//
// FLAT WEIGHT LAYOUT (`get_weights` / `set_weights`)
// [layer0_weights ...] [layer1_weights ...] [layer2_weights ...] [layer3_weights ...]
// [layer0_biases ...] [layer1_biases ...] [layer2_biases ...] [layer3_biases ...]
// Documented in detail near `weight_count()`.
//
// CONCEPT SATISFACTION
// The class satisfies `nisps::MLEngine`:
// set_input, process, outputs, add_example, train (no-arg overload
// returning float), move_weights(speed, spread), draw_weights(spread),
// reset, seed.
// Plus diagnostics required by the broader API (see Stream 2 brief in
// architecture.md): eval_loss, layer_stats, get/set_weights, weight_count,
// infer_batch, loss_history.
#pragma once
#include <array>
#include <cmath>
#include <cstddef>
#include <cstdint>
#include <span>
#include "../core/concepts.hpp"
#include "../core/perf.hpp"
#include "../core/rng.hpp"
#include "activations.hpp"
#include "init.hpp"
#include "loss.hpp"
#include "rl.hpp"
#include "stats.hpp"
#include "training.hpp"
namespace nisps::ml {
// Layer<FanIn, FanOut, Act>. Stores its weights, biases, and the work
// buffers needed for forward + backprop. Header-only, all sizes compile-
// time. Each method is small; the compiler will inline through.
template <std::size_t FanIn, std::size_t FanOut, Activation Act>
struct Layer {
static constexpr std::size_t kFanIn = FanIn;
static constexpr std::size_t kFanOut = FanOut;
static constexpr Activation kAct = Act;
std::array<float, FanIn * FanOut> weights{};
std::array<float, FanOut> biases{};
// Cached during forward(); consumed during backprop().
std::array<float, FanOut> pre_activation{};
std::array<float, FanOut> activation{};
// Gradient accumulators — used per-sample for SGD weight update.
std::array<float, FanIn * FanOut> grad_w{};
std::array<float, FanOut> grad_b{};
NISPS_FORCE_INLINE float& w(std::size_t node, std::size_t in) noexcept {
return weights[node * FanIn + in];
}
NISPS_FORCE_INLINE float w(std::size_t node, std::size_t in) const noexcept {
return weights[node * FanIn + in];
}
// Forward: compute pre_activation and activation given an input span.
NISPS_HOT NISPS_FORCE_INLINE
void forward(std::span<const float, FanIn> input) noexcept {
for (std::size_t node = 0; node < FanOut; ++node) {
const std::size_t row = node * FanIn;
float sum = biases[node];
for (std::size_t j = 0; j < FanIn; ++j) {
sum += weights[row + j] * input[j];
}
pre_activation[node] = sum;
activation[node] = activate<Act>(sum);
}
}
// Compute incoming-error vector for the previous layer:
// delta_in[j] = sum_node (err_signal[node] * w[node, j])
// Where err_signal[node] = upstream_err[node] * d/dpre activation.
// Also accumulates per-weight and per-bias gradients (no LR yet).
NISPS_HOT NISPS_FORCE_INLINE
void backprop_accumulate(std::span<const float, FanIn> input,
std::span<const float, FanOut> upstream_err,
std::span<float, FanIn> delta_in,
float sample_weight) noexcept {
for (std::size_t j = 0; j < FanIn; ++j) delta_in[j] = 0.f;
for (std::size_t node = 0; node < FanOut; ++node) {
const float err_signal =
upstream_err[node] * activate_deriv_pre<Act>(pre_activation[node]) * sample_weight;
const std::size_t row = node * FanIn;
for (std::size_t j = 0; j < FanIn; ++j) {
grad_w[row + j] += err_signal * input[j];
delta_in[j] += err_signal * weights[row + j];
}
grad_b[node] += err_signal;
}
}
// Apply accumulated gradient to weights+biases with clipping. Resets
// the accumulators to zero for the next sample/iteration.
NISPS_FORCE_INLINE
void apply_grad(float lr) noexcept {
for (std::size_t i = 0; i < FanIn * FanOut; ++i) {
const float g = clip_gradient(grad_w[i]);
weights[i] -= lr * g;
grad_w[i] = 0.f;
}
for (std::size_t i = 0; i < FanOut; ++i) {
const float g = clip_gradient(grad_b[i]);
biases[i] -= lr * g;
grad_b[i] = 0.f;
}
}
NISPS_FORCE_INLINE
void clear_grad() noexcept {
for (std::size_t i = 0; i < FanIn * FanOut; ++i) grad_w[i] = 0.f;
for (std::size_t i = 0; i < FanOut; ++i) grad_b[i] = 0.f;
}
};
// MLP<NIn, NHidden1, NHidden2, NHidden3, NOut, NMaxExamples = 128>
//
// MaxIterTrain caps the loss-history buffer; if a caller asks for more
// iterations they will be honored at runtime, but only the first
// kMaxIterTrain are recorded for inspection. 4096 fits the playground's
// upper bound and costs 16 KiB.
template <std::size_t NIn,
std::size_t NHidden1,
std::size_t NHidden2,
std::size_t NHidden3,
std::size_t NOut,
std::size_t NMaxExamples = 128u,
std::size_t NMaxIterTrain = 4096u>
class MLP {
public:
static constexpr std::size_t kInput = NIn;
static constexpr std::size_t kHidden1 = NHidden1;
static constexpr std::size_t kHidden2 = NHidden2;
static constexpr std::size_t kHidden3 = NHidden3;
static constexpr std::size_t kOutput = NOut;
static constexpr std::size_t kMaxExamples = NMaxExamples;
static constexpr std::size_t kMaxIterTrain = NMaxIterTrain;
static constexpr std::size_t kNumLayers = 4u;
using Layer0 = Layer<NIn, NHidden1, Activation::ReLU>;
using Layer1 = Layer<NHidden1, NHidden2, Activation::ReLU>;
using Layer2 = Layer<NHidden2, NHidden3, Activation::ReLU>;
using Layer3 = Layer<NHidden3, NOut, Activation::Sigmoid>;
// ---------------------------------------------------------------
// Lifecycle
// ---------------------------------------------------------------
explicit MLP(std::uint64_t seed) noexcept : rng_(seed) {
// Default-init weights with spread=1 (Xavier-like). The IML
// interface caller is expected to draw_weights() with the
// playground spread before the first inference; this default
// simply gives us a non-degenerate starting state for tests
// that skip an explicit draw.
draw_weights(1.f);
clear_dataset_();
loss_history_count_ = 0u;
}
// ---------------------------------------------------------------
// Inference API (concept: set_input / process / outputs)
// ---------------------------------------------------------------
NISPS_FORCE_INLINE void set_input(std::size_t i, float v) noexcept {
if (i < NIn) input_[i] = v;
}
NISPS_HOT void process() noexcept {
forward_(std::span<const float, NIn>(input_));
// Mirror final activation into the output buffer so callers can
// read a stable span.
const auto& a = layer3_.activation;
for (std::size_t i = 0; i < NOut; ++i) output_[i] = a[i];
}
NISPS_FORCE_INLINE std::span<const float> outputs() const noexcept {
return std::span<const float>(output_.data(), NOut);
}
// ---------------------------------------------------------------
// Dataset / Training (concept: add_example / train)
// ---------------------------------------------------------------
// FIFO ring buffer; oldest example evicted when full. No allocation.
// We don't track logical insertion order during training because SGD
// doesn't care — the iteration order over the buffer is arbitrary.
void add_example(std::span<const float> features,
std::span<const float> labels) noexcept {
if (features.size() < NIn || labels.size() < NOut) return;
std::size_t slot;
if (dataset_count_ < NMaxExamples) {
slot = dataset_count_++;
} else {
// Buffer full: overwrite the slot pointed at by head_ (oldest)
// and advance head_ to the next-oldest.
slot = dataset_head_;
dataset_head_ = (dataset_head_ + 1u) % NMaxExamples;
}
const std::size_t f_off = slot * NIn;
const std::size_t l_off = slot * NOut;
for (std::size_t i = 0; i < NIn; ++i) ds_features_[f_off + i] = features[i];
for (std::size_t i = 0; i < NOut; ++i) ds_labels_ [l_off + i] = labels[i];
}
// Concept-required no-arg overload. Default learning rate matches the
// playground's "sane RL training" knob; max_iter and min_err follow.
float train() noexcept {
return train(1.f, 1000u, 0.001f, std::span<const float>{});
}
// Full SGD training. `sample_weights`, if non-empty, must size to the
// current example count and sum to 1.0 (caller's responsibility — we
// do NOT renormalize).
//
// Returns final epoch loss. Records per-iteration loss in
// `loss_history_` (bounded by kMaxIterTrain).
float train(float lr,
std::size_t max_iter,
float min_err,
std::span<const float> sample_weights = {}) noexcept {
loss_history_count_ = 0u;
if (dataset_count_ == 0u) return 0.f;
const bool weighted = !sample_weights.empty();
const float uniform_w = 1.f / static_cast<float>(dataset_count_);
float epoch_loss = 0.f;
for (std::size_t iter = 0; iter < max_iter; ++iter) {
epoch_loss = 0.f;
// SGD: per-sample forward → loss → backprop+update. The order
// is the dataset insertion order; we do not shuffle (matches
// the legacy `Train()` exactly — `TrainBatch` shuffles, but
// we're not implementing batch yet).
for (std::size_t s = 0; s < dataset_count_; ++s) {
const float w = weighted ? sample_weights[s] : uniform_w;
// Forward pass on sample s.
std::span<const float, NIn> x = sample_features_(s);
forward_(x);
// Per-sample loss (NOT scaled by 1/N — the meml-ues fix).
std::array<float, NOut> deriv{};
const float sample_loss = mse_per_sample(
sample_labels_(s),
std::span<const float>(layer3_.activation.data(), NOut),
std::span<float>(deriv.data(), NOut));
// Aggregate weighted loss.
epoch_loss += w * sample_loss;
// Backprop with the same w as the gradient scaler.
backprop_(x, std::span<const float, NOut>(deriv), w);
// Apply gradient (per-sample, SGD).
layer3_.apply_grad(lr);
layer2_.apply_grad(lr);
layer1_.apply_grad(lr);
layer0_.apply_grad(lr);
}
if (loss_history_count_ < NMaxIterTrain) {
loss_history_[loss_history_count_++] = epoch_loss;
}
if (epoch_loss < min_err) break;
}
return epoch_loss;
}
// ---------------------------------------------------------------
// RL ops (concept: move_weights / draw_weights)
// ---------------------------------------------------------------
void move_weights(float speed, float spread,
std::span<const std::uint8_t> output_pin_mask = {}) noexcept {
move_weights_layer(std::span<float>(layer0_.weights), std::span<float>(layer0_.biases),
Layer0::kFanIn, speed, spread, /*final=*/false, {}, rng_);
move_weights_layer(std::span<float>(layer1_.weights), std::span<float>(layer1_.biases),
Layer1::kFanIn, speed, spread, /*final=*/false, {}, rng_);
move_weights_layer(std::span<float>(layer2_.weights), std::span<float>(layer2_.biases),
Layer2::kFanIn, speed, spread, /*final=*/false, {}, rng_);
move_weights_layer(std::span<float>(layer3_.weights), std::span<float>(layer3_.biases),
Layer3::kFanIn, speed, spread, /*final=*/true, output_pin_mask, rng_);
}
void draw_weights(float spread) noexcept {
draw_weights_layer(std::span<float>(layer0_.weights), std::span<float>(layer0_.biases),
Layer0::kFanIn, spread, rng_);
draw_weights_layer(std::span<float>(layer1_.weights), std::span<float>(layer1_.biases),
Layer1::kFanIn, spread, rng_);
draw_weights_layer(std::span<float>(layer2_.weights), std::span<float>(layer2_.biases),
Layer2::kFanIn, spread, rng_);
draw_weights_layer(std::span<float>(layer3_.weights), std::span<float>(layer3_.biases),
Layer3::kFanIn, spread, rng_);
layer0_.clear_grad();
layer1_.clear_grad();
layer2_.clear_grad();
layer3_.clear_grad();
}
// Concept reset: clear weights, dataset, and loss history. Seed is
// intentionally NOT reset (use `seed()` for that).
void reset() noexcept {
clear_dataset_();
loss_history_count_ = 0u;
// Re-init weights from current rng state with default spread.
draw_weights(1.f);
for (std::size_t i = 0; i < NIn; ++i) input_[i] = 0.f;
for (std::size_t i = 0; i < NOut; ++i) output_[i] = 0.f;
}
void seed(std::uint64_t s) noexcept { rng_.seed(s); }
// ---------------------------------------------------------------
// Diagnostics
// ---------------------------------------------------------------
// Forward pass + MSE on a single (input, label) implied by current
// input_ and the most recent training labels — i.e. "what would the
// loss be on the current input if the label were the current output?"
// For now we report the average loss across the training set without
// updating weights. Useful for non-destructive evaluation.
float eval_loss() const noexcept {
if (dataset_count_ == 0u) return 0.f;
const float inv_n = 1.f / static_cast<float>(dataset_count_);
// We need a non-const forward pass to use the cached buffers; since
// eval_loss is logically const, fork a local computation that
// doesn't touch member buffers. That means recomputing through the
// layer weights against scratch arrays — no allocation, just stack.
float total = 0.f;
for (std::size_t s = 0; s < dataset_count_; ++s) {
std::array<float, NHidden1> a1{};
std::array<float, NHidden2> a2{};
std::array<float, NHidden3> a3{};
std::array<float, NOut> ao{};
const std::size_t f_off = s * NIn;
forward_const_layer<NIn, NHidden1, Activation::ReLU>(
std::span<const float, NIn>(ds_features_.data() + f_off, NIn),
layer0_.weights, layer0_.biases, a1);
forward_const_layer<NHidden1, NHidden2, Activation::ReLU>(
std::span<const float, NHidden1>(a1), layer1_.weights, layer1_.biases, a2);
forward_const_layer<NHidden2, NHidden3, Activation::ReLU>(
std::span<const float, NHidden2>(a2), layer2_.weights, layer2_.biases, a3);
forward_const_layer<NHidden3, NOut, Activation::Sigmoid>(
std::span<const float, NHidden3>(a3), layer3_.weights, layer3_.biases, ao);
float sse = 0.f;
const float inv_o = 1.f / static_cast<float>(NOut);
const std::size_t l_off = s * NOut;
for (std::size_t j = 0; j < NOut; ++j) {
const float d = ds_labels_[l_off + j] - ao[j];
sse += d * d * inv_o;
}
total += sse * inv_n;
}
return total;
}
LayerStats layer_stats(std::size_t layer_idx) const noexcept {
switch (layer_idx) {
case 0: return compute_layer_stats(layer0_.weights, layer0_.biases);
case 1: return compute_layer_stats(layer1_.weights, layer1_.biases);
case 2: return compute_layer_stats(layer2_.weights, layer2_.biases);
case 3: return compute_layer_stats(layer3_.weights, layer3_.biases);
default: return {};
}
}
// Flat layout: layer0 weights, layer1 weights, layer2 weights, layer3
// weights, then layer0..3 biases. The split is `weights first all
// layers, then biases all layers` so callers serializing weights can
// pre-compute offsets without consulting layer-specific tables.
static constexpr std::size_t weight_count() noexcept {
return NIn * NHidden1 + NHidden1 * NHidden2 + NHidden2 * NHidden3 + NHidden3 * NOut
+ NHidden1 + NHidden2 + NHidden3 + NOut;
}
// Returns a span into a member-owned scratch buffer that holds a copy
// of the flat weights+biases. The buffer is regenerated on each call,
// so don't hold onto the span across mutations.
std::span<const float> get_weights() noexcept {
std::size_t k = 0u;
// Weights, layer-major.
for (float v : layer0_.weights) flat_weight_buf_[k++] = v;
for (float v : layer1_.weights) flat_weight_buf_[k++] = v;
for (float v : layer2_.weights) flat_weight_buf_[k++] = v;
for (float v : layer3_.weights) flat_weight_buf_[k++] = v;
// Biases.
for (float v : layer0_.biases) flat_weight_buf_[k++] = v;
for (float v : layer1_.biases) flat_weight_buf_[k++] = v;
for (float v : layer2_.biases) flat_weight_buf_[k++] = v;
for (float v : layer3_.biases) flat_weight_buf_[k++] = v;
return std::span<const float>(flat_weight_buf_.data(), weight_count());
}
void set_weights(std::span<const float> w) noexcept {
if (w.size() < weight_count()) return;
std::size_t k = 0u;
for (float& v : layer0_.weights) v = w[k++];
for (float& v : layer1_.weights) v = w[k++];
for (float& v : layer2_.weights) v = w[k++];
for (float& v : layer3_.weights) v = w[k++];
for (float& v : layer0_.biases) v = w[k++];
for (float& v : layer1_.biases) v = w[k++];
for (float& v : layer2_.biases) v = w[k++];
for (float& v : layer3_.biases) v = w[k++];
}
// Run inference on N points (each NIn-sized) and write N output vectors
// (each NOut-sized) into `outs`. NO heap. Modifies the internal cached
// activations as a side effect.
void infer_batch(std::span<const float> points,
std::span<float> outs) noexcept {
const std::size_t n = points.size() / NIn;
if (outs.size() < n * NOut) return;
for (std::size_t i = 0; i < n; ++i) {
const std::size_t in_off = i * NIn;
for (std::size_t j = 0; j < NIn; ++j) input_[j] = points[in_off + j];
process();
const std::size_t out_off = i * NOut;
for (std::size_t j = 0; j < NOut; ++j) outs[out_off + j] = output_[j];
}
}
std::span<const float> loss_history() const noexcept {
return std::span<const float>(loss_history_.data(), loss_history_count_);
}
std::size_t example_count() const noexcept { return dataset_count_; }
void clear_examples() noexcept { clear_dataset_(); }
private:
// ---------------------------------------------------------------
// Internal helpers
// ---------------------------------------------------------------
NISPS_HOT NISPS_FORCE_INLINE
void forward_(std::span<const float, NIn> in) noexcept {
layer0_.forward(in);
layer1_.forward(std::span<const float, NHidden1>(layer0_.activation));
layer2_.forward(std::span<const float, NHidden2>(layer1_.activation));
layer3_.forward(std::span<const float, NHidden3>(layer2_.activation));
}
// Backprop with sample_weight applied to every error signal (so the
// accumulated gradient is already weighted). No weight update happens
// here — caller does it after each sample.
NISPS_HOT NISPS_FORCE_INLINE
void backprop_(std::span<const float, NIn> input,
std::span<const float, NOut> output_deriv,
float sample_weight) noexcept {
std::array<float, NHidden3> d3{};
std::array<float, NHidden2> d2{};
std::array<float, NHidden1> d1{};
std::array<float, NIn> d0{};
layer3_.backprop_accumulate(
std::span<const float, NHidden3>(layer2_.activation),
output_deriv,
std::span<float, NHidden3>(d3),
sample_weight);
layer2_.backprop_accumulate(
std::span<const float, NHidden2>(layer1_.activation),
std::span<const float, NHidden3>(d3),
std::span<float, NHidden2>(d2),
1.f); // weight already in d3
layer1_.backprop_accumulate(
std::span<const float, NHidden1>(layer0_.activation),
std::span<const float, NHidden2>(d2),
std::span<float, NHidden1>(d1),
1.f);
layer0_.backprop_accumulate(
input,
std::span<const float, NHidden1>(d1),
std::span<float, NIn>(d0),
1.f);
}
// const forward pass for diagnostics. Doesn't touch layer caches.
template <std::size_t Fi, std::size_t Fo, Activation A>
static NISPS_FORCE_INLINE void forward_const_layer(
std::span<const float, Fi> in,
const std::array<float, Fi*Fo>& w,
const std::array<float, Fo>& b,
std::array<float, Fo>& out) noexcept {
for (std::size_t node = 0; node < Fo; ++node) {
const std::size_t row = node * Fi;
float sum = b[node];
for (std::size_t j = 0; j < Fi; ++j) sum += w[row + j] * in[j];
out[node] = activate<A>(sum);
}
}
NISPS_FORCE_INLINE std::span<const float, NIn> sample_features_(std::size_t s) const noexcept {
return std::span<const float, NIn>(ds_features_.data() + s * NIn, NIn);
}
NISPS_FORCE_INLINE std::span<const float> sample_labels_(std::size_t s) const noexcept {
return std::span<const float>(ds_labels_.data() + s * NOut, NOut);
}
void clear_dataset_() noexcept {
dataset_count_ = 0u;
dataset_head_ = 0u;
}
// ---------------------------------------------------------------
// Members
// ---------------------------------------------------------------
Layer0 layer0_{};
Layer1 layer1_{};
Layer2 layer2_{};
Layer3 layer3_{};
std::array<float, NIn> input_{};
std::array<float, NOut> output_{};
std::array<float, NMaxExamples * NIn> ds_features_{};
std::array<float, NMaxExamples * NOut> ds_labels_{};
std::size_t dataset_count_ = 0u;
std::size_t dataset_head_ = 0u;
std::array<float, weight_count()> flat_weight_buf_{};
std::array<float, NMaxIterTrain> loss_history_{};
std::size_t loss_history_count_ = 0u;
Rng rng_;
};
// MLP satisfies the MLEngine concept. We keep a static_assert in the test
// suite (test_mlp_concept_satisfied) — see test_mlp_init.cpp.
} // namespace nisps::ml

112
nisps/ml/rl.hpp Normal file
View file

@ -0,0 +1,112 @@
// nisps/ml/rl.hpp — reinforcement-style weight perturbation primitives.
//
// `move_weights` and `draw_weights` are the two playground RL operations
// applied to network weights:
// - `draw_weights(spread)`: re-randomize all weights using the spread-
// aware scale (see init.hpp). This is the "thumbs-down really hard"
// button.
// - `move_weights(speed, spread, output_pin_mask)`: per-call weight decay
// followed by gaussian noise injection. This is the "thumbs-down a
// little" feedback. The pin mask, if supplied, freezes weights feeding
// specific OUTPUT nodes — only on the FINAL layer.
//
// PARITY CONTRACT WITH LEGACY JS PLAYGROUND
// These match `playground/_archive/js/nisps/mlp.js::moveWeights` and
// `drawWeights` semantically. Specifically:
// * Per-layer scale: layer_scale = (1 - spread) + spread / sqrt(fan_in)
// * Weight decay: each weight gets multiplied by (1 - 0.1 * spread)
// BEFORE noise is added. spread=0 ⇒ no decay,
// spread=1 ⇒ ~10% decay per call.
// * Noise: gaussian via rng.next_float_gaussian(speed *
// layer_scale). The legacy code uses a sum-of-three-
// uniforms gaussian, and our `Rng` matches that
// shape (see core/rng.hpp).
// * Output pin mask: 1 byte per output node. If mask[i]==1, all weights
// AND the bias feeding output node i in the final
// layer are skipped. Hidden layers are unaffected by
// the mask.
// * Biases: included in the noise injection, unaffected by
// decay (matching the legacy JS behavior — decay was
// applied to weights only there).
//
// The functions are layer-scoped; the MLP class (mlp.hpp) iterates through
// its layers, calling these for each.
#pragma once
#include <cmath>
#include <cstddef>
#include <cstdint>
#include <span>
#include "../core/rng.hpp"
#include "init.hpp"
namespace nisps::ml {
// Per-call weight-decay coefficient. spread=0 ⇒ 1.0 (no decay), spread=1 ⇒
// 0.9 (10% decay per call). See architecture spread §3 in CLAUDE.md.
inline float weight_decay_factor(float spread) noexcept {
return 1.f - 0.1f * spread;
}
// Apply RL perturbation to one layer's weights+biases. `is_final_layer`
// gates the application of `output_pin_mask`. `output_pin_mask`, if
// non-empty, must be sized fan_out.
//
// Effects, applied in order per weight w:
// w *= (1 - 0.1 * spread)
// w += gaussian(speed * layer_scale)
// Biases:
// b += gaussian(speed * layer_scale)
// Skip rule on final layer if mask[node]==1: leave w/b untouched.
inline void move_weights_layer(std::span<float> weights,
std::span<float> biases,
std::size_t fan_in,
float speed,
float spread,
bool is_final_layer,
std::span<const std::uint8_t> output_pin_mask,
Rng& rng) noexcept {
const float layer_scale = spread_scale(spread, fan_in);
const float noise_stddev = speed * layer_scale;
const float decay = weight_decay_factor(spread);
const std::size_t fan_out = biases.size();
for (std::size_t node = 0; node < fan_out; ++node) {
const bool skip =
is_final_layer && !output_pin_mask.empty() && output_pin_mask[node] != 0u;
// Bias perturbation.
if (!skip) {
biases[node] += rng.next_float_gaussian(noise_stddev);
}
// Per-weight decay+noise. We always advance the RNG even on skip
// so that the random stream is independent of pin-mask state — this
// makes parity tests deterministic regardless of which outputs are
// pinned. Otherwise toggling a pin would shift the entire
// downstream noise sequence, which would be surprising.
const std::size_t row_off = node * fan_in;
for (std::size_t j = 0; j < fan_in; ++j) {
const float noise = rng.next_float_gaussian(noise_stddev);
if (!skip) {
float& w = weights[row_off + j];
w = w * decay + noise;
}
}
}
}
// Re-randomize one layer's weights using the spread-aware scale. Biases are
// reset to zero, matching `init.hpp::spread_init`.
inline void draw_weights_layer(std::span<float> weights,
std::span<float> biases,
std::size_t fan_in,
float spread,
Rng& rng) noexcept {
spread_init(weights, biases, fan_in, spread, rng);
}
} // namespace nisps::ml

76
nisps/ml/stats.hpp Normal file
View file

@ -0,0 +1,76 @@
// nisps/ml/stats.hpp — per-layer weight-health diagnostics.
//
// Mirrors the `nisps_mlp_get_layer_stats` WASM binding semantics: for each
// layer, return mean(|w|), max(|w|), dead fraction, and saturating fraction.
// These are used by the playground "weight health" panel.
//
// Definitions:
// - mean_abs: average of |w| across all weights+biases of the layer
// - max_abs: max of |w| across all weights+biases
// - dead_frac: fraction of weights with |w| < kDeadThresh
// - saturating_frac: fraction with |w| > kSaturatingThresh
//
// Thresholds are chosen to match the JS engine's heuristic targets:
// dead: 0.001 (a weight smaller than this contributes negligibly)
// sat: 3.0 (weights larger than this typically push sigmoid to its
// rails — a smell during normal training but normal during
// RL exploration)
//
// We intentionally include biases in mean/max but not in dead/sat fractions
// — biases are scarce relative to weights and including them would skew the
// fractions, especially at startup when biases are zero.
#pragma once
#include <cmath>
#include <cstddef>
#include <span>
namespace nisps::ml {
inline constexpr float kDeadThresh = 0.001f;
inline constexpr float kSaturatingThresh = 3.0f;
struct LayerStats {
float mean_abs = 0.f;
float max_abs = 0.f;
float dead_frac = 0.f;
float saturating_frac = 0.f;
};
inline LayerStats compute_layer_stats(std::span<const float> weights,
std::span<const float> biases) noexcept {
LayerStats out;
const std::size_t total = weights.size() + biases.size();
if (total == 0u) return out;
float sum = 0.f;
float max_abs = 0.f;
std::size_t dead = 0u;
std::size_t sat = 0u;
for (float w : weights) {
const float a = std::fabs(w);
sum += a;
if (a > max_abs) max_abs = a;
if (a < kDeadThresh) ++dead;
if (a > kSaturatingThresh) ++sat;
}
for (float b : biases) {
const float a = std::fabs(b);
sum += a;
if (a > max_abs) max_abs = a;
}
out.mean_abs = sum / static_cast<float>(total);
out.max_abs = max_abs;
const float w_count = static_cast<float>(weights.size());
if (w_count > 0.f) {
out.dead_frac = static_cast<float>(dead) / w_count;
out.saturating_frac = static_cast<float>(sat) / w_count;
}
return out;
}
} // namespace nisps::ml

37
nisps/ml/training.hpp Normal file
View file

@ -0,0 +1,37 @@
// nisps/ml/training.hpp — gradient-clipping helper and per-layer SGD update.
//
// The MLP class owns the training loop because it knows the dataset layout
// and weight buffers. This header provides:
// - kGradClip: ±10.0, matches the legacy firmware/Layer.h clamp.
// - clip_gradient(): scalar clipper, applied to each accumulated gradient
// before the weight update.
//
// The full training loop (forward, loss, backprop, weight update) is
// implemented inline in mlp.hpp because everything it touches is either a
// member array or layer-templated. Splitting it across translation units
// would require type-erasing the layers, which we don't want.
//
// Optimizer choice: this MVP ships SGD only. RMSProp is planned (the legacy
// firmware uses it for `TrainBatch`) but the playground was using plain SGD
// until very recently and the XOR-convergence benchmark in the test suite
// is the clearer target. RMSProp can land as a follow-up — see the bd
// issue notes. For now `train()` is SGD with optional sample-weight
// scaling and gradient clipping.
#pragma once
#include "../core/perf.hpp"
namespace nisps::ml {
// Per-element gradient clip threshold. Matches the legacy firmware's
// gradientClipValue in Layer.h::ApplyAccumulatedGradients.
inline constexpr float kGradClip = 10.f;
NISPS_FORCE_INLINE float clip_gradient(float g) noexcept {
if (g > kGradClip) return kGradClip;
if (g < -kGradClip) return -kGradClip;
return g;
}
} // namespace nisps::ml

114
tests/cpp/test_mlp_inference.cpp Normal file
View file

@ -0,0 +1,114 @@
// tests/cpp/test_mlp_inference.cpp — verify that the forward pass computes
// the expected linear+activation chain. We construct an MLP with known
// weights via set_weights(), run inference, and compare against a hand-
// computed result.
#include <array>
#include <cmath>
#include "test_helpers.hpp"
#include "../../nisps/ml/mlp.hpp"
namespace {
using TinyMLP = nisps::ml::MLP<2, 2, 2, 2, 1, 8, 32>;
// Manually compute the forward pass using the same activation rules as the
// MLP. ReLU on hidden layers (leaky 0.01), sigmoid on output.
float manual_forward(float x0, float x1,
const std::array<float, 4>& w0, const std::array<float, 2>& b0,
const std::array<float, 4>& w1, const std::array<float, 2>& b1,
const std::array<float, 4>& w2, const std::array<float, 2>& b2,
const std::array<float, 2>& w3, const std::array<float, 1>& b3) {
auto leaky_relu = [](float v) { return v > 0.f ? v : 0.01f * v; };
auto sig = [](float v) {
if (v > 40.f) return 1.f;
if (v < -40.f) return 0.f;
return 1.f / (1.f + std::exp(-v));
};
// Layer 0: 2 → 2, ReLU
float a0_0 = leaky_relu(b0[0] + w0[0] * x0 + w0[1] * x1);
float a0_1 = leaky_relu(b0[1] + w0[2] * x0 + w0[3] * x1);
// Layer 1: 2 → 2, ReLU
float a1_0 = leaky_relu(b1[0] + w1[0] * a0_0 + w1[1] * a0_1);
float a1_1 = leaky_relu(b1[1] + w1[2] * a0_0 + w1[3] * a0_1);
// Layer 2: 2 → 2, ReLU
float a2_0 = leaky_relu(b2[0] + w2[0] * a1_0 + w2[1] * a1_1);
float a2_1 = leaky_relu(b2[1] + w2[2] * a1_0 + w2[3] * a1_1);
// Layer 3: 2 → 1, Sigmoid
return sig(b3[0] + w3[0] * a2_0 + w3[1] * a2_1);
}
NISPS_TEST(mlp_forward_matches_handcomputed) {
TinyMLP m(0ull);
// Set weights to a known pattern: first weight 0.5, second 0.25, etc.
// (we construct a flat buffer matching the documented layout).
constexpr std::size_t WC = TinyMLP::weight_count();
std::array<float, WC> flat{};
// Weights layer by layer:
// L0: 4 weights; L1: 4; L2: 4; L3: 2 → 14 weights
// Biases: L0: 2; L1: 2; L2: 2; L3: 1 → 7 biases
// Total: 21 (sanity check)
NISPS_EXPECT(WC == 21u);
float val = 0.1f;
for (std::size_t i = 0; i < WC; ++i) {
flat[i] = val;
val += 0.05f;
if (val > 0.7f) val = -0.6f;
}
m.set_weights(std::span<const float>(flat));
// Extract layer slices for the manual check.
std::array<float, 4> w0{flat[0], flat[1], flat[2], flat[3]};
std::array<float, 4> w1{flat[4], flat[5], flat[6], flat[7]};
std::array<float, 4> w2{flat[8], flat[9], flat[10], flat[11]};
std::array<float, 2> w3{flat[12], flat[13]};
std::array<float, 2> b0{flat[14], flat[15]};
std::array<float, 2> b1{flat[16], flat[17]};
std::array<float, 2> b2{flat[18], flat[19]};
std::array<float, 1> b3{flat[20]};
const float x0 = 0.3f, x1 = 0.7f;
m.set_input(0, x0);
m.set_input(1, x1);
m.process();
auto out = m.outputs();
NISPS_EXPECT(out.size() == 1u);
const float expected = manual_forward(x0, x1, w0, b0, w1, b1, w2, b2, w3, b3);
NISPS_EXPECT_NEAR(out[0], expected, 1e-5);
}
NISPS_TEST(mlp_forward_outputs_in_unit_range) {
// Sigmoid output guarantees [0, 1].
using M = nisps::ml::MLP<3, 10, 10, 14, 126, 16, 64>;
M m(42ull);
m.draw_weights(0.5f);
for (int t = 0; t < 100; ++t) {
m.set_input(0, static_cast<float>(t) * 0.013f);
m.set_input(1, static_cast<float>(t) * 0.027f - 0.5f);
m.set_input(2, static_cast<float>(t) * 0.041f);
m.process();
for (float v : m.outputs()) {
NISPS_EXPECT(v >= 0.f);
NISPS_EXPECT(v <= 1.f);
}
}
}
NISPS_TEST(mlp_set_input_out_of_range_silently_ignored) {
TinyMLP m(0ull);
// Should not crash; index 99 simply does nothing.
m.set_input(99u, 1.234f);
m.set_input(0, 0.5f);
m.set_input(1, 0.5f);
m.process();
auto out = m.outputs();
NISPS_EXPECT(out.size() == 1u);
}
} // namespace

98
tests/cpp/test_mlp_init.cpp Normal file
View file

@ -0,0 +1,98 @@
// tests/cpp/test_mlp_init.cpp — exercises MLP construction, deterministic
// seeding, the spread parameter, and the static_assert that the class
// satisfies the MLEngine concept.
#include "test_helpers.hpp"
#include "../../nisps/core/concepts.hpp"
#include "../../nisps/ml/mlp.hpp"
namespace {
// Compact alias used across the test suite.
using SmallMLP = nisps::ml::MLP<2, 8, 8, 8, 4, 16, 64>;
// Hard ground-truth: the class satisfies MLEngine. Compile-time check.
static_assert(nisps::MLEngine<SmallMLP>,
"MLP<...> must satisfy nisps::MLEngine concept");
NISPS_TEST(mlp_same_seed_same_init_weights) {
SmallMLP a(123ull);
SmallMLP b(123ull);
auto wa = a.get_weights();
auto wb = b.get_weights();
NISPS_EXPECT(wa.size() == wb.size());
for (std::size_t i = 0; i < wa.size(); ++i) {
NISPS_EXPECT(wa[i] == wb[i]);
}
}
NISPS_TEST(mlp_diff_seed_diff_init_weights) {
SmallMLP a(1ull);
SmallMLP b(2ull);
auto wa = a.get_weights();
auto wb = b.get_weights();
int distinct = 0;
for (std::size_t i = 0; i < wa.size(); ++i) {
if (wa[i] != wb[i]) ++distinct;
}
// With ~200+ weights, almost all should differ.
NISPS_EXPECT(distinct > static_cast<int>(wa.size()) / 2);
}
NISPS_TEST(mlp_seed_method_resets_rng) {
SmallMLP a(5ull);
a.seed(42ull);
a.draw_weights(0.5f);
SmallMLP b(99ull);
b.seed(42ull);
b.draw_weights(0.5f);
auto wa = a.get_weights();
auto wb = b.get_weights();
for (std::size_t i = 0; i < wa.size(); ++i) {
NISPS_EXPECT(wa[i] == wb[i]);
}
}
NISPS_TEST(mlp_weight_count_matches_topology) {
using M = nisps::ml::MLP<3, 10, 10, 14, 126>;
// Layer fan_in*fan_out: 3*10 + 10*10 + 10*14 + 14*126 = 30+100+140+1764 = 2034
// Biases: 10+10+14+126 = 160
// Total: 2194
NISPS_EXPECT(M::weight_count() == 2194u);
}
NISPS_TEST(mlp_spread_zero_uniform_in_minus_one_one_range) {
// spread=0 → weights drawn from U[-1,1]. The max |w| should be ≤ 1.
SmallMLP m(7ull);
m.draw_weights(0.f);
auto w = m.get_weights();
float maxabs = 0.f;
for (float v : w) {
const float a = v >= 0.f ? v : -v;
if (a > maxabs) maxabs = a;
}
NISPS_EXPECT(maxabs <= 1.f);
}
NISPS_TEST(mlp_spread_one_xavier_smaller_range) {
// spread=1 → weights scaled by 1/sqrt(fan_in). For fan_in≥2, weights
// should be strictly smaller in magnitude than the spread=0 case.
SmallMLP m(7ull);
m.draw_weights(1.f);
auto w = m.get_weights();
float maxabs = 0.f;
for (float v : w) {
const float a = v >= 0.f ? v : -v;
if (a > maxabs) maxabs = a;
}
// 1/sqrt(2) ≈ 0.707 — the smallest fan_in is 2 (input). So max possible
// is ≈ 0.707 (drawn from rng_signed, which is < 1).
NISPS_EXPECT(maxabs < 1.f);
}
} // namespace

125
tests/cpp/test_mlp_loss.cpp Normal file
View file

@ -0,0 +1,125 @@
// tests/cpp/test_mlp_loss.cpp — regression test for the meml-ues
// double-scaling bug.
//
// Setup: one training example (label = current network output exactly).
// Expected loss = 0.
//
// Then: one example with label = output + small_delta.
// Expected loss = (1/NOut) * sum(small_delta^2). NOT scaled additionally
// by 1/N (here N=1 so any extra 1/N would still be 1, but the per-element
// derivative path is what would expose double-scaling).
//
// We compare the value reported by `train()` after a single iteration
// (with lr=0 — so no weight update — but the loss should be reported)
// against the expected formula.
#include <array>
#include <cmath>
#include "test_helpers.hpp"
#include "../../nisps/ml/loss.hpp"
#include "../../nisps/ml/mlp.hpp"
namespace {
// Direct test of mse_per_sample.
NISPS_TEST(mse_per_sample_matches_formula) {
std::array<float, 4> label = {0.5f, 0.25f, 0.75f, 0.f};
std::array<float, 4> pred = {0.4f, 0.30f, 0.70f, 0.1f};
std::array<float, 4> deriv{};
const float l = nisps::ml::mse_per_sample(
std::span<const float>(label),
std::span<const float>(pred),
std::span<float>(deriv));
// (0.1^2 + 0.05^2 + 0.05^2 + 0.1^2) / 4 = (0.01 + 0.0025 + 0.0025 + 0.01) / 4 = 0.025/4
const float expected = (0.01f + 0.0025f + 0.0025f + 0.01f) * 0.25f;
NISPS_EXPECT_NEAR(l, expected, 1e-7);
// Derivatives: -(2/N) * (label - pred). For element 0: -(2/4) * (0.5 - 0.4) = -0.05
NISPS_EXPECT_NEAR(deriv[0], -0.05f, 1e-7);
NISPS_EXPECT_NEAR(deriv[3], 0.05f, 1e-7); // -(2/4) * (0 - 0.1) = +0.05
}
NISPS_TEST(mlp_train_loss_not_double_scaled_single_sample) {
// 1 example. The reported epoch_loss should equal the per-sample MSE
// exactly (with sample_weight=1/N=1.0). If the legacy double-scaling
// bug crept in, we'd see it scaled by an extra 1/1 = 1.0 (invisible)
// — so we use TWO samples and check the relationship.
using M = nisps::ml::MLP<1, 2, 2, 2, 2, 4, 8>;
M m(11ull);
m.draw_weights(0.5f);
std::array<float, 1> f1{0.3f};
std::array<float, 2> l1{0.7f, 0.2f};
std::array<float, 1> f2{0.6f};
std::array<float, 2> l2{0.1f, 0.9f};
m.add_example(std::span<const float>(f1), std::span<const float>(l1));
m.add_example(std::span<const float>(f2), std::span<const float>(l2));
// Compute the expected loss manually: average of per-sample MSE.
// Run inference on each, get pred, compute MSE manually.
auto compute_mse = [&](std::array<float, 1>& f, std::array<float, 2>& lab) {
m.set_input(0, f[0]);
m.process();
auto out = m.outputs();
float sse = 0.f;
for (std::size_t j = 0; j < 2u; ++j) {
const float d = lab[j] - out[j];
sse += d * d;
}
return sse * 0.5f; // (1/NOut) sum of squared diff
};
const float mse1 = compute_mse(f1, l1);
const float mse2 = compute_mse(f2, l2);
const float expected_epoch_loss = 0.5f * (mse1 + mse2); // avg over 2 samples
// train() with lr=0 leaves weights unchanged but reports the per-iter
// loss. We use max_iter=1 to grab exactly one epoch_loss.
const float reported = m.train(/*lr=*/0.f, /*max_iter=*/1u, /*min_err=*/-1.f);
NISPS_EXPECT_NEAR(reported, expected_epoch_loss, 1e-5);
// If the legacy bug were present, reported would be expected/2 (extra
// 1/N multiplication). That would fail at 1e-5 tolerance for any
// non-tiny mse. Sanity-confirm:
NISPS_EXPECT(reported > expected_epoch_loss * 0.9f);
}
NISPS_TEST(mlp_train_with_sample_weights_uses_them) {
using M = nisps::ml::MLP<1, 2, 2, 2, 2, 4, 8>;
M m(13ull);
m.draw_weights(0.5f);
std::array<float, 1> f1{0.1f};
std::array<float, 2> l1{0.5f, 0.5f};
std::array<float, 1> f2{0.9f};
std::array<float, 2> l2{0.5f, 0.5f};
m.add_example(std::span<const float>(f1), std::span<const float>(l1));
m.add_example(std::span<const float>(f2), std::span<const float>(l2));
// Concentrate all weight on sample 0.
std::array<float, 2> sw{1.f, 0.f};
auto compute_mse = [&](std::array<float, 1>& f, std::array<float, 2>& lab) {
m.set_input(0, f[0]);
m.process();
auto out = m.outputs();
float sse = 0.f;
for (std::size_t j = 0; j < 2u; ++j) {
const float d = lab[j] - out[j];
sse += d * d;
}
return sse * 0.5f;
};
const float mse1 = compute_mse(f1, l1);
const float reported = m.train(0.f, 1u, -1.f, std::span<const float>(sw));
// Under uniform weights this would be 0.5*(mse1+mse2). With sw={1,0}
// it should equal mse1.
NISPS_EXPECT_NEAR(reported, mse1, 1e-5);
}
} // namespace

133
tests/cpp/test_mlp_rl.cpp Normal file
View file

@ -0,0 +1,133 @@
// tests/cpp/test_mlp_rl.cpp — verify draw_weights and move_weights
// behavior:
// - draw_weights with spread=0 vs spread=1 produces different weight
// magnitude regimes.
// - move_weights perturbs all weights in unmasked layers.
// - move_weights with output_pin_mask preserves the corresponding
// final-layer weight rows AND biases exactly.
#include <array>
#include <cstdint>
#include "test_helpers.hpp"
#include "../../nisps/ml/mlp.hpp"
namespace {
using SmallMLP = nisps::ml::MLP<2, 4, 4, 4, 6, 8, 32>;
NISPS_TEST(mlp_draw_weights_spread_zero_full_range) {
SmallMLP m(0ull);
m.draw_weights(0.f);
auto w = m.get_weights();
float maxabs = 0.f;
for (float v : w) {
const float a = v >= 0.f ? v : -v;
if (a > maxabs) maxabs = a;
}
// spread=0 ⇒ U[-1,1] (no Xavier scale). Most weights will fall in
// (0.5, 1.0) band — we expect at least one above 0.5.
NISPS_EXPECT(maxabs > 0.5f);
}
NISPS_TEST(mlp_draw_weights_spread_one_xavier_compressed) {
SmallMLP m(0ull);
m.draw_weights(1.f);
auto w = m.get_weights();
float maxabs = 0.f;
for (float v : w) {
const float a = v >= 0.f ? v : -v;
if (a > maxabs) maxabs = a;
}
// Xavier scale: 1/sqrt(fan_in). Smallest fan_in = 2 → scale ~0.707.
// Largest weight magnitude bounded by that.
NISPS_EXPECT(maxabs < 0.71f);
}
NISPS_TEST(mlp_move_weights_changes_unpinned_weights) {
SmallMLP m(99ull);
m.draw_weights(0.5f);
// Snapshot weights.
std::array<float, SmallMLP::weight_count()> before{};
{
auto w = m.get_weights();
for (std::size_t i = 0; i < w.size(); ++i) before[i] = w[i];
}
m.move_weights(0.1f, 0.5f);
auto after = m.get_weights();
int distinct = 0;
for (std::size_t i = 0; i < after.size(); ++i) {
if (before[i] != after[i]) ++distinct;
}
// Most weights should have changed (gaussian noise + decay).
NISPS_EXPECT(distinct > static_cast<int>(after.size()) / 2);
}
NISPS_TEST(mlp_move_weights_pin_mask_skips_final_outputs) {
SmallMLP m(33ull);
m.draw_weights(0.5f);
auto w_before = m.get_weights();
std::array<float, SmallMLP::weight_count()> before{};
for (std::size_t i = 0; i < w_before.size(); ++i) before[i] = w_before[i];
// NOut = 6. Pin nodes 0, 2, 4.
std::array<std::uint8_t, 6> mask{1, 0, 1, 0, 1, 0};
m.move_weights(0.1f, 0.5f, std::span<const std::uint8_t>(mask));
auto after = m.get_weights();
// Layout: weights = [L0(8) L1(16) L2(16) L3(24)] then biases [L0(4) L1(4) L2(4) L3(6)].
constexpr std::size_t L0_W = 2*4; // 8
constexpr std::size_t L1_W = 4*4; // 16
constexpr std::size_t L2_W = 4*4; // 16
constexpr std::size_t L3_W = 4*6; // 24
constexpr std::size_t L3_W_OFF = L0_W + L1_W + L2_W;
constexpr std::size_t BIAS_OFF = L0_W + L1_W + L2_W + L3_W;
constexpr std::size_t L3_B_OFF = BIAS_OFF + 4 + 4 + 4; // 76 + 0 → biases start
// Final layer weights row-major: [node*4 + j] for j ∈ [0,4).
// Pinned nodes 0, 2, 4 → rows 0, 2, 4 should be preserved.
for (std::size_t node : {0u, 2u, 4u}) {
for (std::size_t j = 0; j < 4u; ++j) {
const std::size_t idx = L3_W_OFF + node * 4u + j;
NISPS_EXPECT(before[idx] == after[idx]);
}
// Bias too.
NISPS_EXPECT(before[L3_B_OFF + node] == after[L3_B_OFF + node]);
}
// Unpinned nodes 1, 3, 5 → rows changed (most weights distinct).
int unpinned_changed = 0;
for (std::size_t node : {1u, 3u, 5u}) {
for (std::size_t j = 0; j < 4u; ++j) {
const std::size_t idx = L3_W_OFF + node * 4u + j;
if (before[idx] != after[idx]) ++unpinned_changed;
}
}
// 12 unpinned weights; gaussian noise w/ stddev > 0 → almost all change.
NISPS_EXPECT(unpinned_changed >= 10);
}
NISPS_TEST(mlp_draw_weights_clears_grad_accumulators) {
// After draw_weights, calling train() shouldn't see stale gradients.
// Smoke test: draw → train → loss should drop normally.
using M = nisps::ml::MLP<2, 4, 4, 4, 1, 4, 64>;
M m(0ull);
m.draw_weights(0.5f);
std::array<float, 2> f{0.3f, 0.7f};
std::array<float, 1> l{0.5f};
m.add_example(std::span<const float>(f), std::span<const float>(l));
const float l1 = m.train(0.5f, 50u, -1.f);
m.draw_weights(0.5f);
const float l2 = m.train(0.5f, 50u, -1.f);
// Both should be finite and >= 0. We don't assert about ordering; the
// point is that no NaNs leak through stale grads.
NISPS_EXPECT(l1 >= 0.f && l1 < 100.f);
NISPS_EXPECT(l2 >= 0.f && l2 < 100.f);
}
} // namespace

146
tests/cpp/test_mlp_serialize.cpp Normal file
View file

@ -0,0 +1,146 @@
// tests/cpp/test_mlp_serialize.cpp — round-trip get_weights → set_weights.
//
// Serializing weights, restoring them in a fresh MLP, and running inference
// must produce bit-identical outputs. This is the gate for cross-platform
// state transfer (firmware ↔ browser).
#include <array>
#include <vector>
#include "test_helpers.hpp"
#include "../../nisps/ml/mlp.hpp"
namespace {
using TestMLP = nisps::ml::MLP<3, 8, 8, 10, 5, 16, 32>;
NISPS_TEST(mlp_get_set_weights_roundtrip_preserves_inference) {
TestMLP a(123ull);
a.draw_weights(0.4f);
// Snapshot weights.
auto w = a.get_weights();
std::vector<float> snap(w.begin(), w.end());
TestMLP b(0ull); // different seed → different starting weights
b.set_weights(std::span<const float>(snap.data(), snap.size()));
// Inference on identical input should produce identical outputs.
const float in[3] = {0.2f, -0.3f, 0.5f};
for (std::size_t i = 0; i < 3u; ++i) {
a.set_input(i, in[i]);
b.set_input(i, in[i]);
}
a.process();
b.process();
auto oa = a.outputs();
auto ob = b.outputs();
NISPS_EXPECT(oa.size() == ob.size());
for (std::size_t i = 0; i < oa.size(); ++i) {
NISPS_EXPECT(oa[i] == ob[i]);
}
}
NISPS_TEST(mlp_weight_count_matches_get_weights_size) {
TestMLP m(0ull);
auto w = m.get_weights();
NISPS_EXPECT(w.size() == TestMLP::weight_count());
}
NISPS_TEST(mlp_get_weights_contains_layer_concatenation) {
// Verify the documented flat layout: weights all layers, then biases
// all layers. We do this by setting a known pattern via set_weights
// and reading it back.
TestMLP m(0ull);
constexpr std::size_t WC = TestMLP::weight_count();
std::vector<float> pattern(WC);
for (std::size_t i = 0; i < WC; ++i) {
pattern[i] = static_cast<float>(i) * 0.001f - 0.5f;
}
m.set_weights(std::span<const float>(pattern.data(), pattern.size()));
auto w = m.get_weights();
NISPS_EXPECT(w.size() == WC);
for (std::size_t i = 0; i < WC; ++i) {
NISPS_EXPECT_NEAR(w[i], pattern[i], 1e-7);
}
}
NISPS_TEST(mlp_set_weights_too_short_is_ignored) {
TestMLP m(7ull);
auto before = m.get_weights();
std::vector<float> snap(before.begin(), before.end());
// Pass a too-short buffer.
std::array<float, 2> tiny{1.f, 2.f};
m.set_weights(std::span<const float>(tiny));
auto after = m.get_weights();
// Weights should be unchanged.
for (std::size_t i = 0; i < snap.size(); ++i) {
NISPS_EXPECT(snap[i] == after[i]);
}
}
NISPS_TEST(mlp_eval_loss_does_not_modify_state) {
TestMLP m(5ull);
m.draw_weights(0.5f);
std::array<float, 3> f{0.1f, 0.2f, 0.3f};
std::array<float, 5> l{0.5f, 0.5f, 0.5f, 0.5f, 0.5f};
m.add_example(std::span<const float>(f), std::span<const float>(l));
// Snapshot weights.
auto w_before = m.get_weights();
std::vector<float> snap(w_before.begin(), w_before.end());
const float L = m.eval_loss();
NISPS_EXPECT(L >= 0.f);
auto w_after = m.get_weights();
for (std::size_t i = 0; i < snap.size(); ++i) {
NISPS_EXPECT(snap[i] == w_after[i]);
}
}
NISPS_TEST(mlp_infer_batch_matches_individual_inference) {
TestMLP m(42ull);
m.draw_weights(0.5f);
constexpr std::size_t N = 5;
constexpr std::size_t NI = 3;
constexpr std::size_t NO = 5;
std::array<float, N * NI> points{};
for (std::size_t i = 0; i < N * NI; ++i) {
points[i] = static_cast<float>(i) * 0.07f - 0.3f;
}
std::array<float, N * NO> outs{};
m.infer_batch(std::span<const float>(points), std::span<float>(outs));
// Compare against individual inference.
for (std::size_t i = 0; i < N; ++i) {
for (std::size_t j = 0; j < NI; ++j) {
m.set_input(j, points[i * NI + j]);
}
m.process();
auto o = m.outputs();
for (std::size_t j = 0; j < NO; ++j) {
NISPS_EXPECT(o[j] == outs[i * NO + j]);
}
}
}
NISPS_TEST(mlp_layer_stats_reasonable_after_init) {
TestMLP m(0ull);
m.draw_weights(0.5f);
for (std::size_t L = 0; L < 4; ++L) {
const auto s = m.layer_stats(L);
NISPS_EXPECT(s.mean_abs > 0.f);
NISPS_EXPECT(s.max_abs >= s.mean_abs);
NISPS_EXPECT(s.dead_frac >= 0.f && s.dead_frac <= 1.f);
NISPS_EXPECT(s.saturating_frac >= 0.f && s.saturating_frac <= 1.f);
}
}
} // namespace

92
tests/cpp/test_mlp_training.cpp Normal file
View file

@ -0,0 +1,92 @@
// tests/cpp/test_mlp_training.cpp — convergence test on XOR.
//
// XOR is the classic minimum non-linear problem: a 2-layer linear net
// cannot solve it; an MLP with one hidden layer (and a non-linear
// activation) can. Our 4-layer MLP with sigmoid output is more than enough.
//
// We check loss < 0.01 within a generous iteration budget. If this test
// regresses to taking >1000 iterations, something is wrong with the
// gradient or weight-update path.
#include <array>
#include "test_helpers.hpp"
#include "../../nisps/ml/mlp.hpp"
namespace {
NISPS_TEST(mlp_xor_converges) {
// Modest network: 2 inputs, [4, 4, 4] hidden, 1 output.
using M = nisps::ml::MLP<2, 4, 4, 4, 1, 4, 1024>;
M m(7ull);
m.draw_weights(1.f); // Xavier-ish; needed for sigmoid output to start reasonable
// XOR truth table.
std::array<std::array<float, 2>, 4> X = {{
{0.f, 0.f},
{0.f, 1.f},
{1.f, 0.f},
{1.f, 1.f},
}};
std::array<std::array<float, 1>, 4> Y = {{
{0.f},
{1.f},
{1.f},
{0.f},
}};
for (std::size_t i = 0; i < 4u; ++i) {
m.add_example(std::span<const float>(X[i]), std::span<const float>(Y[i]));
}
NISPS_EXPECT(m.example_count() == 4u);
// Train. Higher LR + more iterations is fine — the test is "did it
// converge AT ALL within a generous budget".
const float final_loss = m.train(/*lr=*/0.5f, /*max_iter=*/2000u, /*min_err=*/0.01f);
NISPS_EXPECT(final_loss < 0.01f);
// Sanity: outputs should be near labels for each input.
for (std::size_t i = 0; i < 4u; ++i) {
m.set_input(0, X[i][0]);
m.set_input(1, X[i][1]);
m.process();
const float o = m.outputs()[0];
NISPS_EXPECT_NEAR(o, Y[i][0], 0.2);
}
// Loss history should record at least one entry.
NISPS_EXPECT(m.loss_history().size() >= 1u);
}
NISPS_TEST(mlp_train_with_no_examples_returns_zero) {
using M = nisps::ml::MLP<2, 4, 4, 4, 1, 4, 8>;
M m(0ull);
const float loss = m.train(0.5f, 100u, 0.001f);
NISPS_EXPECT(loss == 0.f);
}
NISPS_TEST(mlp_clear_examples_works) {
using M = nisps::ml::MLP<2, 4, 4, 4, 1, 4, 8>;
M m(0ull);
std::array<float, 2> f{0.f, 1.f};
std::array<float, 1> l{0.5f};
m.add_example(std::span<const float>(f), std::span<const float>(l));
NISPS_EXPECT(m.example_count() == 1u);
m.clear_examples();
NISPS_EXPECT(m.example_count() == 0u);
}
NISPS_TEST(mlp_dataset_ring_buffer_evicts_oldest) {
// NMaxExamples=4, add 6 examples; oldest 2 should be evicted.
using M = nisps::ml::MLP<1, 2, 2, 2, 1, 4, 8>;
M m(0ull);
for (int i = 0; i < 6; ++i) {
std::array<float, 1> f{static_cast<float>(i)};
std::array<float, 1> l{static_cast<float>(i) * 0.1f};
m.add_example(std::span<const float>(f), std::span<const float>(l));
}
NISPS_EXPECT(m.example_count() == 4u);
}
} // namespace