fastlsq v0.1.5 · scalar & vector-valued estimator

Solve PDEs
in one shot.

FastLSQ skips neural-net training entirely. Random Fourier features have a closed-form cyclic derivative, so a single least-squares call reaches 10⁻⁷ relative L2 error in 0.07 s on linear PDEs. Orders of magnitude faster than PINNs.

$ pip install fastlsq
view source read paper
§ 01 · live demo

Random sinusoids,
fit in one solve.

N features → least-squares → uₙ(x)
N 4 ‖uₙ − u‖L²/‖u‖ · solve ·
target u(x) fit uₙ(x) features φⱼ
features N σ = 8.0
§ 02 · method

Least-squares on a
random Fourier basis.

u(x) ≈ Σⱼ βⱼ · sin(Wⱼ·x + bⱼ)

01 · the basis

The unknown solution is expanded in a fixed basis of random sinusoids: φⱼ(x) = sin(Wⱼ·x + bⱼ), with Wⱼ and bⱼ sampled once and never updated.

02 · why it works

Averaging enough random Fourier features recovers the Gaussian RBF kernel (Bochner). Smooth solutions live in that RKHS, so a few hundred sinusoids cover the space.

03 · the fit

Find β by a single least-squares call. Scalar or vector-valued targets: components stack into the same block-LSQ, one solve.

§ 03 · API

Compose any operator.
One formula handles it.

helmholtz_2d.py
helmholtz_2d.py · k = 10, [0, 1]² 
import numpy as np
import fastlsq as fl
from fastlsq.basis import SinusoidalBasis, Op

# Helmholtz 2D:  Δu + k² u = f   on  [0, 1]²
basis = SinusoidalBasis(d=2, N=1500, sigma=10.0)
L     = Op.laplacian(d=2) + (10.0**2) * Op.identity(d=2)

u, info = fl.solve(
    L      = L,
    f      = lambda x: 0.0,
    bc     = lambda x: np.sin(np.pi * x[..., 0]),
    domain = fl.Box(d=2),
    basis  = basis,
)

print(info.l2_error)   # 1.9e-6
print(info.wall_clock) # 0.08 s
status solved ✓ rel L² 1.9 × 10⁻⁶ wall-clock 0.08 s features 1500 scikit-fem P2 baseline 4.0 × 10⁻⁵
§ 04 · examples

Two physical demos.
One basis, one solve, two very different problems.

synchrotron radiation · 9-mirror beam shaper
Animation: two point charges orbiting a ring at β=0.70, with Liénard–Wiechert scalar potential and Doppler-beamed emission lobes
Liénard–Wiechert potential · Doppler-beamed lobes · past-emission wavefronts β = 0.70

Synchrotron radiation.
The whole field, in one solve.

Two point charges orbit a ring in opposite directions at 70% the speed of light. Each emits light in every direction, but because the charge is moving, the light is bunched into the forward direction. The bright lobes are this Doppler beaming; the faint spirals are photon shells emitted at earlier times, expanding outward at c.

The hard part is the retarded time: an observer at (r, t) sees the charge not at its present position but where it was when the light left. That gives an implicit equation with no closed form, normally solved per pixel by a root finder.

fastlsq instead fits the entire scalar field τ(r, t) as one random-Fourier expansion. Coefficients β come from one least-squares call against 25 000 root-finder samples. After that, τ and all of its space–time derivatives are available in closed form at any query point, so evaluating the Liénard–Wiechert potential on a dense grid is essentially free.

Animation: 9-mirror cascade reshaping an aberrated Gaussian into a ring, joint Gauss–Newton convergence
9-mirror cascade · per-mirror phases φk · RMS convergence 40 GN iters · RMS 8.6e-4

Inverse-design beam shaper.
Sculpt nine mirrors at once.

Nine corner mirrors are folded into an L/staircase layout. An aberrated Gaussian enters; the goal is a clean ring at the output. Each mirror has a known textbook starting shape (axicon + parabolic relays) plus a free-form correction and two tilt angles.

The free-form correction lives in fastlsq's random-Fourier feature subspace, so every recovered surface is smooth and manufacturable by construction, no spline knots, no Zernike truncation. Tilts ride on the same coefficient vector by appending two columns to each mirror's basis matrix.

Differentiating the FFT cascade analytically gives a block Jacobian for every mirror in one batched FFT. Stack the blocks, hand to solve_lstsq with Levenberg–Marquardt damping, and one joint Gauss–Newton step updates every degree of freedom on every mirror.

§ 05 · cite

The paper.

Full method, proofs, ablations, and 17-PDE benchmark are described in Sulc, A. Solving PDEs in One Shot via Fourier Features with Exact Analytical Derivatives, arXiv:2602.10541.

Code, examples, and the inverse-problem demos are available in the reference implementation on GitHub.

bibtex 
@misc{sulc2026solvingpdesshotfourier,
  title         = {Solving PDEs in One Shot via Fourier Features
                   with Exact Analytical Derivatives},
  author        = {Antonin Sulc},
  year          = {2026},
  eprint        = {2602.10541},
  archivePrefix = {arXiv},
  primaryClass  = {math.NA},
  url           = {https://arxiv.org/abs/2602.10541},
}