Module 3 of 14 · one week
A kernel integral layer costs \(O(J^2)\) on \(J\) grid points, which is unaffordable at any usable resolution. GNO, LNO, MGNO and FNO are four ways to pay less, and each one is a numerical linear algebra idea that predates deep learning by decades: Nyström sampling, low rank, fast multipole, and the FFT.
Module 2 derived the Fourier layer by making one assumption, that the kernel depends only on \(x-y\) on a periodic box, and it flagged at the time that this was a trade, not a free win. This module is where the trade gets priced. The JMLR paper offers four parameterizations of the same general layer, and FNO is the one with no learned kernel network at all, which is exactly why it is the fastest and why it is the worst of the four at small parameter budgets.
There is also a practical reason not to learn only FNO. The geometry architectures of Module 6, and GINO and Transolver in particular, are built out of the graph and low-rank ideas, not the Fourier one, because their meshes are not uniform grids and the FFT has no uniform grid to run on. If the only layer you know is the spectral one, those papers read as unrelated inventions instead of as the obvious next move.
Module 2 wrote down the Green's function of the destination example and checked it by its slope jump: for \(-u''=f\) on \((0,1)\) with zero boundary values,
\[u(x) = \int_0^1 g(x,y)f(y)\,dy, \qquad g(x,y) = \min(x,y)\,\big[1 - \max(x,y)\big].\]
What that module did not do is look at the operator's spectrum, and the spectrum is what decides whether any of the four cost-saving tricks below can work at all. Write \(G\) for the map \(f \mapsto u\). Module 2 used \(G(x,y)\) for the kernel; this module writes the kernel \(g(x,y)\) so that the operator can have the capital letter.
Claim 1 (the eigenvalues). \(G\big[\sin(n\pi x)\big] = (n\pi)^{-2}\sin(n\pi x)\) for every \(n \ge 1\).
Proof. \(-\partial_{xx}\sin(n\pi x) = (n\pi)^2\sin(n\pi x)\), and \(\sin(n\pi x)\) vanishes at both endpoints, so it is in the domain where \(G\) inverts \(-\partial_{xx}\). Applying \(G\) to both sides and dividing by \((n\pi)^2\) gives the claim. \(\square\)
The functions \(\sqrt2\sin(n\pi x)\) are an orthonormal basis of \(L^2(0,1)\), so \(G\) is diagonal in that basis with eigenvalues \(\lambda_n = (n\pi)^{-2}\): \(0.1013\) at \(n=1\), \(1.013\times10^{-3}\) at \(n=10\), \(1.013\times10^{-5}\) at \(n=100\).
Claim 2 (truncation has a uniform error, unlike the identity). Let \(G_m\) keep only the first \(m\) modes. Then \[\|G - G_m\|_{L^2\to L^2} = \frac{1}{\pi^2(m+1)^2}.\]
Proof. \(G - G_m\) is diagonal with entries \(\lambda_n\) for \(n > m\) and zero otherwise. The norm of a diagonal operator is the supremum of the absolute values of its entries, which here is the largest surviving eigenvalue, \(\lambda_{m+1}\). \(\square\)
Consequence, and the contrast with Module 0. Fact 3 showed that truncating the identity has \(\|I-P_m\|=1\) for every \(m\), so no mode count gives a uniform guarantee on the unit ball. Truncating \(G\) has \(\|G-G_m\| = \pi^{-2}(m+1)^{-2}\), which does. Keeping ten modes bounds the error by \(8.4\times10^{-4}\) for every input of norm one, with no compactness assumption and no restriction on the input set at all. The difference is that \(G\) smooths and the identity does not. This also proves \(G\) is compact, since it is the norm limit of the finite-rank operators \(G_m\).
That contrast is the load-bearing fact of this module, and it is easiest to hold as a design rule. Low-rank and truncation methods work when the operator's spectrum decays, and elliptic solution operators have decaying spectra because they integrate. They degrade on operators that do not smooth, and the degradation shows up in robustness before it shows up in clean accuracy. The JMLR paper's Table 9 (p.47) makes the distinction sharply. On clean data FNO handles one-dimensional advection with a random square-wave input as well as anything: training error 0.002, test error 0.002. Add 10 percent input noise and that test error becomes 0.094, a factor of 47, while Darcy moves only from 0.011 to 0.012 and Navier-Stokes from 0.024 to 0.039. The paper's own reading is that FNO "is less robust on the advection equation, which is not smoothing, and on Burgers equation which, whilst smoothing also forms steep fronts". So a headline accuracy number on clean data does not tell you whether the operator's spectrum decays. Module 1 gave the theoretical half of the same point.
The paper also checks all of this on the operator itself: a neural operator with zero hidden layers learns \(G\) to relative test error \(10^{-7}\), and its Figure 7 puts the learned kernel next to the analytic \(g\) (research_notes/pdf-fno-jmlr.md:537-545). Every remaining operator in the course has a kernel nobody can write down, so this is the last time the answer is checkable in closed form. Use it as the calibration point.
Module 2 introduced the layer in its simplest form. The paper's Equation 6 on page 10 is the same object with the pieces named, and the naming matters because the four answers below each modify a different piece:
\[v_0 = P(a), \qquad v_{t+1}(x) = \sigma\Big(W_t\,v_t(x) + \int_D \kappa_t(x,y)\,v_t(y)\,dy + b_t(x)\Big), \qquad G_\theta(a) = Q(v_T).\]
\(P\) and \(Q\) are pointwise lifting and projection maps. \(W_t\) is a matrix mixing channels at a single point, so it is local. The integral is the only term that moves information between locations, so it is the only term any of the four answers touch. \(b_t\) is a bias function of \(x\), not a vector, which is what lets a layer know where it is on a non-periodic domain. And \(\sigma\) acts pointwise.
The paper allows the kernel to depend on more than the two positions, and the three forms it lists (Equations 7 to 9) need keeping straight because they have different consequences.
Attention is a special case, and the proof says which one. JMLR Proposition 6, on page 28, states that "the attention mechanism in transformer models is a special case of a neural operator layer". The proof picks the third kernel form, sets the bias to zero and \(W\) to the identity, chooses a softmax kernel, and discretizes the integral by Monte Carlo, at which point the layer becomes a pre-normalized single-headed attention block as in Vaswani et al. This is a real containment rather than a slogan, and it also shows what a transformer is missing: the Monte Carlo discretization is over the tokens you happen to have, so nothing in the construction makes the layer converge as the token set refines. That is the gap Module 6's Transolver has to close.
Claim 3. If \(\kappa \in L^2(D\times D)\) and \((Kv)(x) = \int_D\kappa(x,y)v(y)\,dy\), then \(\|Kv\|_2 \le \|\kappa\|_{L^2(D\times D)}\|v\|_2\).
Proof. Fix \(x\) and apply Cauchy-Schwarz in \(y\): \(|Kv(x)|^2 \le \int_D|\kappa(x,y)|^2dy \cdot \int_D|v(y)|^2dy\). Integrate over \(x\); the first factor integrates to \(\|\kappa\|^2_{L^2(D\times D)}\) and the second is a constant. \(\square\)
Consequence. With a Lipschitz activation, one layer obeys \(\|L(v)-L(w)\|_2 \le \mathrm{Lip}(\sigma)\big(\|W\| + \|\kappa\|_{L^2}\big)\|v-w\|_2\), and \(T\) layers multiply those constants. The bound is loose, but it names the amplification correctly: depth times local gain plus kernel gain. It also makes two limits visible. Because \(W\) may contain an identity component, the layer as a whole need not be compact even though \(K\) is; and the bound says nothing about accuracy, only about stability. Module 5 reuses this decomposition when it states the discretization theorem.
Discretize the integral on a mesh of \(J\) points and it becomes \(\sum_j \kappa(x_i,y_j)v(y_j)\Delta y\) for every one of the \(J\) values of \(i\): a dense \(J\times J\) matrix times a vector, per layer, per channel pair (JMLR p.13, Equation 14).
The size of that cost matters before reading the four repairs. At \(J = 10^6\), which is a \(1024\times1024\) grid and modest by simulation standards, the dense kernel matrix has \(10^{12}\) entries. Stored in single precision that is 4 terabytes, for one channel pair in one layer, before any arithmetic happens. For comparison, Lab 03 runs a complete four-layer FNO forward pass at \(1024^2 = 1{,}048{,}576\) points in 143 ms using 1.02 GiB of peak memory on a laptop GPU (../labs/03_scale/results_scale.json). The gap between those two numbers is what this module is about.
Subsample \(J' \ll J\) points and approximate the kernel matrix by the rank-\(J'\) factorization \(K \approx K_{JJ'}K_{J'J'}^{-1}K_{J'J}\). Or truncate the integral to a ball of radius \(r\) around each \(x\), costing \(O(r^dJ^2)\). Truncation loses long-range interaction, but composing \(L\) truncated layers with \(2^{L-1}r \ge 1\) recovers the full reach at cost \(O(Lr^dJ^2)\). Implemented as message passing on the mesh graph, where averaging over neighbours is exactly the Monte Carlo estimate of the truncated integral. Pages 14 to 16.
Impose \(\kappa(x,y)=\sum_{j=1}^r\varphi_j(x)\psi_j(y)\). Then \(u(x)=\sum_j\langle\psi_j,v\rangle\varphi_j(x)\), and the \(r\) inner products do not depend on \(x\), so they are computed once: cost \(O(rJ)\). This is a truncated SVD of a rank-\(r\) operator, and Claim 2 above is exactly the statement that it works when the spectrum decays. It is also DeepONet's structure made discretization-consistent, which Module 4 takes up. Pages 16 to 17.
Split \(K=K_1+\dots+K_L\) by interaction range. The short-range part is sparse and full rank; the long-range part is dense and low rank. Build a hierarchy of coarser meshes \(J_1\ge\dots\ge J_L\) with growing radii and run a V-cycle down and back up, as multigrid does. Choosing \(J_l\) and \(r_l\) so that \(\sum_l J_l^2r_l^d = O(J)\) makes the whole thing linear. Pages 17 to 20.
Impose \(\kappa(x,y)=\kappa(x-y)\) on a periodic box. Then \(Kv=\mathcal F^{-1}(\hat\kappa\cdot\hat v)\), and the multipliers \(\hat\kappa(k)\) for \(|k|\le k_{\max}\) are the parameters themselves. Cost \(O(J\log J)\), dominated by the FFT, with no kernel network to evaluate. The price is a uniform grid and a periodic domain, both of which Module 6 has to undo. Pages 21 to 24; derived in full in Module 2.
Why the CNN is not a fifth answer. A convolution with a kernel of fixed radius \(r\) in physical units is a legitimate GNO layer, and nothing here rules it out. A standard CNN fixes the kernel in grid points, so as \(J\to\infty\) its physical support shrinks to a point and the layer becomes local. The only repair is to grow the filter with \(J\), and then the parameter count diverges. That is JMLR Section 4.1, pages 15 to 16 (research_notes/pdf-fno-jmlr.md:372-379), and Module 1's Lab 01 run C measures it: the same weights that match an FNO exactly at 64 points give a relative error of 1.46 at 1024, worse than predicting zero.
All from JMLR Section 7.4, on a single V100, with \(N=1000\), 4 layers and 500 epochs (research_notes/pdf-fno-jmlr.md:527-535, :627-661). The last two columns are accuracy on the two benchmark problems; the timing column is what one epoch of Darcy costs.
| Method | Idea (Table 6) | Complexity (Table 7) | s/epoch | Darcy s=141 (Table 2) | Burgers s=1024 (Table 3) |
|---|---|---|---|---|---|
| GNO | Nyström, kernel network, graph-based | \(O(J'^2r^2)\) | 4 | 0.0332 | 0.0651 |
| LNO | low rank, kernel network | \(O(J)\) | 20 | 0.0461 | 0.0217 |
| MGNO | multi-level graphs on GNO | \(\sum_l O(J_l^2r_l^2) \sim O(J)\) | 8 | 0.0428 | 0.0374 |
| FNO | convolution theorem, no kernel network | \(O(J\log J)\) | 4 | 0.0109 | 0.0018 |
Three things in that table matter more than the ranking. LNO is \(O(J)\) on paper and the slowest of the four in wall clock, five times slower than FNO, because it carries a kernel network that has to be evaluated; complexity classes do not settle wall clock when the constants differ by an order of magnitude. The FNO Burgers entry is the JMLR paper's 0.0018 rather than the ICLR paper's 0.0160 at the same resolution, for the reason Module 2 documents, and the other three entries are identical in both papers. And the honest sentence in the paper's own conclusion is that "the full \(O(J^2)\) integration method still outperforms the FNO by about 40%, albeit at greater cost" (p.58, research_notes/pdf-fno-jmlr.md:719-721). All four answers buy speed with accuracy. FNO buys the most of both.
The most useful table in the paper is also the one most often skipped, because it is the only one where FNO loses.
| Method | 10³ params | 10⁴ | 10⁵ | 10⁶ |
|---|---|---|---|---|
| GNO | 0.075 | 0.065 | 0.060 | 0.035 |
| LNO | 0.080 | 0.070 | 0.060 | 0.040 |
| MGNO | 0.070 | 0.050 | 0.040 | 0.030 |
| FNO | 0.200 | 0.035 | 0.020 | 0.015 |
At a budget of one thousand parameters FNO is the worst of the four by a factor of nearly three, 0.200 against 0.070 to 0.080. From ten thousand upward it is the best, and by \(10^6\) it is twice as accurate as the runner-up. The paper explains the shape rather than hiding it: the other three "have the kernel networks, the slope of their error rates are flat: they can work with a very small number of parameters. On the other hand, FNO does not have the sub-network. It needs at a larger magnitude of parameters to obtain an acceptable error rate" (p.46).
Do not read fine structure into these sixteen numbers. The caption states that "the errors above are approximated value roundup to 0.05", and that each entry is the lowest test error the model reached subject to a parameter cap. So the difference between GNO's 0.075 and MGNO's 0.070 at \(10^3\) carries no information, and neither does GNO 0.060 against LNO 0.060 at \(10^5\). Exactly one comparison in the table is larger than the stated rounding: FNO's 0.200 at \(10^3\) against everything else in that column. That single number is the finding. Everything else is a tie within the paper's own error bars, and any course that ranks the four methods off this table is over-reading it.
LAB Part 1, not yet run. Reproduce two rows of the refinability table with neuraloperator on the bundled Darcy 16 data: an FNO with hidden_channels and n_modes tuned to land near \(10^3\), \(10^4\) and \(10^5\) parameters, measured with count_model_params, against GNO-style layers through neuralop.layers.gno_block.GNOBlock. Plot error against parameter count on log axes. The bundled files darcy_train_16.pt, darcy_test_16.pt and darcy_test_32.pt live in the library's GitHub data folder rather than the pip package, so download them into ../data/ the way Lab 01 did. The prediction to test is the crossover between \(10^3\) and \(10^4\) parameters, and the thing to watch is whether it survives at all once the rounding in the warning above is taken seriously.
Part 2, reading. Swapping the transform is a design family rather than a trick, and each member picks a basis suited to the physics. Laplace Neural Operator replaces Fourier with Laplace to get poles and residues, which matters for transients that a periodic basis represents badly (arXiv 2303.10528). Wavelet Neural Operator uses wavelets to buy spatial localization back (arXiv 2205.02191). Spherical FNO uses spherical harmonics because the atmosphere is a sphere and a periodic box is the wrong domain for it (arXiv 2306.03838). Read the three abstracts through research_notes/lineage.md:148-168 and :82-91, and write one paragraph placing each of the three against the four answers above. All three are variations on the fourth answer, and saying precisely why is the exercise.
research_notes/pdf-fno-jmlr.md:352-426.21-1524.pdf p.28; notes at research_notes/pdf-fno-jmlr.md:453-463.21-1524.pdf; notes at research_notes/pdf-fno-jmlr.md:627-661.research_notes/pdf-fno-jmlr.md:719-721.research_notes/pdf-fno-jmlr.md:537-545. Advection robustness, Table 9, p.46.../labs/03_scale/results_scale.json. CNN failure at unseen grids: ../labs/01_burgers/results_c_rf49.json.research_notes/lineage.md:38-58, :82-91, :148-168.You can now say why truncation works on an elliptic solution operator and not on the identity, place any new operator architecture into one of four cost strategies, and read the refinability table without over-reading it.