fem-surrogate-toolkit

A small, dependency-light toolkit for building fast surrogate models of expensive simulations: Latin-hypercube design of experiments, POD/PCA output compression, ridge regression with polynomial features, and an end-to-end Surrogate that wires them together. Pure NumPy.

The idea is the one behind every practical simulation surrogate: run the expensive solver at a well-chosen handful of designs, learn the map from parameters to response, and then evaluate thousands of new designs in milliseconds. This toolkit is the honest baseline version of that pipeline, the thing to build and beat before reaching for a neural network.

The pipeline

parameters ──▶ StandardScaler ──▶ Ridge(+polynomial) ──▶ POD⁻¹ ──▶ response

• latin_hypercube — stratified DoE sampling of the parameter box, so a few hundred solver runs cover the space evenly.
• POD — proper orthogonal decomposition (PCA via SVD) that compresses a whole field or spectrum to a handful of mode coefficients.
• Ridge — closed-form ridge regression, optionally on polynomial features, mapping parameters to those coefficients.
• Surrogate — the whole thing behind one fit / predict.

Install

pip install -e .

Python 3.9+, NumPy. Matplotlib only for the example.

Use

from femsurrogate import latin_hypercube, Surrogate, r2_score

bounds = [[0.5, 2.0], [0.02, 0.2], [0.2, 0.8]]   # one [low, high] per parameter
X = latin_hypercube(bounds, n_samples=250, seed=0)
Y = expensive_solver(X)                           # (n_samples, n_outputs)

sm = Surrogate(alpha=1e-6, degree=4, pod_variance=0.999).fit(X, Y)
Y_new = sm.predict(new_designs)                   # milliseconds for thousands

examples/spectrum_surrogate.py builds a surrogate for a 128-point response curve. On unseen designs it retains ~24 POD modes, scores R² ≈ 0.95, and predicts at about 3 µs per design:

Notice the surrogate nails the broad response but rounds off the sharp resonance peaks. That is the honest signature of a linear-basis surrogate: sharp, parameter-shifting features are exactly where it struggles, and exactly where a physics-informed or convolutional model earns its keep. Knowing that boundary is the point.

Tests check the maths, not just execution

python -m pytest tests -q

• Latin hypercube puts exactly one sample in every 1/n stratum on every axis, is reproducible from its seed, and centres on stratum midpoints when asked
• the scaler round-trips and survives a constant (held-fixed) parameter column
• POD reconstructs genuinely low-rank data to machine precision, and its variance criterion picks the minimal rank
• ridge recovers a known linear map as alpha → 0, and its penalty shrinks the coefficients as alpha grows
• the end-to-end surrogate generalises to held-out designs at R² > 0.98 on a smooth analytic response, and uses POD to hit the same bar on a vector field

Limits, stated plainly

• Linear-basis surrogate. Ridge on polynomial features plus a linear POD basis. It handles smooth response surfaces well and sharp, translating features (moving resonance peaks, discontinuities) poorly. It is a baseline, not a universal approximator.
• Dense, in-memory. Designed for hundreds to low-thousands of training designs and up to a few thousand output points. Not a big-data trainer.
• No uncertainty estimate. It returns a point prediction, not error bars. Gaussian-process or ensemble methods are the next step if you need those.
• No adaptive sampling. The DoE is one-shot; active/sequential design is not implemented here.
• NumPy only. No PyTorch/ONNX export in the core. Wiring the fitted linear maps into an ONNX graph is straightforward and left as an extension.

Licence

MIT.