Project: Formula Discovery Without Search

The symbolic-regression engine finds formulas by trying four million of them. This page finds them without trying any: one fixed three-parameter function family contains $\log (1+x)$ and $x/(1+x)$ as specific parameter settings, so "which function is this?" becomes gradient descent on three numbers. That is the trick behind symbolic metamodeling (Alaa & van der Schaar, NeurIPS 2019), built here on Meijer G-functions. Same Knot → C → WebAssembly pipeline as the SCF; the gradients are digamma functions, hand-rolled, and the whole engine is 58 KB.

One family, many functions

The Meijer G-function $G_{2,2}^{1,2}(a,a;a,b\mid cx)$ is a two-line contour integral in general, but the engine only ever needs its truncated power series — four terms whose coefficients are ratios of Gamma functions of $a-b$:

f(a, b, c; x) = (cx)^a [ 1/Γ(a−b+1) − cx/Γ(a−b+2) + 2(cx)²/Γ(a−b+3) − 6(cx)³/Γ(a−b+4) ]

f(1, 0, 1; x)  =  x − x²/2 + x³/3 − x⁴/4   =  log(1 + x)   (truncated)
f(1, 1, 1; x)  =  x (1 − x + x² − x³)      =  x/(1 + x)    (truncated)

Read those two rows again: the same expression, at two different parameter points, reproduces the series of two entirely different elementary functions. The family is a curved surface through function space with familiar functions sitting at special coordinates. Fitting data with it turns symbolic discovery into ordinary optimization — no expression trees, no enumeration, no combinatorics. The discrete question is deferred to the end, where you ask which special coordinate the fit landed near.

Gradients through the Gamma function

Descending in $(a,b,c)$ needs $\partial f/\partial b$, and $b$ lives inside $\Gamma (a-b+k+1)$. The derivative of $\log \Gamma$ is the digamma function $\psi$, so the gradient formulas are the value formulas with $\psi$ sprinkled where $\Gamma$ stood — closed form, no autodiff:

# d/db of the basis -- digamma where the value has Gamma, because
# d/dz log Gamma(z) = psi(z). Transliterated from the paper's basis_grad,
# verified against scipy to 1e-8 (see the tests at the bottom of the source).
out[2] = p * (0.0 - 6.0 * cx3 * K4 / G4 + 2.0 * cx2 * K3 / G3 - cx * K2 / G2 + K1 / G1)

Knot has neither $\Gamma$ nor $\psi$, so the source (tab ① above) builds both from scratch: Lanczos approximation for Gamma, shift-and-asymptotic-series for digamma, each verified against scipy references in the tests at the bottom of the file. The optimizer is Adam, written out as six lines of moment bookkeeping — the loss surface has a long flat valley where $b$ and $c$ trade off against each other, and raw gradient descent crawls in it.

Identification by snapping, not by reading off parameters

That flat valley is why the final answer is not read off the fitted parameters. On a short interval many $(a,b,c)$ triples draw nearly the same curve, so the fit for $\log (1+x)$ ends somewhere like $(1.00,0.12,0.96)$ — near the catalog point, not on it. The identifiable thing is function-space distance: snap the parameters to each catalog point and re-measure the error. Snapping to $(1,0,1)$ keeps the fit tight; snapping to $(1,1,1)$ makes it 100× worse. Verdict: it's the log. And when the target is $\arctan (x)$ — a function this family does not contain — every snap ruins the fit, and the engine says so. That is the method's failure mode, reported instead of hidden: a fitted G-function that lands away from every reducible point is a hypergeometric with no familiar name.

The two engines make a matched pair. Brute-force enumeration searches a discrete space and decides which function during the search; metamodeling fixes one continuous family and defers the which-function question to a snap at the end. Discrete search pays in combinatorics; continuous relaxation pays in hypothesis-class rigidity. Every trick in the two papers is a way of managing one cost or the other.

Reproduce it

scripts/build_mg_wasm.sh rebuilds the pipeline; the Knot source carries its own test suite (Gamma and digamma against scipy references, the gradient transliteration against the paper's code at a spot-check point, and the two catalog identities in exact arithmetic). Every number above comes out of the wasm running on your machine.