Module 11 of 14 · one week
A simulator can make training data on demand, in a curriculum you control. An equation can score an output with no labels at all. The seed episode calls the second of these "dense feedback", and it is a real mechanism with published numbers. It is also an optimization problem whose conditioning gets worse exactly where you would want to use it, and this module derives how much worse.
Between 15:00 and 18:30 the seed episode says the model improves itself because physics laws give dense feedback, and the company site says the same in different words: known laws validate outputs and provide "directional feedback", reducing the need for labels (research_notes/company-episode.md:112).
The mechanism behind that is public and small. A PDE residual is a number you can compute from any candidate output, with no reference solution, so it is a loss you can evaluate anywhere. That much is not in dispute. What is in dispute is whether it scales, and the answer turns out to depend on a property of the residual loss that has nothing to do with neural networks: how badly conditioned it is. This module gives the published forms of both uses, the numbers they reach, and then works out the conditioning for this course's own thermal problem so the caveat is concrete.
Krishnapriyan et al.'s fix from Module 4 is the simplest form. Train at an easy coefficient, then warm-start the next harder one from those weights. On convection at \(\beta=30\), a standard PINN gives 8.97e-1 relative error and the curriculum gives 2.02e-2, a factor of 44 (web-papers.md:620-626).
Notice what the "data engine" is here: the simulator is trivial, since the exact solution \(\sin(x-\beta t)\) is known in closed form, and the curriculum runs over the equation rather than over the data. That is the same structure Anandkumar describes as composing learned sub-phenomena, heat first, then stretching, then the coupled problem (earlier episode 01:14:20 to 01:15:48). The idea does not need a solver at all; it needs a family of problems ordered by difficulty and a way to move along it.
George, Zhao, Kossaifi, Li and Anandkumar (arXiv 2211.15188) run the curriculum along a different axis: grow the number of Fourier modes and the training resolution during training instead of fixing them. They report 10 percent lower test error with 20 percent fewer modes and 30 percent faster training than a fixed FNO (lineage.md:104-113). In neuraloperator this is IncrementalFNOTrainer in neuralop/training/incremental.py, with an example at examples/training/plot_incremental_FNO_darcy.py (repos.md:63, :99).
Module 2's frequency-doubling argument explains why growing modes late is cheap. Early in training the model's error lives in the low modes, where the signal's energy is; the high modes carry a small fraction of the total and the nonlinearity can regenerate part of them anyway. Spending parameters on modes the model cannot yet use is spending them on the wrong thing. Lab 03's measurement in Module 9 says the same from the cost side: going from 16 to 64 modes multiplies the parameters by 14.4 while adding 6 percent to activation memory, so modes are almost pure storage cost, and deferring that cost saves real time.
Both labs in this course are data engines. Lab 01 generates 1200 Burgers trajectories in 120 seconds; Lab 02 generates 1200 thermal dies in 14 seconds. Each one runs a self-check that halts the script if the solver is wrong, and that check is the entire difference between a data engine and a noise engine.
PDEBench's Appendix H is the standing warning, and Module 8 gives it in full: a 128×128 "ground truth" solver carried discretization error of 1.8e-1 against a 512×512 reference, which is larger than the FNO's 8.1e-3 error against the same reference (pdf-mgn-pdebench.md:286). If the data engine is coarser than the model it feeds, the model is being graded against the wrong answer, and every downstream number is about the solver, not the model.
The loss is \(\mathcal J = \mathcal J_{\text{data}} + \lambda\mathcal J_{\text{pde}}\), with the residual computed by exact Fourier differentiation of the operator's output on a grid where no data exists (web-papers.md:466-479). Module 4 has the full numbers; the two that matter here are that PINO on Burgers reaches 1.22 percent with data and 1.50 percent with no paired data at all, against DeepONet's 6.97 percent, and that on Darcy an 11×11 supervision grid plus a 61×61 residual reaches 1.56 percent at four times the training resolution against 9.46 percent for the same data without the residual. The residual is doing real work, and it is doing it precisely where labels are missing.
Majumdar et al. (2023) close the loop: iterative pseudo-labelling of a PINO with the physics loss as the only signal. Their result is reported as a ratio to fully supervised training, and it is close: 1.07 times the supervised error on Burgers, 1.02 times on Darcy (foundation.md:156-159). That is the published form of "self-improvement", and it works. It is also two small PDEs, with nothing at foundation-model scale, which is the gap Module 12 grades.
The recent work reframes the residual as a reward signal instead of a loss term. PIRF (arXiv 2509.20570) uses local finite-difference residuals to fine-tune diffusion models for PDE generation, updating only the high-resolution layers because the reward is local. arXiv 2509.02846 uses a reward model to guide inference-time compute for stochastic PDE predictors, motivated explicitly by the admission that pretrained models "struggle with auto-regressive rollout performance, especially in out-of-distribution cases". And arXiv 2603.15431 fine-tunes PDE foundation models with residual and boundary losses, reporting that a hybrid data-plus-physics strategy generalizes best under minimal data. Exact numbers were not extracted from any of the three; these are abstract-level readings, graded B (foundation.md:161-167).
Every residual term inherits the Krishnapriyan result, and the derivation beats a quotation here, because it says exactly how the problem scales with the grid.
The setting. Take Lab 02's operator: \(-kt\nabla^2T = q\) on a square with \(T=0\) on the edge. Discretize with the five-point stencil on an \(N\times N\) grid of spacing \(h\) and call the resulting matrix \(L\). A pure residual loss is \(\mathcal J(T) = \|LT-q\|_2^2\), whose Hessian in \(T\) is \(2L^\top L\). The conditioning of that Hessian is what a gradient method has to fight.
The eigenvalues of \(L\). Start in one dimension, where the Dirichlet second-difference matrix has eigenvalues \[\lambda_j = \frac{4}{h^2}\sin^2\!\Big(\frac{j\pi}{2(N+1)}\Big), \qquad j = 1,\dots,N .\] At \(j=1\) the sine argument is small, so \(\sin^2x\approx x^2\) gives \(\lambda_{\min}\approx\pi^2\) with \(h=1/(N+1)\). At \(j=N\) the argument approaches \(\pi/2\), so \(\lambda_{\max}\to4/h^2 = 4(N+1)^2\). The two-dimensional five-point stencil is the sum of two such operators acting on different indices, so its eigenvalues are \(\lambda_i+\lambda_j\) and both the smallest and the largest double. The ratio is therefore unchanged: \[\mathrm{cond}(L) = \frac{\lambda_{\max}}{\lambda_{\min}} \approx \frac{4(N+1)^2}{\pi^2} \approx 0.405\,N^2 ,\] in one dimension and in two. The conditioning of the loss. Squaring the operator squares the condition number: \(\mathrm{cond}(L^\top L) \approx 0.164\,N^4\).
The numbers for this course's grids. The asymptotic form gives \(\mathrm{cond}(L)\approx 1.7\times10^3\) and \(\mathrm{cond}(L^\top L)\approx 2.8\times10^6\) at \(N=64\), and \(6.6\times10^3\) and \(4.4\times10^7\) at \(N=128\). The exact ratio of the two sine values, without the small-angle approximation, gives \(1.71\times10^3\) and \(2.93\times10^6\) at \(N=64\), and \(6.74\times10^3\) and \(4.55\times10^7\) at \(N=128\), so the asymptotic form is low by about 6 percent at 64 and 3 percent at 128. Both sets are computed on this page from the closed-form eigenvalues; no matrix was built and no lab file exists for them. Every doubling of the grid multiplies the loss's condition number by 16.
Set that against the promise. "Dense feedback" is dense in space and in time, and it is genuinely label-free. It is also a quadratic form whose condition number grows like the fourth power of the resolution, which means the finer the grid and the harder the physics, the worse the optimization problem, and a foundation model targets exactly the fine-grid, hard-physics regime. Krishnapriyan et al. report the same scaling from the other direction, \((\beta N)^2\) for convection and \((\nu N^2)^2\) for diffusion (web-papers.md:613-616), and the \(N^4\) above is the second of those with \(\nu\) absorbed.
What the field does about it, and what it does not. PINO's answer is to put an operator in front, so the residual only fine-tunes an already-good solution and never has to find one from scratch. Wang, Wang and Perdikaris' answer is a hand-tuned weight, swept over \(\{1,5,10,20,50,100\}\) and set to 20 (:565-573). Neither is a theory of the weight, and the JMLR paper lists "physics-informed operator learning" as a direction and not a result (pdf-fno-jmlr.md:719-724). The library ships the pieces, BurgersEqnLoss, PoissonEqnLoss and WeightedSumLoss, and it also ships SoftAdapt and Relobralo in neuralop/losses/ (repos.md:57-58). Those last two are adaptive weighting schemes, and the fact that a library needs two of them is itself evidence about how settled the question is.
One more caveat, carried over from Module 4 and sharpened by Module 8. Automatic differentiation of a model's output is not evidence of physical fidelity. A residual can be small on a wrong field if the field is smooth enough, because smoothing suppresses exactly the derivatives the residual measures. PDEBench's cRMSE and bRMSE exist because bulk error and residuals both miss conservation and boundary violations, and Lab 02's own boundary error of 0.651 K against a true boundary mean of 0.377 K is what that looks like in practice.
LAB Not yet run. Lab 02's operator is \(q\mapsto T\) with \(-kt\nabla^2T=q\) and \(T=0\) on the edge, and its residual is cheap to compute: a five-point Laplacian of the model output, times \(kt\), minus \(q\), normalized by \(\|q\|\). The boundary term is \(\|T\|\) on the edge cells. Both terms use the same normalization as the data term, so \(\lambda=1\) means genuinely equal weight and not an arbitrary scale.
../labs/02_thermal/run_thermal.py and add \(\lambda_r\|\text{residual}\|_2/\|q\|_2 + \lambda_b\|T_{\text{edge}}\|_2\) to the existing relative \(L^2\) data loss.Two predictions to record before running, because a prediction written afterwards is not a prediction. From PINO's Darcy result, expect a gain at the small data fraction and none at full data. From the conditioning derivation above, expect instability at \(\lambda_r = 10\) unless the learning rate drops, and expect it to be worse if you repeat the sweep at 128² than at 64², by roughly the factor of 16 the \(N^4\) scaling predicts.
research_notes/web-papers.md:585-653, :613-626.research_notes/lineage.md:104-113; library support at research_notes/repos.md:63, :99.web-papers.md:452-519. Physics-informed DeepONet weight sweep: :565-573.research_notes/foundation.md:154-167.research_notes/pdf-mgn-pdebench.md:286. JMLR physics-informed direction: pdf-fno-jmlr.md:719-724.research_notes/company-episode.md:107-112.You can now name the two distinct ways physics substitutes for data, quote the best published numbers for each, and derive the reason the second one gets harder, not easier as the problem gets more interesting.
neuraloperator, and say what their existence tells you about the state of the theory.