AMStraMGRAM : Adaptive Multi-cutoff Strategy Modification for ANaGRAM

## Abstract

Abstract Recent works have shown that natural gradient methods can significantly outperform standard optimizers when training physics-informed neural networks (PINNs). In this paper, we analyze the training dynamics of PINNs optimized with ANaGRAM, a natural-gradient-inspired approach employing singular value decomposition with cutoff regularization. Building on this analysis, we propose a multi-cutoff adaptation strategy that further enhances ANaGRAM’s performance. Experiments on benchmark PDEs validate the effectiveness of our method, which allows to reach machine precision on some experiments. To provide theoretical grounding, we develop a framework based on spectral theory that explains the necessity of regularization and extend previous shown connections with Green’s functions theory.

## 1 Introduction

Physics-informed neural networks (PINNs) have recently emerged as a promising alternative for the numerical solution of partial differential equations (PDEs) (Raissi et al., 2019). By leveraging neural networks as universal function approximators (Leshno et al., 1993), PINNs replace traditional mesh-based discretizations with sampling-based collocation methods, enabling a straightforward extension to high-dimensional domains. This mesh-free formulation not only circumvents the “curse of dimensionality” inherent in grid-based approaches, but also allows continuous evaluation of the solution throughout the domain without explicit mesh generation (Cuomo et al., 2022).

Despite these advantages, achieving low training error with PINNs remains a major challenge (Wang et al., 2023; Urbán et al., 2025; Kiyani et al., 2025; De Ryck et al., 2024).
Open questions include how to select and distribute collocation points, how to balance the PDE residual against boundary-condition penalties, and which optimization strategies most effectively minimize the composite loss (Krishnapriyan et al., 2021; Wang et al., 2021; McClenny & Braga-Neto, 2022).

A different line of research has recently reexamined PINNs from the perspective of functional geometry (Müller & Zeinhofer, 2023, 2024; Jnini et al., 2024), providing a mathematically principled view of the training dynamics.
In this vein, the ANaGRAM algorithm (Schwencke & Furtlehner, 2025) applies a natural-gradient update (Amari, 1998; Ollivier, 2015), based on a reinterpretation and generalization of the neural tangent kernel (NTK; Jacot et al. (2018)) as the kernel of the projection onto the neural network’s tangent space.
This leads to a notion of the empirical natural gradient that projects the true functional gradient onto the empirical tangent space, yielding significantly faster convergence and lower errors compared to standard optimizers on PDE benchmarks.

Nevertheless, while ANaGRAM improves over standard optimizers, it still falls short of the accuracy attained by classical mesh-based methods, such as the finite element method (Grossmann et al., 2024).
Moreover, its final performance is highly affected by the way the pseudo-inverse of the feature matrix is computed. In particular, ANaGRAM sets a fixed level of *cutoff*: a value below which the singular values of the feature matrix are ignored, i.e. it controls how much loss signal is incorporated into an update. ANaGRAM’s cutoff is currently chosen manually, as no automatic selection procedure has been proposed.

In this paper, we study the performance and training dynamics of ANaGRAM, with a particular focus on the role of the chosen cutoff. Typically, the training loss of ANaGRAM exhibits the slow convergence at the early iterations followed by a sudden drop at the end of the training – similar behavior is shown by the eNGD method (Müller & Zeinhofer, 2023). We discover that it is closely connected to what we further refer as the *flattening phenomenon*, which we define and characterize using the *reconstruction error*: a novel metric that measures how much of the loss signal is lost by different choices of cutoffs. Relying on the adaptive multi-cutoff strategy, our new algorithm AMStraMGRAM manages to capitalize on this phenomenon, resulting in a significant improvement (of several orders of magnitude) on various PDE benchmarks. To complement our empirical findings, we also present a functional-analytic view linking cutoff (and ridge regularization) to (generalized) Green operator theory, clarifying why cutoff regularization is essential and not just a mere fix to stabilize training.

## 2 Problem Statement

### 2.1 Differential Operators and Physics-Informed Neural Networks (PINNs)

Let $\Omega\subset\mathbb{R}^{d}$ be a domain. We introduce two operators, $D$ and $B$, defined on a Hilbert space $\mathcal{H}$ of real-valued functions, acting respectively on $\Omega$ and on its boundary $\partial\Omega$:

$$ $\displaystyle D:$ $\displaystyle\left\{\begin{array}[]{lll}\mathcal{H}&\to&\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)\\ u&\mapsto&D[u]\end{array}\right.,$ $\displaystyle B:$ $\displaystyle\left\{\begin{array}[]{lll}\mathcal{H}&\to&\textrm{L}^{2}(\partial\Omega\to\mathbb{R},\sigma)\\ u&\mapsto&B[u]\end{array}\right..$ (5) $$

Here, $D$ denotes a differential operator, while $B$ represents a boundary operator.
A function $u\in\mathcal{H}$ is said to be a *classical solution* to the *Partial Differential Equation* (PDE) associated with $D$ and $B$ if it satisfies

$$ $\begin{cases}D(u)=f\in\textrm{L}^{2}(\Omega\to\mathbb{R},\mu),&\text{in }\Omega,\\ B(u)=g\in\textrm{L}^{2}(\partial\Omega\to\mathbb{R},\sigma),&\text{on }\partial\Omega,\end{cases}$ (6) $$

A *physics-informed neural network* (PINN) approximates the solution $u$ by a parametric model $u_{{\bm{\theta}}}$, where $u_{{\bm{\theta}}}$ is a neural network with parameters ${\bm{\theta}}\in\mathbb{R}^{P}$. The learning objective is to minimize the empirical loss

$$ $\ell_{D,B}({\bm{\theta}}):=\frac{1}{2S_{D}}\sum_{i=1}^{S_{D}}\left(D[u_{{\bm{\theta}}}](x_{i}^{D})-f(x_{i}^{D})\right)^{2}+\frac{1}{2S_{B}}\sum_{i=1}^{S_{B}}\left(B[u_{{\bm{\theta}}}](x_{i}^{B})-g(x_{i}^{B})\right)^{2}.$ (7) $$

### 2.2 PINNs Optimizers

Training PINNs is notoriously challenging. Issues such as spectral bias, where networks struggle to learn high-frequency components, and the difficulty of balancing residual and boundary loss terms—often with vastly different magnitudes— result in unsatisfactory performance of standard deep learning optimizers (Wang et al., 2021; De Ryck et al., 2024; Krishnapriyan et al., 2021; Liu et al., 2024).

To mitigate these challenges, researchers have proposed various strategies. These include adaptive sampling approaches that focus on regions with high error (Krishnapriyan et al., 2021), dynamic loss weighting schemes (McClenny & Braga-Neto, 2022), and architectural modifications (Wang et al., 2024). Another promising line of research has focused on modifying the optimizers. In particular, two main branches of optimization approaches for PINNs have emerged:

- (i)
Second-Order Methods. These methods, based on Quasi-Newton techniques, particularly the BFGS algorithm (Nocedal & Wright, 1999, Chapter 6) and its memory-efficient approximation L-BFGS (Liu & Nocedal, 1989), address some of the training difficulties by considering the curvature of the loss landscape. This curvature arises from the non-linearities of both the neural network and the differential operators (Rathore et al., 2024). Recently, Urbán et al. (2025) extended this approach by modifying the self-scaled BFGS (SSBFGS; Al-Baali, 1998) and self-scaled Broyden (SSBroyden; Al-Baali & Khalfan, 2005), along with other computational enhancements such as point resampling (Wu et al., 2023) and boundary condition enforcement (Wang et al., 2023), achieving state-of-the-art results (Kiyani et al., 2025).
- (ii)
Natural Gradient Methods. In contrast to second-order methods, natural gradient methods are first-order techniques(^1^11contrary to a widespread misconception, which arises from their analogy in the context of information theory) that provide a principled way to incorporate the geometry and metric structure of the problem space. Initially introduced in the context of information geometry by Amari (1998) and later extended by Ollivier (2015), these methods were introduced for PINNs by Müller & Zeinhofer (2023). In subsequent work, Schwencke & Furtlehner (2025) connected these methods to kernel methods, yielding an efficient implementation they linked to Green’s function theory (Duffy, 2015).

### 2.3 Natural Gradient Methods for PINNs

As a preliminary observation highlighted in Schwencke & Furtlehner (2025, Section 4.1), PINNs can be interpreted as a quadratic regression problem. This viewpoint arises naturally once the parametric model $u_{{\bm{\theta}}}$ is replaced with the following compound model:

$$ $(D,B)\circ u:\left\{\begin{array}[]{lllll}\mathbb{R}^{P}&\to&\mathcal{H}&\to&\textrm{L}^{2}(\Omega,\mu)\times\textrm{L}^{2}(\partial\Omega,\sigma)\\ {\bm{\theta}}&\mapsto&u_{{\bm{\theta}}}&\mapsto&(D[u_{{\bm{\theta}}}],B[u_{{\bm{\theta}}}])\end{array}\right..$ (8) $$

For ease of exposition, and without loss of generality, we restrict attention to regression in $\textrm{L}^{2}(\Omega,\mu)$.
Given $f\in\textrm{L}^{2}(\Omega,\mu)$, we define the associated empirical loss

$$ $\ell({\bm{\theta}}):=\frac{1}{2S}\sum_{i=1}^{S}\left(u_{{\bm{\theta}}}(x_{i})-f(x_{i})\right)^{2},$ (9) $$

which can be seen as a discretization of the functional loss

$$ $\mathcal{L}(u):=\tfrac{1}{2}\left\|u-f\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)},\qquad u\in\textrm{L}^{2}(\Omega,\mu).$ (10) $$

The natural gradient approach seeks to compute the optimal update direction in function space and then pull it back to parameter space.
A single Fréchet derivative of the functional loss [Equation˜10](#S2.E10) yields $\nabla\mathcal{L}_{|u}=u-f$.
The key insight is that admissible updates are constrained to the tangent space of the parametric model,

$$ $T_{{\bm{\theta}}}\mathcal{M}:=\mathrm{Im}(\mathrm{d}u_{{\bm{\theta}}})=\operatorname*{Span}\left(\partial_{p}u_{{\bm{\theta}}}:1\leq p\leq P\right)\subset\mathcal{H},$ (11) $$

where $\mathcal{M}:=\mathrm{Im}(u)=\{u_{{\bm{\theta}}}:{\bm{\theta}}\in\mathbb{R}^{P}\}\subset\mathcal{H}$ is the manifold of functions parametrized by ${\bm{\theta}}$.
Thus, the optimal update in function space is the projection of $\nabla\mathcal{L}_{|u}$ onto the tangent space (cf. [Figure˜4](#A1.F4)),

$$ $\displaystyle u_{{\bm{\theta}}_{t+1}}$ $\displaystyle\leftarrow u_{{\bm{\theta}}_{t}}-\eta\,\Pi_{T_{{\bm{\theta}}_{t}}\mathcal{M}}\!\left(\nabla\mathcal{L}_{u_{{\bm{\theta}}_{t}}}\right)\,;$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,\mathrm{d}u_{{\bm{\theta}}_{t}}^{\dagger}\!\left(\Pi_{T_{{\bm{\theta}}_{t}}\mathcal{M}}\!\left(\nabla\mathcal{L}_{u_{{\bm{\theta}}_{t}}}\right)\right),$ (12) $$

where the second equation is simply the pullback of the functional update to parameter space. We prove in [Section˜H.1](#A8.SS1) that this update is equivalent to the Gram–matrix formulation:

$$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,G_{{\bm{\theta}}_{t}}^{\dagger}\nabla\ell({\bm{\theta}}_{t})\,;$ $\displaystyle{G_{{\bm{\theta}}_{t}}}_{p,q}$ $\displaystyle:=\left\langle\partial_{p}u_{{\bm{\theta}}_{t}}\,,\,\partial_{q}u_{{\bm{\theta}}_{t}}\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}.$ (13) $$

### 2.4 ANaGRAM: Empirical Natural Gradient

The $O(P^{3})$ complexity of matrix inversion in [Equation˜13](#S2.E13) renders a direct implementation prohibitively expensive.
ANaGRAM (Schwencke & Furtlehner, 2025) circumvents this by exploiting a motivated approximation.
The key observation is that the update can be expressed in terms of the empirical feature matrix $\widehat{\phi}\in\mathbb{R}^{S\times P}$ and the empirical functional residuals $\widehat{\mathcal{L}_{{\bm{\theta}}}}\in\mathbb{R}^{S}$:

$$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,\widehat{\phi}^{\dagger}\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}_{t}};$ $\displaystyle\widehat{\phi}_{i,p}$ $\displaystyle:=\partial_{p}u_{{\bm{\theta}}}(x_{i});$ $\displaystyle\left(\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}}\right)_{i}$ $\displaystyle:=u_{{\bm{\theta}}}(x_{i})-f(x_{i}).$ (14) $$

Here, the pseudo-inverse is computed via singular value decomposition (SVD):
$\widehat{\phi}^{\dagger}=\widehat{U}\widehat{\Delta}^{\dagger}\widehat{V}^{T}$ with $\widehat{\phi}=\widehat{V}\widehat{\Delta}\widehat{U}^{T}$,
where $\widehat{U}\in\mathbb{R}^{P\times\textrm{r}_{\text{svd}}}$, $\widehat{\Delta}\in\mathbb{R}^{\textrm{r}_{\text{svd}}\times\textrm{r}_{\text{svd}}}$, $\widehat{V}\in\mathbb{R}^{S\times\textrm{r}_{\text{svd}}}$, and $\textrm{r}_{\text{svd}}=\min(P,S)$.
This reduces computational cost to $O(\min(PS^{2},P^{2}S))$, which is tractable in practice.
A comparable complexity was later obtained by Guzmán-Cordero et al. (2025) using a Cholesky factorization approach.

For further details on the derivation of the empirical natural gradient, we refer to Schwencke & Furtlehner (2025).
In what follows, we adopt a slight abuse of notation by omitting the explicit dependence on ${\bm{\theta}}$ whenever it is clear from context.
When iteration indices matter, we explicitly write $t$ to emphasize the connection to ${\bm{\theta}}_{t}$.

### 2.5 Regularization

As discussed in [Section˜G.1](#A7.SS1), the type of problem we consider is ill-conditioned, which necessitates the use of regularization.
We distinguish between two main regularization schemes: (i) *ridge regression*, which consists in adding a factor $\alpha^{2}I_{d}$ (or, according to conventions, $\alpha^{-2}I_{d}$) to the Gram matrix $G_{\bm{\theta}}$ in [Equation˜13](#S2.E13) (or its approximation $\widehat{\mathcal{G}}_{{\bm{\theta}}}$), thereby making it invertible or (ii) *cutoff regularization*, a scheme that applies a binary threshold (used in ANaGRAM):

$$ $\widehat{\Delta}_{t,i}^{\dagger}=\begin{cases}\widehat{\Delta}_{t,i}^{-1},&\text{if }\widehat{\Delta}_{t,i}\geq\alpha,\\ 0,&\text{otherwise}.\end{cases}$ (15) $$

Here $\alpha$ denotes the cutoff threshold. This regularization is the focus of our analysis in [Section˜3](#S3). For completeness, we provide a geometric interpretation of each scheme in [Appendix˜G](#A7).
We further show that cutoff regularization extends previously established connections between natural gradient methods and Green’s function theory (Schwencke & Furtlehner, 2025). In particular, we obtain:

###### Theorem 1 .

The generalized Green’s function of the operator $D$ in the regularized space $\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$ is given, for all $x,y\in\Omega$, by

$$ $g_{D}(x,y):=D[k_{D}(x,\cdot)](y),$ (16) $$

where $\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$ is a regularized space with reproducing kernel $k_{D}$, defined in [Section˜G.4](#A7.SS4).

###### Theorem 1 .

## 3 Insights on ANaGRAM’s Training Dynamics

In this section, we will look at relevant quantities of interest to understand this empirical phenomenon.

### 3.1 Reconstruction Error of Functional Gradient

Let ${\bm{\theta}}\in\mathbb{R}^{P}$, the empirical feature matrix $\widehat{\phi}\in\mathbb{R}^{S\times P}$, and the empirical functional gradient $\widehat{\nabla\mathcal{L}}\in\mathbb{R}^{S}$ as defined in [Equation˜14](#S2.E14).
Let us consider various empirical tangent spaces formed by taking different ranges of right singular vectors of $\widehat{\phi}=\widehat{U}\widehat{\Delta}\widehat{V}^{T}$, i.e. $\widehat{T^{M}_{N}\mathcal{M}}=\operatorname*{Span}(\widehat{V}_{t,{i}}\,:\,M\leq i\leq N)$.
For $1\leq N\leq\textrm{r}_{\text{svd}}$, reconstruction error measures how much information from the functional gradient signal is lost when considering only first $N$ components in SVD (the error caused by the projection onto the empirical tangent space $\widehat{T^{0}_{N}\mathcal{M}}$) is defined as follows

$$ $\textrm{RCE}^{S}_{N}:=\frac{1}{\sqrt{S}}\left\|\widehat{V}\Pi^{0}_{N}\widehat{V}^{T}\widehat{\nabla\mathcal{L}}-\widehat{\nabla\mathcal{L}}\right\|_{\mathbb{R}^{S}}=\frac{1}{\sqrt{S}}\|\Pi^{\bot}_{\widehat{T^{0}_{N}\mathcal{M}}}\widehat{\nabla\mathcal{L}}-\widehat{\nabla\mathcal{L}}\|,$ (17) $$

where we define $\Pi^{M}_{N}\in\mathbb{R}^{\textrm{r}_{\text{svd}}\times\textrm{r}_{\text{svd}}}$ as a projection operator onto $\widehat{T^{M}_{N}\mathcal{M}}$:

$$ $\Pi^{M}_{N}=\sum_{p=M+1}^{N}{\bm{e}}^{(p)}{\bm{e}}^{(p)^{T}},$ (18) $$

with $(\mathbf{e}^{(p)})_{1\leq p\leq\textrm{r}_{\text{svd}}}$ being the canonical basis of $\mathbb{R}^{\textrm{r}_{\text{svd}}}$.

###### Proposition 1 .

$\textrm{RCE}^{S}_{N}$ is a non-increasing function of $N$, i.e. for all $1\leq M,N\leq\textrm{r}_{\text{svd}}$:

$$ $M\leq N\implies\textrm{RCE}^{S}_{M}\geq\textrm{RCE}^{S}_{N}.$ (19) $$

Furthermore, assuming that $(x_{i})_{i=1}^{S}$ are i.i.d sampled from $\mu$, we have $\mu$-almost surely

$$ $\displaystyle\lim\limits_{S\to\infty}\textrm{RCE}^{S}_{N}=\left\|\nabla\mathcal{L}_{u_{{\bm{\theta}}}}-\Pi^{\bot}_{T^{0}_{N}\mathcal{M}}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|_{\textrm{L}^{2}(\Omega,\mu)}=\left\|\Pi^{\bot}_{\left[{T^{0}_{N}\mathcal{M}}\right]^{\bot}}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|_{\textrm{L}^{2}(\Omega,\mu)},$ (20) $$

where $T^{M}_{N}\mathcal{M}=\operatorname*{Span}(V_{t,{i}}\,:\,M\leq i\leq N)$,
while $\left(V_{t,{i}}\right)_{1\leq i\leq\textrm{r}_{\text{svd}}}$ are the right singular-vectors of the differential $\textrm{d}u_{{\bm{\theta}}}$ ordered in a decreasing order according to their associated singular values.

###### Remark 1 .

Note that
$\widehat{V}_{t,{i}}\in\mathbb{R}^{S}$ for $i\in 1,\dots,N$, the right singular vectors of $\widehat{\phi}$, can be seen as discretized versions of $V_{t,{i}}$ from [˜1](#ThmProp1). Indeed, a weak convergence holds, i.e. $\forall h\in\mathcal{H}$, $\frac{1}{S}\sum_{j=1}^{S}\widehat{V}_{t,{i,j}}h_{j}=\frac{1}{S}\sum_{j=1}^{S}V_{t,{i}}(x_{j})h(x_{j})\stackrel{{\scriptstyle S\to\infty}}{{\to}}\left\langle V_{t,{i}},h\right\rangle_{L^{2}}$.

Proof of [˜1](#ThmProp1) can be found in [Section˜H.3](#A8.SS3).
From [˜1](#ThmProp1) RCE is related to the concept of *expressivity bottleneck* illustrated in Verbockhaven et al. (2024), and measures what part of the learning signal is not captured by truncating at $N$ components for natural gradient computation.
Therefore, this metric allows us to explicitly estimate and compare different cutoff choices. Note that this metric incurs no additional computational cost since ANaGRAM already computes the required SVD.

###### Proposition 1 .

###### Remark 1 .

### 3.2 Empirical Observations: Flattening

Here we illustrate the evolution of training loss and reconstruction error, where [Figure˜1](#S3.F1) schematically outlines key stages of ANaGRAM’s training dynamics. The plot of a real experiment is provided in [Appendix˜E](#A5).

Figure: (a) Early iterations, RCE at intersection with singular values is above the desired precision threshold.
Refer to caption: /html/2510.15998/assets/figures/ANAGRAM_legend.png

Let $\alpha$ is a cutoff level (also referred to as precision) and $\textrm{r}_{\text{cutoff}}$ denote the number of components retained by the cutoff, i.e., $\textrm{r}_{\text{cutoff}}(t)=\max\{j:\widehat{\Delta}_{t,j}\geq\alpha\}$.
In [Figure˜1](#S3.F1), we observe different stages of the training.
First, the reconstruction error is above the wanted precision ([Figure˜5(a)](#A5.F5.sf1)). As the training progresses, the training loss drops and the reconstruction error drops until reaching the cutoff precision ([Figure˜5(b)](#A5.F5.sf2)). Eventually, the reconstruction error drops below the cutoff threshold ([Figure˜5(c)](#A5.F5.sf3)). During this phase, the training loss (corresponding to the RCE for 0 component (green line in the figure)) is not decreasing a lot.

Then, a phenomenon that we call "flattening" occurs: once the reconstruction error is small compared to the cutoff precision value, reconstruction error
flattens over the interval $[N_{\text{flat}},\textrm{r}_{\text{cutoff}}]$, where $N_{\text{flat}}$ is the smallest number such as

$$ $\textrm{RCE}^{S}_{N_{\text{flat}}}-\textrm{RCE}^{S}_{\textrm{r}_{\text{cutoff}}}\approx 0.$ (21) $$

Eventually, the phenomenon propagates toward low numbers of retained components ([Figure˜5(e)](#A5.F5.sf5)) and $N_{\text{flat}}=0$. Reconstruction error is now constant for all retained components and the training ends with training loss at cutoff precision. We refer a reader to [Section˜H.3](#A8.SS3) to have a more theoretical insight on what is happening during the flattening.

###### Remark 2 .

This phenomenon sheds light on the sharp drop in training loss observed near the end of optimization, as reported in Schwencke & Furtlehner (2025). By combining [Equations˜9](#S2.E9), [17](#S3.E17) and [14](#S2.E14) and using that $\Pi^{0}_{0}=0$, we obtain

$$ ${\textrm{RCE}^{S}_{0}}^{2}\stackrel{{\scriptstyle\ref{eqn:reconstruction-loss-orthogonal}}}{{=}}\frac{1}{S}\left\|\widehat{V}\Pi^{0}_{0}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}-\widehat{\nabla\mathcal{L}}\right\|^{2}_{\mathbb{R}^{S}}=\frac{1}{S}\left\|\widehat{\nabla\mathcal{L}}\right\|^{2}_{\mathbb{R}^{S}}\stackrel{{\scriptstyle\ref{eqn:eng-computation-trick}}}{{=}}\frac{1}{S}\sum_{i=1}^{S}\left(u_{{\bm{\theta}}}(x_{i})-f(x_{i})\right)^{2}\stackrel{{\scriptstyle\ref{eqn:empirical-quadratic-loss}}}{{=}}\ell({\bm{\theta}}).$ (22) $$

Thus, the last iteration of flattening is directly responsible for the sudden drop of train loss at the end of the training.

###### Remark 3 .

We see that for higher precision than the cutoff value ($N>\textrm{r}_{\text{cutoff}}$), the RCE is still decreasing as we increase the number of components kept. This indicates that there is still information to capture in the functional eigenspace composed of components associated to lower eigenvalues, see also [Section˜H.3](#A8.SS3).

The final interesting observation is that

$$ $\textrm{RCE}^{S}_{0}-\textrm{RCE}^{S}_{\textrm{r}_{\text{cutoff}}}\simeq 0\qquad\Leftrightarrow\qquad\Pi^{0}_{\textrm{r}_{\text{cutoff}}}\widehat{V}^{T}\widehat{\nabla\mathcal{L}}\approx 0.$ (23) $$

Thus, the flattening phenomenon means that the projection of the signal onto the first $\textrm{r}_{\text{cutoff}}$ components retained by the cutoff is negligible. In other words, the optimization has extracted all the *usable* signal from these components at this cutoff level.

###### Remark 2 .

###### Remark 3 .

### 3.3 Incomplete Flattening and Adaptive Strategies

In practice, for some experiments we observe that the flattening may remain incomplete with $\lim_{t\to\infty}N_{\text{flat}}=N_{\text{flat}}^{\infty}>0$: the system remains in a state similar to that shown in [Figure˜1(c)](#S3.F1.sf3) and never (at least not within a reasonable number of iterations) reaches the configuration illustrated in [Figure˜1(d)](#S3.F1.sf4).
A natural question arises: *what happens if we adjust the cutoff to retain exactly $N_{\text{flat}}^{\infty}$ components?*

If we try this trick in practice (see [Figure˜6](#A5.F6)), then a single natural gradient step with an adjusted cutoff can be enough to get immediate and complete flattening ($N_{\text{flat}}=0$) and eventually dramatically reduce training loss.
This abrupt flattening when restricting cutoff to low number of feature is typically accompanied by a learning rate found by the line search to be very close to one. A possible explanation is that this may represent an iteration in the *lazy training* regime (NTK and the feature matrix are nearly constant), where we regress linearly (and thus fast) based on learned features. This hypothesis should be further explored in future work.

This empirical insight motivates the use of an adaptive algorithm: by dynamically adjusting cutoffs, we can hope to accelerate convergence and achieve higher precision.

## 4 Algorithmic Design: Exploiting Flattening

Building upon the empirical analysis presented in [Section˜3](#S3), we develop a principled algorithm that controls and exploits the flattening phenomenon identified in ANaGRAM’s training dynamics. Our approach is based on tracking the relationship between reconstruction error and singular values to automatically determine well-adapted cutoff in order to reach the target precision (error) $\epsilon$ at the end of the training. This well-adapted cutoff should vary from one iteration to another to adjust to the currently learned weights and training dynamics in such a way to avoid early flattening (if flattening happens too early, the training stagnates at higher values of losses) and when intersection between RCE and singular values goes below the target precision $\epsilon$, we enforce the flattening, so that the final training loss also drops to $\epsilon$.

### 4.1 Adaptive Cutoff Strategy

In what follows, we suggest an adaptive cutoff rank $\textrm{r}_{\text{cutoff}}$ that indicates how much components of $\widehat{\Delta}$ are retained for the next update of ANaGRAM.
Our algorithm operates by dynamically selecting cutoff ranks based on the relationship between reconstruction error and singular values:

$$ $\textrm{r}_{\text{cutoff}}(t)=\begin{cases}\textrm{r}_{\text{int}}(t):=\max\left\{j:\textrm{RCE}^{S}_{j}(t)\leq\widehat{\Delta}_{t,j}\right\}&\text{if }\textrm{RCE}^{S}_{\textrm{r}_{\text{int}}(t)}(t)>\epsilon\text{ (intersection rank)},\\ \textrm{r}_{\epsilon}(t):=\max\left\{j:\textrm{RCE}^{S}_{j}(t)\geq\epsilon\right\}&\text{if }\textrm{RCE}^{S}_{\textrm{r}_{\text{int}}(t)}(t)\leq\epsilon\text{ (precision rank)}.\end{cases}$ (24) $$

The algorithm terminates when $\textrm{r}_{\epsilon}(t)=0$, indicating that the reconstruction error $\textrm{RCE}^{S}_{0}$ that is equal to the training error is indeed below the predefined precision threshold.

For ease of presentation, we provide only the core elements of AMStraMGRAM in Algorithm [1](#algorithm1) consisting in adaptively choosing, which $\textrm{r}_{\text{cutoff}}$ to apply for $\widehat{\Delta}$ at each update of ANaGRAM. The final algorithm is explained in [Appendix˜C](#A3). Final Algortihm [2](#algorithm2) addresses some irregularities observed in evolution of RCE and singular values that we explain in more details in [Section˜C.4](#A3.SS4).

Figure: (a) Early iterations ($\textrm{r}_{\text{cutoff}}=\textrm{r}_{\text{int}}$).
Refer to caption: /html/2510.15998/assets/figures/AMSTRAMGRAM_legend_horizontal.png

### 4.2 Geometrical Interpretation of the Adaptive Strategy

The algorithm exploits the geometric relationship between the empirical tangent space and the functional gradient. By tracking the intersection, we maximize the projection of the functional gradient onto the empirical tangent space while staying out of flattening. Once the intersection reach the precision level, we exploit the flattening phenomenon to achieve prescribed precision.

According to [˜1](#ThmProp1), the reconstruction error $\textrm{RCE}^{S}_{N}$ measures how much of the functional gradient signal remains to be captured by the first $N$ components. The intersection point thus represents the good balance between signal capture and phase transition.

## 5 Experimental Results

We first compare in [Table˜1](#S5.T1) our method implemented in JAX(^2^22[https://anonymous.4open.science/r/AMStraMGRAM-8D1B/](https://anonymous.4open.science/r/AMStraMGRAM-8D1B/)) with the ANaGRAM method (Schwencke & Furtlehner, 2025) on the benchmark problems presented in their paper, with modified datasets. As we see, for every equation, we perfom better.

**Table 1: Performance comparison between AMStraMGRAM (our method) and ANaGRAM Schwencke & Furtlehner (2025). The adaptive strategy demonstrates significant improvements across all benchmark problems, with L2 error improvements of up to 8 orders of magnitude.**
| Experiment | Mean Squared Error (MSE) | $L_{2}$ Error |  |  |
| --- | --- | --- | --- | --- |
|  | Ours | ANaGRAM | Ours | ANaGRAM |
| Heat Equation | 6.29e-29 $\pm$ 6.78e-30 | 8.56e-11 $\pm$ 7.05e-11 | 2.32e-14 $\pm$ 1.14e-14 | 1.28e-06 $\pm$ 1.75e-06 |
| Laplace 2D | 1.46e-28 $\pm$ 1.87e-29 | 4.27e-13 $\pm$ 4.66e-13 | 2.24e-15 $\pm$ 2.52e-16 | 3.49e-09 $\pm$ 3.58e-09 |
| Laplace 5D | 2.04e-08 $\pm$ 1.16e-08 | 6.37e-08 $\pm$ 7.01e-08 | 2.12e-05 $\pm$ 8.15e-06 | 4.00e-05 $\pm$ 2.93e-05 |
| Allen–Cahn | 3.19e-11 $\pm$ 2.37e-11 | 2.19e-04 $\pm$ 4.16e-04 | 5.87e-05 $\pm$ 6.25e-06 | 4.32e-03 $\pm$ 5.93e-03 |

We then compare our method with the baseline methods from Urbán et al. (2025) on the benchmark problems presented in their paper.
Note that in our case we do not need to enforce boundary constraints.
The methodology of sampling is also sighltly different, as we sample the data from a fixed grid, following the methodology of Schwencke & Furtlehner (2025), while in Urbán et al. (2025) they perform batching of randomly sampled points.

**Table 2: Performance comparison between AMStraMGRAM (our method) and baseline Urbán et al. (2025) methods. Our method demonstrates improvements across benchmark problems, without requiring enforcement of boundary constraints.**
| Experiment | Mean Squared Error (MSE) | $L_{2}$ Error |  |  |
| --- | --- | --- | --- | --- |
|  | Ours | SSBroyden* | Ours | SSBroyden* |
| One-dimensional Burgers (1DB) | 2.99e-12 $\pm$ 9.26e-13 | 2.92e-10 $\pm$ 1.45e-10 | 1.5e-06 $\pm$ 9.43e-7 | 1.59e-06 $\pm$ 1.02e-6 |
| Non-Linear Poisson (k=1) | 8.51e-24 $\pm$ 2.24e-24 | 3.03e-16 $\pm$ 3.82e-16 | 6.81e-10 $\pm$ 1.41e-09 | 9.29e-12 $\pm$ 5.85e-12 |
| Allen–Cahn (AC) | 3.19e-11 $\pm$ 2.37e-11 | 6.42e-12 $\pm$ 5.52e-12 | 5.87e-05 $\pm$ 6.25e-06 | 3.94e-06 $\pm$ 1.72e-06 |

## 6 Limitations

Despite its effectiveness, AMStraMGRAM can exhibit overfitting, particularly in problems with sharp features like the Allen–Cahn equation. The algorithm drives the training error to machine precision on the sampled points, but the learned function may develop high-frequency oscillations between them, especially in regions of high curvature where the approximation is the most challenging. These artifacts, visible as “overfitting lines” in [Figure˜3](#S6.F3), are an imprint of the sampling lattice (see regions around $x=\pm 0.5$). They arise because the SVD cutoff effectively projects the update onto a low-rank subspace of the tangent space. This subspace is often aligned with the grid axes, leading to anisotropic smoothing that perfectly fits the data on the grid lines but interpolates poorly in the under-sampled regions between them. Once the flattening phase begins, the training enters a quasi-linear regime that can “lock in” these geometric artifacts.

This phenomenon highlights that while our method significantly improves on ANaGRAM, the quality of the final solution remains fundamentally limited by the sampling strategy. Mitigating such overfitting requires co-designing the sampler and the optimizer. Potential remedies include adaptive sampling, where new collocation points are added in regions of high reconstruction error, or curriculum-based approaches that progressively refine the sampling grid.

Figure: Figure 3: Allen–Cahn overfitting: residual lines align with sampling lines. Low-rank (post-cutoff) tangent projections fit exactly on sampled fibers while interpolation between them inherits weakly constrained oscillations in regions of steep interface curvature.
Refer to caption: /html/2510.15998/assets/figures/loss-PDE-test.png

## 7 Conclusion

In this work, we have introduced AMStraMGRAM, an adaptive multi-cutoff strategy that enhances the ANaGRAM natural gradient method for training PINNs. Our work provides an analytical framework to explain ANaGRAM’s convergence behavior, uncovering a *flattening* phenomenon that clarifies its training dynamics. The proposed algorithm automatically adjusts cutoff regularization. Notably, AMStraMGRAM exhibits “overfitting” as demonstrated in Allen-Cahn experiments. These results underscore the potential of natural gradient optimization for PINNs while highlighting the critical role of sampling strategies in realizing their full accuracy.

Future research will focus on integrating residual-based methods to further stabilize training, establishing rigorous convergence guarantees for our adaptive cutoff scheme, and extending the approach to higher-dimensional PDEs and complex geometries. Exploring the interplay between network architecture and optimization—as well as further developing sampling techniques—will be essential to address the fundamental challenge of balancing optimization power with data representation. Ultimately, our findings suggest that with careful algorithmic design, PINNs can achieve the precision required for practical scientific computing, paving the way for mesh-free methods in computational science.

## References

- Adcock & Huybrechs (2019)
Ben Adcock and Daan Huybrechs.
Frames and Numerical Approximation.
*SIAM Review*, 61(3):443–473, January 2019.
ISSN 0036-1445, 1095-7200.
doi:[10.1137/17M1114697](https://doi.org/10.1137/17M1114697).
- Adcock & Huybrechs (2020)
Ben Adcock and Daan Huybrechs.
Frames and numerical approximation II: Generalized sampling, July 2020.
- Al-Baali (1998)
M. Al-Baali.
Numerical Experience with a Class of Self-Scaling Quasi-Newton Algorithms.
*Journal of Optimization Theory and Applications*, 96(3):533–553, March 1998.
ISSN 0022-3239, 1573-2878.
doi:[10.1023/A:1022608410710](https://doi.org/10.1023/A:1022608410710).
- Al-Baali & Khalfan (2005)
Mehiddin Al-Baali and Humaid Khalfan.
Wide interval for efficient self-scaling quasi-Newton algorithms.
*Optimization Methods and Software*, 20(6):679–691, December 2005.
ISSN 1055-6788, 1029-4937.
doi:[10.1080/10556780410001709448](https://doi.org/10.1080/10556780410001709448).
- Amari (1998)
Shun-Ichi Amari.
Natural gradient works efficiently in learning.
*Neural computation*, 10(2):251–276, 1998.
- Berezansky et al. (1996)
Yurij M. Berezansky, Zinovij G. Sheftel, and Georgij F. Us.
*Functional Analysis. Vol. II*, volume 86 of *Operator Theory Advances and Applications*.
Birkhäuser, 1996.
- Cuomo et al. (2022)
Salvatore Cuomo, Vincenzo Schiano Di Cola, Fabio Giampaolo, Gianluigi Rozza, Maziar Raissi, and Francesco Piccialli.
Scientific Machine Learning Through Physics–Informed Neural Networks: Where we are and What’s Next.
*Journal of Scientific Computing*, 92(3):88, July 2022.
ISSN 1573-7691.
doi:[10.1007/s10915-022-01939-z](https://doi.org/10.1007/s10915-022-01939-z).
- De Ryck et al. (2024)
Tim De Ryck, Florent Bonnet, Siddhartha Mishra, and Emmanuel de Bézenac.
An operator preconditioning perspective on training in physics-informed machine learning, May 2024.
- Duffy (2015)
Dean G. Duffy.
*Green’s Functions with Applications*.
Chapman and Hall/CRC, 2015.
- Grossmann et al. (2024)
Tamara G Grossmann, Urszula Julia Komorowska, Jonas Latz, and Carola-Bibiane Schönlieb.
Can physics-informed neural networks beat the finite element method?
*IMA Journal of Applied Mathematics*, 89(1):143–174, January 2024.
ISSN 0272-4960.
doi:[10.1093/imamat/hxae011](https://doi.org/10.1093/imamat/hxae011).
- Guzmán-Cordero et al. (2025)
Andrés Guzmán-Cordero, Felix Dangel, Gil Goldshlager, and Marius Zeinhofer.
Improving Energy Natural Gradient Descent through Woodbury, Momentum, and Randomization, May 2025.
- Jacot et al. (2018)
Arthur Jacot, Franck Gabriel, and Clément Hongler.
Neural tangent kernel: Convergence and generalization in neural networks.
*Advances in neural information processing systems*, 31, 2018.
- Jnini et al. (2024)
Anas Jnini, Flavio Vella, and Marius Zeinhofer.
Gauss-Newton Natural Gradient Descent for Physics-Informed Computational Fluid Dynamics, February 2024.
- Kirsch (2021)
Andreas Kirsch.
*An Introduction to the Mathematical Theory of Inverse Problems*, volume 120 of *Applied Mathematical Sciences*.
Springer International Publishing, Cham, 2021.
ISBN 978-3-030-63342-4 978-3-030-63343-1.
doi:[10.1007/978-3-030-63343-1](https://doi.org/10.1007/978-3-030-63343-1).
- Kiyani et al. (2025)
Elham Kiyani, Khemraj Shukla, Jorge F. Urbán, Jérôme Darbon, and George Em Karniadakis.
Which Optimizer Works Best for Physics-Informed Neural Networks and Kolmogorov-Arnold Networks?, April 2025.
- Kress (2014)
Rainer Kress.
*Linear Integral Equations*, volume 82 of *Applied Mathematical Sciences*.
Springer, New York, NY, 2014.
ISBN 978-1-4614-9592-5 978-1-4614-9593-2.
doi:[10.1007/978-1-4614-9593-2](https://doi.org/10.1007/978-1-4614-9593-2).
- Krishnapriyan et al. (2021)
Aditi Krishnapriyan, Amir Gholami, Shandian Zhe, Robert Kirby, and Michael W Mahoney.
Characterizing possible failure modes in physics-informed neural networks.
In *Advances in Neural Information Processing Systems*, volume 34, pp. 26548–26560. Curran Associates, Inc., 2021.
- Leshno et al. (1993)
Moshe Leshno, Vladimir Ya Lin, Allan Pinkus, and Shimon Schocken.
Multilayer feedforward networks with a nonpolynomial activation function can approximate any function.
*Neural networks*, 6(6):861–867, 1993.
- Liu & Nocedal (1989)
Dong C. Liu and Jorge Nocedal.
On the limited memory BFGS method for large scale optimization.
*Mathematical Programming*, 45(1):503–528, August 1989.
ISSN 1436-4646.
doi:[10.1007/BF01589116](https://doi.org/10.1007/BF01589116).
- Liu et al. (2024)
Songming Liu, Chang Su, Jiachen Yao, Zhongkai Hao, Hang Su, Youjia Wu, and Jun Zhu.
Preconditioning for Physics-Informed Neural Networks, February 2024.
- McClenny & Braga-Neto (2022)
Levi McClenny and Ulisses Braga-Neto.
Self-Adaptive Physics-Informed Neural Networks using a Soft Attention Mechanism, April 2022.
- Müller & Zeinhofer (2023)
Johannes Müller and Marius Zeinhofer.
Achieving high accuracy with PINNs via energy natural gradient descent.
In *International Conference on Machine Learning*, pp. 25471–25485. PMLR, 2023.
- Müller & Zeinhofer (2024)
Johannes Müller and Marius Zeinhofer.
Position: Optimization in SciML Should Employ the Function Space Geometry.
In *Forty-First International Conference on Machine Learning*, February 2024.
- Nocedal & Wright (1999)
Jorge Nocedal and Stephen J. Wright.
*Numerical Optimization*.
Springer, 1999.
- Ollivier (2015)
Yann Ollivier.
Riemannian metrics for neural networks I: Feedforward networks.
*Information and Inference: A Journal of the IMA*, 4(2):108–153, 2015.
- Paulsen & Raghupathi (2016)
Vern I. Paulsen and Mrinal Raghupathi.
*An Introduction to the Theory of Reproducing Kernel Hilbert Spaces*, volume 152.
Cambridge university press, 2016.
- Raissi et al. (2019)
M. Raissi, P. Perdikaris, and G.E. Karniadakis.
Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations.
*Journal of Computational Physics*, 378:686–707, February 2019.
ISSN 00219991.
doi:[10.1016/j.jcp.2018.10.045](https://doi.org/10.1016/j.jcp.2018.10.045).
- Rathore et al. (2024)
Pratik Rathore, Weimu Lei, Zachary Frangella, Lu Lu, and Madeleine Udell.
Challenges in Training PINNs: A Loss Landscape Perspective, February 2024.
- Schwencke & Furtlehner (2025)
Nilo Schwencke and Cyril Furtlehner.
ANaGRAM: A natural gradient relative to adapted model for efficient PINNs learning.
In *The Thirteenth International Conference on Learning Representations*, 2025.
- Urbán et al. (2025)
Jorge F. Urbán, Petros Stefanou, and José A. Pons.
Unveiling the optimization process of physics informed neural networks: How accurate and competitive can PINNs be?
*Journal of Computational Physics*, 523:113656, February 2025.
ISSN 0021-9991.
doi:[10.1016/j.jcp.2024.113656](https://doi.org/10.1016/j.jcp.2024.113656).
- Verbockhaven et al. (2024)
Manon Verbockhaven, Sylvain Chevallier, and Guillaume Charpiat.
Growing tiny networks: Spotting expressivity bottlenecks and fixing them optimally.
*arXiv preprint arXiv:2405.19816*, 2024.
- Wang et al. (2021)
Sifan Wang, Yujun Teng, and Paris Perdikaris.
Understanding and mitigating gradient flow pathologies in physics-informed neural networks.
*SIAM Journal on Scientific Computing*, 43(5):A3055–A3081, 2021.
- Wang et al. (2023)
Sifan Wang, Shyam Sankaran, Hanwen Wang, and Paris Perdikaris.
An Expert’s Guide to Training Physics-informed Neural Networks, August 2023.
- Wang et al. (2024)
Sifan Wang, Bowen Li, Yuhan Chen, and Paris Perdikaris.
PirateNets: Physics-informed Deep Learning with Residual Adaptive Networks.
*Journal of Machine Learning Research*, 25(402):1–51, 2024.
ISSN 1533-7928.
- Wu et al. (2023)
Chenxi Wu, Min Zhu, Qinyang Tan, Yadhu Kartha, and Lu Lu.
A comprehensive study of non-adaptive and residual-based adaptive sampling for physics-informed neural networks.
*Computer Methods in Applied Mechanics and Engineering*, 403:115671, 2023.

## Appendix A Illustration of Natural Gradient

Figure: Figure 4: Illustration of the orthogonal projection of the functional gradient onto the tangent space. While the ideal update direction would be the functional gradient $\nabla\mathcal{L}_{|{u_{{\bm{\theta}}}}}$ (shown in blue), our model constrains us to follow directions within the tangent space $T_{{\bm{\theta}}}\mathcal{M}$ (shown as a green plane). The optimal feasible direction is thus the orthogonal projection $\Pi^{\perp}_{T_{{\bm{\theta}}}\mathcal{M}}\left(\nabla\mathcal{L}_{|{u_{{\bm{\theta}}}}}\right)$ (shown in red).
Refer to caption: /html/2510.15998/assets/x1.png

## Appendix B Our vocabulary

- •
Domain ($\Omega$).
- •
Boundary ($\partial\Omega$).
- •
Differential operators ($D,B$).
- •
Cutoff ($\alpha_{t}$). A threshold below which the components of the matrix $\widehat{\Delta}$ are truncated, i.e. $\widehat{\Delta}\leftarrow\begin{cases}\widehat{\Delta}&\text{if }\widehat{\Delta}\geq\alpha_{t},\\
0&\text{else.}\end{cases}$
- •
Full rank ($\textrm{r}_{\text{svd}}$). A full rank of feature matrix $\widehat{\phi}$ that we assume, without loss of generality, to be equal to $\min(P,S)$.
- •
Rank ($\textrm{r}_{\text{cutoff}}$). A number of $\widehat{\Delta}$ components that are retained when computing a pseudo-inverse of $\widehat{\Delta}$ in ANaGRAM. Depending on a current regime of the training and a desired effect, it can be set at $\textrm{r}_{\text{int}}$ or $\textrm{r}_{\epsilon}$.
- •
Flattening. The phenomenon described in [Section˜3.2](#S3.SS2), when reconstruction error starts to stabilize for a range of possible ranks.
- •
Flat cutoff ($N_{\text{flat}}$). A number of components that corresponds to the beginning of flattening in reconstruction error curve.
- •
Feature matrix ($\widehat{\phi}\in\mathbb{R}^{P\times S}$). It is defined by a jacobian $\partial_{p}u_{{\bm{\theta}}}(x_{i})$, which is used in an ANaGRAM’s update to "project" a functional gradient onto parameter space of ${\bm{\theta}}$.
- •
Precision ($\epsilon$). A hyperparameter of AMStraMGRAM that prescribes a target error level that the algorithm should achieve.
- •
Intersection rank ($\textrm{r}_{\text{int}}$). Defined in [Equation˜24](#S4.E24), roughly speaking it corresponds to a number of components at which reconstruction error and singular values curves are intersecting.
- •
Precision rank ($\textrm{r}_{\epsilon}$). Defined in [Equation˜24](#S4.E24), it corresponds to a number of components at which reconstruction error curve and precision level are intersecting.
- •
Functional gradient ($\nabla\mathcal{L}$). A Frechet derivative of squared $L^{2}$ loss $\mathcal{L}$, its negative gives the "ideal" update direction in non-parametric case.
- •
Empirical functional gradient ($\widehat{\nabla\mathcal{L}}\in\mathbb{R}^{S}$). A vector obtain by evaluating $\nabla\mathcal{L}$ on some finite number of samples $x_{i}\in\Omega$, for $i\in 1,\dots,S$.
- •
Parametric model ($u_{{\bm{\theta}}}$). A function parametrized with ${\bm{\theta}}$ that serves to approximate a solution to a problem (regression or PDE). Typically, it is a neural network, where ${\bm{\theta}}$ are its full set of weights.
- •
Differential of the model ($d\,u_{{\bm{\theta}}}$). Defined as $d\,u_{{\bm{\theta}}}(h)=\sum_{p=1}^{P}h_{p}\frac{\partial u}{\partial{\bm{\theta}}_{p}}=\lim\limits_{\varepsilon\to 0}\frac{u_{|{\bm{\theta}}+\varepsilon h}-u_{{\bm{\theta}}}}{\varepsilon}$. It measures how much $u_{{\bm{\theta}}}$ changes in a given direction $h$.
- •
Tangent space ($T_{{\bm{\theta}}}\mathcal{M}$). Image of a differential of the model, giving a space of possible updates for a model $u_{{\bm{\theta}}}$.
- •
SVD components of $\widehat{\phi}$ ($\widehat{U}$, $\widehat{\Delta}$, $\widehat{V}$). In particular, $\widehat{\phi}=\widehat{U}\widehat{\Delta}\widehat{V}^{T}$, where $\widehat{U}\in\mathbb{R}^{P\times S}$ is a left singular vector matrix, $\widehat{\Delta}\in\mathbb{R}^{\textrm{r}_{\text{svd}}\times\textrm{r}_{\text{svd}}}$ is a diagonal matrix with singular values on a diagonal ordered in a decreasing order and $\widehat{V}$ is a right singular vector matrix.
- •
Functional singular vectors ($V_{t,{i}}$). Right singular vectors of the differential $du_{{\bm{\theta}}}$.
- •
Empirical tangent space ($T^{M}_{N}\mathcal{M}$). A subspace of tangent space $T_{{\bm{\theta}}}\mathcal{M}$, restricted to a span of the right functional singular vectors $V_{t,{i}}$ corresponding to a range of components from $M$ to $N$, i.e. $\operatorname*{Span}(V_{t,{i}}\,:\,1\leq M\leq N\leq N)$.
- •
Discretized empirical tangent space ($\widehat{T^{M}_{N}\mathcal{M}}$). A version of $T^{M}_{N}\mathcal{M}$ discretized on a set of samples $\{x_{i}\}_{i=1}^{S}$ coming from $\Omega$.
- •
Reconstruction error ($\textrm{RCE}^{S}_{N}$). A measure identifying the portion of the functional gradient signal that is lost when restricting $\widehat{\nabla\mathcal{L}}$ to $\widehat{T^{0}_{N}\mathcal{M}}$.
- •
Feature development phase. The early phase in the training, during which high volatility is observed in both quantities of interest with high sensitivity to the choice of $\textrm{r}_{\text{cutoff}}$.
- •
Flattening phase. The later phase in the training, during which reconstruction error starts to flatten for some values of $N$, at the same time singular values dominate over reconstruction error for all retained components, resulting in a drop of training loss.

## Appendix C Practical Implementation Considerations

While the principled algorithm discussed in the main paper and summarized in Algorithm [1](#algorithm1) provides a sound framework, empirical observations reveal that additional mechanisms are necessary for robust performance across diverse PDE problems. This section describes additional modifications to make the algorithm more practical.

### C.1 The Dual Cutoff Strategy: Addressing Empirical Challenges

Our experiments reveal that the single cutoff approach, while theoretically elegant, suffers from numerical instabilities and incomplete convergence in practice. We observed three critical issues:

- 1.
Ignition failure: The intersection between reconstruction error and singular values sometimes fails to evolve, preventing the algorithm from reaching lower error values.
- 2.
Retreating dynamics: The intersection rank may decrease during training, disrupting convergence.
- 3.
Incomplete flattening: Without additional stabilization, the flattening phenomenon may not complete, leading to suboptimal final accuracy.

To address these challenges, we introduce a dual cutoff strategy inspired by the staged design of rocket launches:

### C.2 Three-Phase Training Dynamics

#### C.2.1 Ignition Phase

We initialize two cutoffs:

- •
Minimum cutoff ($\textrm{r}_{\rm min}$): Set at the intersection point $\textrm{r}_{\text{int}}(t)$
- •
Maximum cutoff ($\textrm{r}_{\rm max}$): Set at the "elbow" of the singular value curve (see algorithm [3](#algorithm3))

The algorithm performs two natural gradient steps per iteration, one with each cutoff. If the intersection position remains static after both updates, we increment $\textrm{r}_{\rm max}$ by one to promote exploration of additional gradient components.

This phase ends when $\textrm{r}_{\rm min}$ reaches $\textrm{r}_{\rm max}$—an event we term liftoff.

#### C.2.2 Ascent Phase

During ascent, both cutoffs track the moving intersection, but with a stability mechanism:

$$ $\textrm{r}_{\rm max}(t)=\max(\textrm{r}_{\rm max}(t-1),\textrm{r}_{\text{int}}(t)).$ (25) $$

This monotonicity constraint prevents the intersection rank from falling to zero, which would disrupt training dynamics.

#### C.2.3 Stage Separation and Precision Locking

When $\textrm{RCE}^{S}_{\textrm{r}_{\text{int}}(t)}(t)\leq\epsilon$, we trigger stage separation:

- •
$\textrm{r}_{\rm min}$ is fixed at the precision level: $\textrm{r}_{\rm min}=\textrm{r}_{\epsilon}(t)$
- •
$\textrm{r}_{\rm max}$ continues tracking the intersection to maintain stability

The algorithm continues until $\textrm{r}_{\rm min}=0$ (booster return), indicating complete convergence.
The final algorithm that combines all three stages is mentioned in Algorithm [2](#algorithm2).

Figure: Algorithm 1 Sketch of the Adaptative MultiCutoff Strategy for ANaGRAM (AMStraMGRAM)

### C.3 Complete Practical Algorithm

Figure: Algorithm 2 AMStraMGRAM : Adaptive Multicutoff Strategy Modification for ANaGRAM

### C.4 Empirical Justification for Design Choices

The dual cutoff strategy addresses specific empirical challenges we observed:

Without the second cutoff, training dynamics sometimes stagnate. The dual approach provides both stability (via $\textrm{r}_{\rm min}$) and exploration (via $\textrm{r}_{\rm max}$).

The elbow point marks where singular values cease contributing meaningful signal, providing a natural upper bound for exploration.

Prevents catastrophic retreat of the intersection point, which we observed in complex equations like Allen-Cahn.

Triggered precisely when the intersection error drops below target precision, ensuring optimal utilization of the flattening phenomenon.

We see in the next section how this practical algorithm successfully improve empirical robustness.

## Appendix D Algorithmic details

Figure: Algorithm 3 Find elbow

Figure: Algorithm 4 Reconstruction Errors

## Appendix E Empirical example of Anagram Training Dynamics

In [Figure˜5](#A5.F5), we analyze ANaGRAM’s training on the heat equation with a fixed cutoff threshold $\alpha=10^{-3}$ and line search for the learning rate. The training loss coincides with $\left\|\widehat{\nabla\mathcal{L}}\right\|^{2}$. We can see the flattening phenomenon to occur on Iteration 120 and completed at 150. As discussed in the main paper, sometimes the flattening can be incomplete, and for many iterations remain without any further progress ($N_{\text{flat}}$ never reaching zero). In this case, changing a cutoff threshold results in an immediate and complete flattening for all first components up to $\textrm{r}_{\text{cutoff}}$, which is demonstrated in [Figure˜6](#A5.F6) for Iteration 120 of [Figure˜5](#A5.F5).

Figure: (a) Iteration 0: intersection point between singular values and reconstruction error lies before cutoff.
Refer to caption: /html/2510.15998/assets/figures/anagram_dynamic/heat_0001_0.png

Figure: (a) Same as [Figure˜5(d)](#A5.F5.sf4): iteration 120 of ANaGRAM with a fixed cutoff at $10^{-3}$.
Refer to caption: /html/2510.15998/assets/figures/instant_flattening/legend_crop.png

## Appendix F Deep dive on selected experiments

In this section we look at curves of training and estimations obtained with AMStraMGRAM on benchmark of PDEs.

### F.1 One Dimensional Burgers Equation

Figure: Figure 7: Training metrics for the One-Dimensional Burgers equation, showing convergence behavior with our adaptive multi-cutoff strategy.
Refer to caption: /html/2510.15998/assets/figures/1db_training_metrics.png

Figure: (a) True solution
Refer to caption: /html/2510.15998/assets/figures/burgers_true_solution_cutoff_1e-06.png

### F.2 Heat Equation

Figure: Figure 9: Convergence results for the Heat equation showing the $L_{2}$ error over iterations. Our method (AMStraMGRAM) converges faster and reaches a lower final error than ANaGRAM and baselines. Variability across runs is due to differing feature development speed from the random initialization.
Refer to caption: /html/2510.15998/assets/figures/AMSTRAMGRAM_HEAT_convergence_cutoff_1e-12_median_q1_q3.png

Figure: (a) True solution
Refer to caption: /html/2510.15998/assets/figures/HEAT_true_solution_cutoff_1e-14.png

### F.3 Laplace Equations (L2D and L5D)

For the Laplace equation in 2D, our method also demonstrates remarkable performance improvements over the baselines. The convergence is both faster and reaches a significantly lower error plateau.

Figure: Figure 11: Convergence results for the Laplace 2D problem, showing the $L_{2}$ error over iterations. Our method (AMStraMGRAM) achieves both faster convergence and lower final error compared to ANaGRAM and other baseline methods. The observed variance between runs can be explained by different speed of convergence depending on the initialization.
Refer to caption: /html/2510.15998/assets/figures/AMSTRAMGRAM_AC_convergence_cutoff_1e-12_median_q1_q3.png

Figure: (a) True solution
Refer to caption: /html/2510.15998/assets/figures/L2D_true_solution_cutoff_1e-14.png

Figure: Figure 13: Convergence results for the Laplace 5D problem, showing the $L_{2}$ error over iterations. Our method (AMStraMGRAM) achieves faster convergence but not lower final error compared to ANaGRAM and other baseline methods. We see that seeds change the speed of convergence of the algorithm
Refer to caption: /html/2510.15998/assets/figures/l5d_loss_l2_error.png

### F.4 Non Linear Poisson Equation

To compare ourselves with Urbán et al. (2025), we select (K=1).

Figure: Figure 14: Convergence results for the Non Linear Poisson equation, showing the $L_{2}$ error over iterations. Our method (AMStraMGRAM) achieves both faster convergence and lower final error compared to ANaGRAM and other baseline methods.
Refer to caption: /html/2510.15998/assets/figures/nlp_loss_l2_error.png

Figure: (a) True solution
Refer to caption: /html/2510.15998/assets/figures/NLP_true_solution_cutoff_1e-04.png

### F.5 Allen-Cahn Equation

Figure: Figure 16: Training curves for the Allen-Cahn equation, showing the evolution of loss and error over iterations.
Refer to caption: /html/2510.15998/assets/figures/AC_training_curves.png

Figure: (a) True solution
Refer to caption: /html/2510.15998/assets/figures/AC_solution-profile-cutoff-None-seed-42.png

## Appendix G Geometrical interpretation of regularizations

### G.1 Why Regularization is Necessary

We recall that our goal is to solve the operator equation $D[u]=f$ by minimizing the squared residual

$$ $\left\|D[u]-f\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}.$ (26) $$

For simplicity, assume $D$ is linear. Then the mapping

$$ $u\in\mathcal{C}^{\infty}(\Omega)\longmapsto\left\|D[u]\right\|_{\textrm{L}^{2}(\Omega,\mu)}$ (27) $$

defines a semi-norm on $\mathcal{C}^{\infty}(\Omega)$.
We can “upgrade” this semi-norm into a true norm by introducing the following generalized Sobolev norm:

$$ $\left\|\cdot\right\|_{\widetilde{\mathcal{H}}_{D}}:\left\{\begin{array}[]{lll}\mathcal{C}^{\infty}(\Omega\to\mathbb{R})&\to&\mathbb{R}^{+}\\ u&\mapsto&\sqrt{\left\|u\right\|^{2}_{\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)}+\left\|D{[u]}\right\|^{2}_{\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)}}\end{array}\right.$ (28) $$

Clearly, for any $u$,

$$ $\left\|u\right\|_{\textrm{L}^{2}(\Omega,\mu)}\leq\left\|u\right\|_{\widetilde{\mathcal{H}}_{D}},$ (29) $$

which guarantees that $\left\|\cdot\right\|_{\widetilde{\mathcal{H}}_{D}}$ is *definite*, i.e. $\left\|u\right\|_{\widetilde{\mathcal{H}}_{D}}=0\iff u=0$.

Completing $\mathcal{C}^{\infty}(\Omega)$ with respect to $\left\|\cdot\right\|_{\widetilde{\mathcal{H}}_{D}}$ yields a generalized Sobolev space $\big(\mathcal{H}_{D},\left\|\cdot\right\|_{\mathcal{H}_{D}}\big)$.
This Hilbert space is the largest subspace of $\textrm{L}^{2}(\Omega,\mu)$ on which $D$ is continuous.
Indeed, for every $u\in\mathcal{H}_{D}$,

$$ $\left\|D[u]\right\|_{\textrm{L}^{2}(\Omega,\mu)}\leq\left\|u\right\|_{\mathcal{H}_{D}}.$ (30) $$

Since our goal is to solve $D[u]=f$, we need $D$ to be continuously invertible.
That is, we need the reverse inequality of [Equation˜30](#A7.E30) to hold (up to a constant $\alpha>0$).
Formally, if $D$ were algebraically invertible (bijective as a mapping), this condition would read:

$$ $\displaystyle\begin{split}&\left(\exists\alpha>0,\,\forall u\in\mathcal{H}_{D},\,\left\|u\right\|_{\mathcal{H}_{D}}\leq\alpha\left\|D{[u]}\right\|_{\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)}\right)\\ \iff\ &\left(\exists\alpha>0,\,\forall u\in\mathcal{H}_{D},\,\left\|D^{-1}{[D{[u]}]}\right\|_{\mathcal{H}_{D}}\leq\alpha\left\|D{[u]}\right\|_{\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)}\right)\\ \iff\ &\left(\exists\alpha>0,\,\forall f\in\textrm{L}^{2}(\Omega\to\mathbb{R},\mu),\,\left\|D^{-1}{[f]}\right\|_{\mathcal{H}_{D}}\leq\alpha\left\|f\right\|_{\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)}\right)\end{split}.$ (31) $$

##### Operator ill-conditioning.

Even if $D$ is bijective, [Equation˜31](#A7.E31) may fail to hold, i.e. $D$ can be ill-conditioned.
Suppose there exists a subspace $\mathcal{H}_{K}\subset\mathcal{H}_{D}$ such that $D$ acts compactly on $\mathcal{H}_{K}$ with infinite rank.
Then $D$ admits a singular value decomposition (Kress, 2014, Theorem 15.16): for $u\in\mathcal{H}_{K}$,

$$ $D[u]=\sum_{n\in\mathbb{N}}e_{n}\lambda_{n}\left\langle v_{n}\,,\,u\right\rangle_{\mathcal{H}_{D}},$ (32) $$

with $(v_{n})$ orthonormal in $\mathcal{H}_{D}$, $(e_{n})$ orthonormal in $\textrm{L}^{2}(\Omega,\mu)$, and $\lambda_{n}\to 0$ as $n\to\infty$.

For [Equation˜31](#A7.E31) to hold, we would need $\inf_{n}\lambda_{n}>0$, contradicting $\lambda_{n}\to 0$.
This is exactly the classical inverse problem setting: $D$ is bijective but ill-conditioned, and regularization is unavoidable.
Among the many schemes developed, Tikhonov regularization is the canonical example (Kirsch, 2021).

##### Non-bijectivity.

If $D$ is not bijective, two additional issues may occur.

###### Non-surjectivity.

If $\operatorname{Im}D$ is a closed subspace, we can still obtain a solution by replacing the target $f$ with its projection $\Pi_{\operatorname{Im}D}f$.
Note that minimizing $\left\|D[u]-f\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}$ yields precisely this least-squares solution.

###### Non-injectivity.

The lack of injectivity is a much more subtle issue.
Since $D$ is linear and continuous, its null space $\operatorname{Ker}D$ is a closed subspace of $\mathcal{H}_{D}$.
In principle, one could restrict the domain of $D$ to $\operatorname{Ker}D^{\perp}$ to make it injective.
The problem, however, is that identifying $\operatorname{Ker}D$ is typically just as hard as solving the original problem itself, since it amounts to characterizing all $u\in\mathcal{H}_{D}$ such that $D[u]=0$.
Therefore, unless one can rely on theoretical results that explicitly describe $\operatorname{Ker}D$, or construct a subspace $\mathcal{H}_{0}\subset\mathcal{H}_{D}$ for which $\operatorname{Ker}D\cap\mathcal{H}_{0}$ is explicitly known (so that $D$ can be restricted to $\mathcal{H}_{0}$), it is generally impossible to “get rid of” $\operatorname{Ker}D$ in practice.

On the other hand, if we do not filter out $\operatorname{Ker}D$, this has the unwanted consequence of introducing “spurious” low-energy signals.
To be concrete, suppose we approximate our solution in a space $\mathcal{H}_{K}$ with orthonormal basis $(u_{n})_{n\in\mathbb{N}}$.
Assume there exists a subsequence $(u_{n}^{S})\notin\operatorname{Ker}D$ converging towards $\operatorname{Ker}D$.
Since $\operatorname{Ker}D$ is closed (by continuity of $D$), this means

$$ $\lim\limits_{n\to\infty}\left\|\Pi_{\operatorname{Ker}D}u_{n}^{S}-u_{n}^{S}\right\|_{\mathcal{H}_{D}}^{2}=0.$ (33) $$

Equivalently, after extraction, this can be rewritten for all $n\in\mathbb{N}$ as

$$ $\frac{\left\|\Pi_{\operatorname{Ker}D}u_{n}^{S}\right\|_{\mathcal{H}_{D}}^{2}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}^{2}}\;\geq\;1-2^{-n}.$ (34) $$

Now consider normalized vectors $u_{n}^{S}/\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}$.
We have

$$ $\displaystyle\begin{split}0<\|D\!\left[\tfrac{u_{n}^{S}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}}\right]\|_{\mathcal{H}_{D}}^{2}&=\left\|D\!\left[\tfrac{\Pi_{\operatorname{Ker}D^{\perp}}u_{n}^{S}+\Pi_{\operatorname{Ker}D}u_{n}^{S}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}}\right]\right\|_{\mathcal{H}_{D}}^{2}\\ &=\Bigg\|D\!\left[\tfrac{\Pi_{\operatorname{Ker}D^{\perp}}u_{n}^{S}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}}\right]+\underbracket{D\!\left[\tfrac{\Pi_{\operatorname{Ker}D}u_{n}^{S}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}}\right]}_{=0}\Bigg\|_{\mathcal{H}_{D}}^{2}\\ &=\left\|D\!\left[\tfrac{\Pi_{\operatorname{Ker}D^{\perp}}u_{n}^{S}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}}\right]\right\|_{\mathcal{H}_{D}}^{2}\\ &\stackrel{{\scriptstyle(\ref{eqn:fondamental-inequality-2})}}{{\leq}}\left\|\tfrac{\Pi_{\operatorname{Ker}D^{\perp}}u_{n}^{S}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}}\right\|_{\mathcal{H}_{D}}^{2}\\ &=1-\frac{\|\Pi_{\operatorname{Ker}D}u_{n}^{S}\|_{\mathcal{H}_{D}}^{2}}{\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}^{2}}\\ &\stackrel{{\scriptstyle(\ref{eqn:convenient-approaching-kernel-characterization})}}{{\leq}}2^{-n}.\end{split}$ (35) $$

In particular, if (for simplicity) the normalized $(u_{n}/\left\|u_{n}\right\|_{\mathcal{H}_{D}})$ are right singular vectors of $D$, then the vectors $\big(u_{n}^{S}/\left\|u_{n}^{S}\right\|_{\mathcal{H}_{D}}\big)$ will correspond to singular values vanishing at least as fast as $(2^{-n})$.
Crucially, however, these vanishing singular values do not reflect an intrinsic ill-conditioning of $D$, but rather an *artificial* ill-conditioning induced by the choice of approximation space $\mathcal{H}_{K}$.
In other words, the spurious instability arises from how we approximate the operator, not from the operator itself.
For more details on this approximation-induced phenomenon, see Adcock & Huybrechs (2019, 2020).

These remarks highlight the *inevitable need for regularization* in practice.
In the next section, we will provide a geometric interpretation of the two regularization schemes introduced in [Section˜2.5](#S2.SS5), emphasizing how fundamentally different they are in nature.

###### Remark 4 .

The above discussion becomes even more critical when we restrict ourselves to a finite-dimensional approximation space $\mathcal{H}_{\text{app}}\subset\mathcal{H}_{D}$.
In this case, the restriction $D_{\text{app}}$ is automatically compact, since it is of finite rank.
As a consequence, both types of ill-conditioning described above may occur simultaneously.
This highlights once again why regularization is not merely convenient but *unavoidable* in numerical practice.

###### Remark 4 .

### G.2 Ridge-regression

Returning to the definition given in [Section˜2.5](#S2.SS5), recall that Ridge regression amounts to adding $\alpha^{2}I_{d}$ (for some $\alpha>0$) to the Gram matrix $G_{\bm{\theta}}$ introduced in [Equation˜13](#S2.E13):

$$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,G_{{\bm{\theta}}_{t}}^{\dagger}\nabla\ell({\bm{\theta}}_{t})\,;$ $\displaystyle{G_{{\bm{\theta}}_{t}}}_{p,q}$ $\displaystyle:=\left\langle\partial_{p}u_{{\bm{\theta}}_{t}}\,,\,\partial_{q}u_{{\bm{\theta}}_{t}}\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}.$ ( 13 ) $$

We can reformulate this observation in the following way: given our model

$$ $u:\mathbb{R}^{P}\to\textrm{L}^{2}(\Omega\to\mathbb{R},\mu),$ (36) $$

consider the *regularized model*

$$ $u^{\alpha}:\left\{\begin{array}[]{lll}\mathbb{R}^{P}&\to&\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}\\ {\bm{\theta}}&\mapsto&(u_{{\bm{\theta}}},\alpha\,{\bm{\theta}}).\end{array}\right.$ (37) $$

The Gram matrix of this regularized model is exactly $G_{{\bm{\theta}}}+\alpha^{2}I_{d}$.
Suppose further that regression is performed with respect to some function $f\in\textrm{L}^{2}(\Omega,\mu)$. Then we must adapt the objective to the regularized model, replacing $f$ with the pair

$$ $(f,\alpha\,{\bm{\theta}})\in\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}.$ (38) $$

A straightforward computation shows that, for all $1\leq p\leq\min(P,S)$,

$$ $\displaystyle\begin{split}\left\langle\partial_{p}u_{{\bm{\theta}}}^{\alpha}\,,\,(f,\alpha\,{\bm{\theta}})-u^{\alpha}_{{\bm{\theta}}}\right\rangle_{\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}}&=\left\langle\partial_{p}u_{{\bm{\theta}}}\,,\,f-u_{{\bm{\theta}}}\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}+\underbracket{\alpha\left\langle{\bm{e}}^{(p)}\,,\,{\bm{\theta}}-{\bm{\theta}}\right\rangle_{\mathbb{R}^{P}}}_{=0}\\ &=\left\langle\partial_{p}u_{{\bm{\theta}}}\,,\,f-u_{{\bm{\theta}}}\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}.\end{split}$ (39) $$

Thus, regression of $(f,\alpha\,{\bm{\theta}})$ with the regularized model is exactly equivalent to Ridge regression.
Equivalently, Ridge regression corresponds to replacing the original model $u$ by the regularized model $u^{\alpha}$, and replacing the objective $f$ by $(f,\alpha\,{\bm{\theta}})$.
From this point of view, the choice of $\alpha\,{\bm{\theta}}$ as the secondary target may be interpreted as a *default assumption* in the absence of prior information on the parameters: one simply uses the current parameters as a reference target.

We can now extract several fundamental facts:

- 1.
As $\alpha\to 0$, the regularized model $u^{\alpha}$ tends in operator norm to the unregularized model $(u,0)$ (i.e. $u$ by abuse of notation). Indeed,
$\sup_{\left\|{\bm{\theta}}\right\|_{\mathbb{R}^{P}}=1}\left\|\textrm{d}u_{{\bm{\theta}}}^{\alpha}-(\textrm{d}u_{{\bm{\theta}}},0)\right\|_{\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}}=\alpha\sup_{\left\|{\bm{\theta}}\right\|_{\mathbb{R}^{P}}=1}\left\|{\bm{\theta}}\right\|_{\mathbb{R}^{P}}=\alpha.$
(40)
- 2.
The model $u^{\alpha}$ is injective and continuous. Since $\textrm{d}u_{{\bm{\theta}}}$ is continuous (as $\mathbb{R}^{P}$ is finite-dimensional), the only possible source of non-injectivity is $\operatorname{Ker}\textrm{d}u^{\alpha}_{{\bm{\theta}}}$. But
$\operatorname{Ker}\textrm{d}u^{\alpha}_{{\bm{\theta}}}=\operatorname{Ker}\textrm{d}u_{{\bm{\theta}}}\cap\operatorname{Ker}(\alpha I_{\mathbb{R}^{P}})\subset\operatorname{Ker}(\alpha I_{\mathbb{R}^{P}})=\{0\},$
(41)
hence injectivity. Restricting $u^{\alpha}$ to its image makes it algebraically bijective, and the inverse is continuous since
$\alpha\left\|{\bm{\theta}}\right\|_{\mathbb{R}^{P}}\leq\left\|\textrm{d}u^{\alpha}_{\bm{\theta}}\right\|_{\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}}.$
(42)
By the equivalence stated in [Equation˜31](#A7.E31), this implies that $\left(\textrm{d}u^{\alpha}_{\bm{\theta}}\right)^{-1}$ is continuous. Consequently, $\operatorname{Im}\textrm{d}u^{\alpha}_{\bm{\theta}}$ is closed in $\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}$, since it is the inverse image of a closed set under $\left(\textrm{d}u^{\alpha}_{\bm{\theta}}\right)^{-1}$. Therefore least-squares solution is well-defined.
- 3.
The least-squares solution of $u^{\alpha}=(f,0)$ is influenced by $\alpha$ as follows: $(f,0)$ is projected onto
$\operatorname{Im}\textrm{d}u^{\alpha}_{{\bm{\theta}}}=\operatorname*{Span}\big((\partial_{p}u_{{\bm{\theta}}},\alpha e^{(p)}):1\leq p\leq P\big).$
(43)
In particular, even if $f\in\operatorname{Im}\textrm{d}u_{{\bm{\theta}}}$ and $f\neq 0$, we still have $(f,0)\not\in\operatorname{Im}\textrm{d}u^{\alpha}_{{\bm{\theta}}}$ (since $\textrm{d}u_{{\bm{\theta}}}(0)=0$). Consequently,
$\left(\Pi^{\bot}_{\operatorname{Im}\textrm{d}u^{\alpha}_{{\bm{\theta}}}}(f,0)\right)_{1}\neq f,$
(44)
where the subscript $1$ denotes projection onto the first component in $\textrm{L}^{2}(\Omega,\mu)\times\mathbb{R}^{P}$.

We illustrate these phenomena in [Figure˜18(a)](#A7.F18.sf1).

Building on the above analysis, we now show that Ridge regression can be extended to the functional setting.
To this end, let us reconsider the operator $D:\mathcal{H}_{D}\to\textrm{L}^{2}(\Omega,\mu)$ introduced in [Section˜G.1](#A7.SS1).
Analogously to what we did for the parametric model $u$, we define the *regularized operator* at level $\alpha>0$ as

$$ $D^{\alpha}:\left\{\begin{array}[]{lll}\mathcal{H}_{D}&\to&\textrm{L}^{2}(\Omega,\mu)\times\mathcal{H}_{D}\\ u&\mapsto&(D[u],\alpha u)\end{array}\right..$ (45) $$

The corresponding target becomes the *regularized objective* $(f,\alpha u)$.

At this level of generality, the equivalence with Gram-matrix regularization no longer holds, since we are dealing with infinite-dimensional operators for which no direct Gram-matrix representation exists.
Nevertheless, the fundamental properties remain valid, namely:

- 1.
When $\alpha\to 0$, the regularized operator $D^{\alpha}$ converges to $(D,0)$ in the operator-norm sense, i.e. to $D$ by a mild abuse of notation. Indeed, we have
$\sup_{\left\|u\right\|_{\mathcal{H}_{D}}=1}\left\|D^{\alpha}{[u]}-(D,0)\right\|_{\textrm{L}^{2}(\Omega,\mu)\times\mathcal{H}_{D}}=\alpha\sup_{\left\|u\right\|_{\mathcal{H}_{D}}=1}\left\|u\right\|_{\mathcal{H}_{D}}=\alpha.$
(46)
- 2.
The operator $D^{\alpha}$ is injective and continuous.
Indeed, $D$ is continuous by the very construction of $\mathcal{H}_{D}$ (see [Section˜G.1](#A7.SS1)), and injectivity follows since
$\operatorname{Ker}D^{\alpha}=\operatorname{Ker}D\cap\operatorname{Ker}(\alpha I_{\mathcal{H}_{D}})\subseteq\operatorname{Ker}(\alpha I_{\mathcal{H}_{D}})=\{0\}.$
(47)
Restricting $D^{\alpha}$ to its image makes it algebraically bijective, and the inverse is continuous: we have $\alpha\left\|u\right\|_{\mathcal{H}_{D}}\leq\left\|D^{\alpha}{[u]}\right\|_{\textrm{L}^{2}(\Omega,\mu)\times\mathcal{H}_{D}}$, which by the equivalence in [Equation˜31](#A7.E31) implies that $\big(D^{\alpha}\big)^{-1}$ is continuous.
Consequently, $\operatorname{Im}D^{\alpha}$ is closed in $\textrm{L}^{2}(\Omega,\mu)\times\mathcal{H}_{D}$, since it is the inverse image of a closed set under $\left(D^{\alpha}\right)^{-1}$. Therefore least-squares solution is well-defined.
- 3.
Least-squares solutions of the regularized problem $D^{\alpha}[u]=(f,0)$ are impacted by $\alpha$ in the following way: we are projecting $(f,0)$ onto
$\operatorname{Im}D^{\alpha}=\operatorname*{Span}\Big(\,(D[h],\alpha h)\,:\,h\in\mathcal{H}_{D}\,\Big).$
(48)
In particular, even if $f\in\operatorname{Im}D$ with $f\neq 0$, we have $(f,0)\notin\operatorname{Im}D^{\alpha}$ (since $D[0]=0$), and hence
$\big(\Pi^{\bot}_{\operatorname{Im}D^{\alpha}}(f,0)\big)_{1}\neq f,$
(49)
where the subscript $1$ denotes the first coordinate in $\textrm{L}^{2}(\Omega,\mu)\times\mathcal{H}_{D}$.

We illustrate these phenomena in [Figure˜18(b)](#A7.F18.sf2).

In summary, Ridge regression can be interpreted as a modification of the operator $D$, rendering it injective and continuously invertible on its image.
However, this comes at a price: the regularized solutions are *never* exact solutions of the original equation $D[u]=f$, even when $\alpha$ is arbitrarily small, since we are in fact solving a different operator equation.
This marks a fundamental distinction from cutoff regularization, which instead acts directly on the approximation space, as we shall see in the next section.

Figure: (a) Illustration of parametric Ridge regression. The green region represents the solution space, while the blue regions denote the target spaces. As $\alpha\to 0$, the regularized graph $\Gamma_{\textrm{d}u_{\bm{\theta}}^{\alpha}}$ of $\textrm{d}u_{\bm{\theta}}^{\alpha}$ approaches the graph $\Gamma_{\textrm{d}u_{\bm{\theta}}}$ of $\textrm{d}u_{\bm{\theta}}$, with the angle between them vanishing at rate $\arctan(\alpha)$. The key consequence is that the projection of the objective $f$ onto $\operatorname{Im}\textrm{d}u_{\bm{\theta}}^{\alpha}$ follows a non-linear path as $\alpha\to 0$, coinciding with $\Pi_{\operatorname{Im}\textrm{d}u_{\bm{\theta}}}f$ only asymptotically.
Refer to caption: /html/2510.15998/assets/x2.png

### G.3 Cutoff regression

As in [Section˜G.2](#A7.SS2), let us return to the setting of [Section˜2.5](#S2.SS5).
In [Equation˜15](#S2.E15), we introduced cutoff regularization from the SVD perspective: given the differential $\textrm{d}u_{{\bm{\theta}}}$ of the model $u$, at the point ${\bm{\theta}}$, and its singular value decomposition $\textrm{d}u_{{\bm{\theta}}}=V_{\bm{\theta}}\Delta_{\bm{\theta}}U_{\bm{\theta}}^{T}$, the cutoff-regularized pseudo-inverse $\textrm{d}u^{\dagger_{\alpha}}_{{\bm{\theta}}}$ at level $\alpha>0$ is defined as

$$ $\displaystyle\textrm{d}u^{\dagger_{\alpha}}_{{\bm{\theta}}}$ $\displaystyle:=U_{\bm{\theta}}\Delta_{\bm{\theta}}^{\dagger_{\alpha}}V_{\bm{\theta}}^{T}\,;$ $\displaystyle\Delta_{{\bm{\theta}},{p}}^{\dagger_{\alpha}}$ $\displaystyle:=\begin{cases}\Delta_{{\bm{\theta}},{p}}^{-1}&\text{if }\Delta_{{\bm{\theta}},{p}}\geq\alpha\\ 0&\text{otherwise}\end{cases},\,1\leq p\leq P.$ (50) $$

Let us reinterpret this construction.
Denote by $N_{\alpha}\in\mathbb{N}$ the number of singular values larger than $\alpha$.
Equivalently, assuming $(\Delta_{{\bm{\theta}},{p}})_{1\leq p\leq P}$ is non-increasing,

$$ $N_{\alpha}:=\operatorname*{arg\,max}_{p\in\mathbb{N}}\{\,\Delta_{{\bm{\theta}},{p}}\geq\alpha\,\}.$ (51) $$

Define

$$ $\Theta_{\alpha}:=\operatorname*{Span}\{U_{{\bm{\theta}},{p}}:1\leq p\leq N_{\alpha}\},\qquad T^{0}_{N_{\alpha}}\mathcal{M}:=\operatorname*{Span}\{V_{t,{p}}:1\leq p\leq N_{\alpha}\},$ (52) $$

so that $T^{0}_{N_{\alpha}}\mathcal{M}=\textrm{d}u_{\bm{\theta}}\left(\Theta_{\alpha}\right)$.
We then have

$$ $\left(\textrm{d}u_{{\bm{\theta}}_{|\Theta_{\alpha}}}^{|T^{0}_{N_{\alpha}}\mathcal{M}}\right)^{-1}=\textrm{d}u^{\dagger_{\alpha}}_{{\bm{\theta}}},$ (53) $$

meaning that the restriction $\textrm{d}u_{{\bm{\theta}}}^{\alpha}:=\textrm{d}u_{{\bm{\theta}}_{|\Theta_{\alpha}}}$ of $\textrm{d}u_{{\bm{\theta}}}$ to the domain $\Theta_{\alpha}$ becomes invertible once its codomain is restricted to its image $T^{0}_{N_{\alpha}}\mathcal{M}$, with inverse given precisely by the cutoff pseudo-inverse $\textrm{d}u^{\dagger_{\alpha}}_{{\bm{\theta}}}$.
Moreover, for any $h\in\Theta_{\alpha}$,

$$ $\left\|\textrm{d}u_{{\bm{\theta}}}(h)\right\|_{\textrm{L}^{2}(\Omega,\mu)}=\left\|V_{\bm{\theta}}\Delta_{\bm{\theta}}U_{\bm{\theta}}^{T}h\right\|_{\textrm{L}^{2}(\Omega,\mu)}\stackrel{{\scriptstyle V_{\bm{\theta}}\text{ unitary}}}{{=}}\left\|\Delta_{\bm{\theta}}U_{\bm{\theta}}^{T}h\right\|_{\mathbb{R}^{P}}\stackrel{{\scriptstyle h\in\Theta_{\alpha}}}{{\geq}}\alpha\left\|U_{\bm{\theta}}^{T}h\right\|_{\mathbb{R}^{P}}\stackrel{{\scriptstyle U_{\bm{\theta}}\text{ unitary}}}{{=}}\alpha\left\|h\right\|_{\mathbb{R}^{P}}.$ (54) $$

In other words, [Equation˜31](#A7.E31) is satisfied by $\textrm{d}u_{{\bm{\theta}}}^{\alpha}$.

Thus, while ridge regularization modifies the model itself, cutoff regularization instead restricts the domain of the model so that, on this restricted domain, [Equation˜31](#A7.E31) holds and the model becomes invertible.
We summarize the fundamental properties:

- 1.
We have
$\bigcap_{\alpha>0}\left(\mathbb{R}^{P}\backslash\Theta_{\alpha}\right)=\operatorname{Ker}\textrm{d}u_{{\bm{\theta}}},$
(55)
that is, $\lim_{\alpha\to 0}\mathbb{R}^{P}\backslash\Theta_{\alpha}=\operatorname{Ker}\textrm{d}u_{{\bm{\theta}}}$, since for all $\alpha>\beta$ we have $\Theta_{\alpha}\subset\Theta_{\beta}$ and then $\mathbb{R}^{P}\backslash\Theta_{\beta}\subset\mathbb{R}^{P}\backslash\Theta_{\alpha}$.
Similarly, $\lim_{\alpha\to 0}T^{0}_{N_{\alpha}}\mathcal{M}=\operatorname{Im}\textrm{d}u_{{\bm{\theta}}}$.
Moreover, for each $\alpha>0$, the restriction $\textrm{d}u_{{\bm{\theta}}}^{\alpha}$ coincides with $\textrm{d}u_{{\bm{\theta}}}$ on $\Theta_{\alpha}$.
- 2.
By [Equation˜54](#A7.E54), $\textrm{d}u_{{\bm{\theta}}}^{\alpha}$ is injective and continuous.
Restricting it to its image $T^{0}_{N_{\alpha}}\mathcal{M}$ makes it bijective and bicontinuous, with inverse exactly the cutoff pseudo-inverse $\textrm{d}u^{\dagger_{\alpha}}_{{\bm{\theta}}}$.
In particular $\textrm{d}u\left(\Theta_{\alpha}\right)$ is closed in $\textrm{L}^{2}(\Omega,\mu)$, since it is the inverse image of a closed set under $\textrm{d}u^{\dagger_{\alpha}}_{{\bm{\theta}}}$. Therefore least-squares solution is well-defined.
- 3.
Solving the least-squares problem $\textrm{d}u_{{\bm{\theta}}}^{\alpha}=f$ is now altered in the following way: the target $f$ is first projected onto $T^{0}_{N_{\alpha}}\mathcal{M}=\operatorname{Im}\textrm{d}u_{{\bm{\theta}}}^{\alpha}$.
In particular, if for some $\alpha>0$ we already have $f\in\operatorname{Im}\textrm{d}u_{{\bm{\theta}}}^{\alpha}$, then the regularized least-squares formulation recovers an *exact solution* to the problem. This stands in sharp contrast with Ridge regression, where such exact recovery can only occur *asymptotically* in the limit $\alpha\to 0$.

As in [Section˜G.2](#A7.SS2), we now need to reinterpret the cutoff regularization in order to extend it to the functional setting.
Let us return once more to the operator $D:\mathcal{H}_{D}\to\textrm{L}^{2}(\Omega,\mu)$ introduced in [Section˜G.1](#A7.SS1).
In general, one cannot define an SVD for such an operator (except when it is compact).
We must therefore appeal to the spectral theorem for bounded self-adjoint operators, which relies on the notion of a *projection-valued measure* (also called a resolution of the identity).
For our purposes, it will be sufficient to simply state the definition.

###### Definition 1 (Projection-valued measure) .

Let $(X,\mathcal{A})$ be a measurable space, where $\mathcal{A}$ denotes its $\sigma$-algebra, and let $\mathcal{H}$ be a Hilbert space.
A *projection-valued measure* (PVM) is a map

$$ $\pi:\mathcal{A}\to\mathcal{L}_{b}(\mathcal{H}\to\mathcal{H}),$ $$

where $\mathcal{L}_{b}(\mathcal{H}\to\mathcal{H})$ denotes the set of bounded operators on $\mathcal{H}$,
such that for every $A\in\mathcal{A}$, $\pi(A)$ is an orthogonal projection on $\mathcal{H}$, and the following properties hold:

- 1.
$\pi(\emptyset)=0$ and $\pi(X)=I_{\mathcal{H}}$, where $I_{\mathcal{H}}$ is the identity operator on $\mathcal{H}$;
- 2.
$\pi(A\cap B)=\pi(A)\pi(B)$ for all $A,B\in\mathcal{A}$;
- 3.
For every countable family $(A_{i})_{i=1}^{\infty}$ of disjoint sets in $\mathcal{A}$,
$\pi\!\left(\bigcup_{i=1}^{\infty}A_{i}\right)=\sum_{i=1}^{\infty}\pi(A_{i}),$
where the series converges in the strong operator topology.

Since projection-valued measures are measures, one can define integrals with respect to them.
We refer to (Berezansky et al., 1996, Chapter 13) for details.
We may now state the spectral theorem.

###### Theorem 2 .

Let $\mathcal{H}$ be a Hilbert space and let $A:\mathcal{H}\to\mathcal{H}$ be a self-adjoint operator.
Then there exists a projection-valued measure $\pi$ on the Borel $\sigma$-algebra of $\mathbb{R}$ such that

$$ $A=\int_{\mathbb{R}}\lambda\,\pi(\textrm{d}\lambda)=\int_{\sigma(A)}\lambda\,\pi(\textrm{d}\lambda),$ (56) $$

where $\sigma(A)$ denotes the spectrum of $A$.

A proof can be found in (Berezansky et al., 1996, Theorem 4.1, Section 4.1, Chapter 13).
In particular, since $\pi$ is a projection-valued measure, we have by [˜1](#ThmDef1):

$$ $I_{\mathcal{H}}=\int_{\mathbb{R}}\pi(\textrm{d}\lambda).$ (57) $$

Since $D$ is continuous, we can define its adjoint $D^{*}:\textrm{L}^{2}(\Omega,\mu)\to\mathcal{H}_{D}$, and hence the self-adjoint operator $D^{*}D:\mathcal{H}_{D}\to\mathcal{H}_{D}$.
Applying [˜2](#ThmTh2), we obtain a projection-valued measure $\pi_{D}$ on $\mathbb{R}$ endowed with its Borel $\sigma$-algebra, such that

$$ $\displaystyle D^{*}D$ $\displaystyle=\int_{\mathbb{R}_{+}}\lambda\,\pi_{D}(\textrm{d}\lambda),$ $\displaystyle I_{\mathcal{H}_{D}}$ $\displaystyle=\int_{\mathbb{R}_{+}}\pi_{D}(\textrm{d}\lambda),$ (58) $$

where the integration is restricted to $\mathbb{R}_{+}$ since $D^{*}D$ is a positive operator.
We can then define

$$ $\Pi_{D}^{\alpha}:=\int_{\alpha^{2}}^{+\infty}\pi_{D}(\textrm{d}\lambda),$ (59) $$

which is an orthogonal projection in $\mathcal{H}_{D}$ since $\pi_{D}$ is a projection-valued measure.
We then define the regularized space $\mathcal{H}_{D}^{\alpha}$ at level $\alpha>0$ by

$$ $\mathcal{H}_{D}^{\alpha}:=\operatorname{Im}\Pi_{D}^{\alpha}\subset\mathcal{H}_{D}.$ (60) $$

For any $u\in\mathcal{H}_{D}^{\alpha}$, we compute

$$ $\displaystyle\begin{split}\left\|D{[u]}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}&=\left\langle D[u]\,,\,D[u]\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}=\left\langle u\,,\,D^{*}D[u]\right\rangle_{\mathcal{H}_{D}}\\ &=\left\langle u\,,\,\int_{\mathbb{R}_{+}}\lambda\pi_{D}(\textrm{d}\lambda)u\right\rangle_{\mathcal{H}_{D}}\\ &\stackrel{{\scriptstyle u\in\mathcal{H}_{D}^{\alpha}}}{{=}}\left\langle u\,,\,\int_{\mathbb{R}_{+}}\lambda_{1}\pi_{D}(\textrm{d}\lambda_{1})\int_{\alpha^{2}}^{+\infty}\pi_{D}(\textrm{d}\lambda_{2})u\right\rangle_{\mathcal{H}_{D}}\\ &\stackrel{{\scriptstyle\pi_{D}\text{ PVM}}}{{=}}\left\langle u\,,\,\int_{\alpha^{2}}^{+\infty}\lambda\pi_{D}(\textrm{d}\lambda)u\right\rangle_{\mathcal{H}_{D}}\\ &=\int_{\alpha^{2}}^{+\infty}\lambda\left\langle u\,,\,\pi_{D}(\textrm{d}\lambda)u\right\rangle_{\mathcal{H}_{D}}\\ &\geq\alpha^{2}\int_{\alpha^{2}}^{+\infty}\left\langle u\,,\,\pi_{D}(\textrm{d}\lambda)u\right\rangle_{\mathcal{H}_{D}}\stackrel{{\scriptstyle u\in\mathcal{H}_{D}^{\alpha}}}{{=}}\alpha^{2}\left\langle u\,,\,u\right\rangle_{\mathcal{H}_{D}}=\alpha^{2}\left\|u\right\|_{\mathcal{H}_{D}}^{2}.\end{split}$ (61) $$

That is,

$$ $\left\|D{[u]}\right\|_{\textrm{L}^{2}(\Omega,\mu)}\;\geq\;\alpha\left\|u\right\|_{\mathcal{H}_{D}},$ (62) $$

so that [Equation˜31](#A7.E31) is verified on $\mathcal{H}_{D}^{\alpha}$.
We denote

$$ $D^{\alpha}:=D_{|\mathcal{H}_{D}^{\alpha}},$ (63) $$

the restriction of $D$ to the domain $\mathcal{H}_{D}^{\alpha}$.
We can now list the fundamental properties:

- 1.
We have
$\bigcap_{\alpha>0}\left(\mathcal{H}_{D}\backslash\mathcal{H}_{D}^{\alpha}\right)=\operatorname{Ker}D$
(64)
that is, $\lim_{\alpha\to 0}\mathcal{H}_{D}\backslash\mathcal{H}_{D}^{\alpha}=\operatorname{Ker}D$, since for all $\alpha>\beta$, $\mathcal{H}_{D}^{\alpha}\subset\mathcal{H}_{D}^{\beta}$ by Property 3 of [˜1](#ThmDef1).
Moreover, by continuity of $D$, we also have $\lim_{\alpha\to 0}D\left[\mathcal{H}_{D}^{\alpha}\right]=\operatorname{Im}D$.
Finally, for each $\alpha>0$, $D^{\alpha}$ coincides with $D$ on $\mathcal{H}_{D}^{\alpha}$ by construction.
- 2.
As established by [Equation˜31](#A7.E31), $D^{\alpha}$ is injective and continuous. When restricted to its image, it is therefore bijective and bicontinuous, hence invertible.
In particular $D\left[\mathcal{H}_{D}^{\alpha}\right]$ is closed in $\textrm{L}^{2}(\Omega,\mu)$, since it is the inverse image of a closed set under $\left(D^{\alpha}\right)^{-1}$. Therefore least-squares solution is well-defined.
- 3.
The least-squares solution of $D^{\alpha}=f$ is now modified as follows: one projects $f$ onto $\operatorname{Im}D^{\alpha}$.
In particular, if for some $\alpha>0$ we already have $f\in\operatorname{Im}D^{\alpha}$, then the regularized least-squares formulation recovers an *exact solution* to the problem $D[u]=f$. This stands in sharp contrast with Ridge regression, where such exact recovery can only occur *asymptotically* in the limit $\alpha\to 0$.

###### Definition 1 (Projection-valued measure) .

###### Theorem 2 .

### G.4 Connection to Green’s Function

To further highlight the difference between the two regularization schemes, we now reinterpret them through the lens of Green’s functions of the operator $D$.
Schwencke & Furtlehner (2025, Theorem 2) established in the finite-dimensional case a connection between the natural gradient for PINNs and Green’s functions. Their proof relies on Schwencke & Furtlehner (2025, Proposition 3), which will be our starting point. We restate the relevant definitions and results for completeness.

###### Definition 2 ( Schwencke & Furtlehner, 2025 , Definition 9 : generalized Green’s function) .

Let $\mathcal{H}$ be an Hilbert space, $D:\mathcal{H}\to\textrm{L}^{2}(\Omega,\mu)$ be a linear differential operator, $\mathcal{H}_{0}\subset\mathcal{H}$ a subspace isometrically embedded in $\mathcal{H}$ and $f\in\textrm{L}^{2}(\Omega,\mu)$. A generalized Green’s function of $D$ on $\mathcal{H}_{0}$ is then any kernel function $g:\Omega\times\Omega\to\mathbb{R}$ such that the operator:

$$ $R_{\mathcal{H}_{0}}:\left\{\begin{array}[]{lll}\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)&\to&\mathcal{H}\\ f&\mapsto&\Big(x\in\Omega\mapsto\int_{\Omega}g(x,s)f(s)\mu(\textrm{d}s)\Big)\end{array}\right.,$ $$

verifies the equation:

$$ $D\circ R_{\mathcal{H}_{0}}=\Pi^{\bot}_{D[\mathcal{H}_{0}]}$ (65) $$

###### Proposition 2 ( Schwencke & Furtlehner, 2025 , Proposition 3 ) .

Let $D:\mathcal{H}\to\textrm{L}^{2}(\Omega,\mu)$ be a linear differential operator, and $\mathcal{H}_{0}:=\operatorname*{Span}(u_{p}:1\leq p\leq P)\subset\mathcal{H}$ a subspace isometrically embedded in $\mathcal{H}$. Then the generalized Green’s function of $D$ on $\mathcal{H}_{0}$ is given by: $\textrm{for all }x,y\in\Omega$

$$ $g_{\mathcal{H}_{0}}(x,y):=\sum_{1\leq p,q\leq P}u_{p}(x)\,G^{\dagger}_{p,q}D[u_{q}](y),$ (66) $$

with: $\textrm{for all }1\leq p,q\leq P$,

$$ $G_{p,q}:=\left\langle D[u_{p}]\,,\,D[u_{q}]\right\rangle_{\textrm{L}^{2}(\Omega\to\mathbb{R},\mu)}.$ (67) $$

###### Definition 2 ( Schwencke & Furtlehner, 2025 , Definition 9 : generalized Green’s function) .

###### Proposition 2 ( Schwencke & Furtlehner, 2025 , Proposition 3 ) .

##### Our goal.

We aim to

- (i)
generalize Schwencke & Furtlehner (2025, Proposition 3) to arbitrary Reproducing Kernel Hilbert Spaces;
- (ii)
establish a direct connection to the regularization framework introduced earlier. This will provide a novel reinterpretation of the Green’s function in the regularized operator setting.

##### Operator framework.

Consider the operator $D:\mathcal{H}_{D}\to\textrm{L}^{2}(\Omega,\mu)$ from [Section˜G.1](#A7.SS1), and assume that there exists an RKHS $\mathcal{H}_{0}$ isometrically embedded in $\mathcal{H}_{D}$ (for instance, any finite-dimensional RKHS, see Schwencke & Furtlehner, 2025, Corollary 1).
For Schwencke & Furtlehner (2025, Definition 9) to be well-posed, the range $D[\mathcal{H}_{0}]$ must be a closed subspace of $\textrm{L}^{2}(\Omega,\mu)$. As argued earlier, this is guaranteed if $D$ is continuously invertible: indeed, in this case

$$ $D[\mathcal{H}_{0}]=(D^{-1})^{-1}[\mathcal{H}_{0}],$ (68) $$

and the inverse image of a closed subspace under a continuous operator is closed.

##### Key observation.

Thus, to generalize Schwencke & Furtlehner (2025, Proposition 3), we require $D$ to be continuously invertible. Conveniently, this is precisely the property enforced by the regularization schemes we introduced earlier.

In what follows, we first focus on the cutoff regularization, which offers the clearest interpretation in terms of Green’s functions. We then briefly revisit the case of Ridge regression.
Before delving further into our main goal, let us first establish two general facts.

###### Lemma 1 .

Let $\big(\mathcal{H}_{0},\left\|\cdot\right\|_{\mathcal{H}_{0}}\big)$ be an RKHS on a set $X$ with reproducing kernel $k$. Suppose that $\left\|\cdot\right\|_{bis}$ is a norm equivalent to $\left\|\cdot\right\|_{\mathcal{H}_{0}}$. Then $\big(\mathcal{H}_{0},\left\|\cdot\right\|_{bis}\big)$ is also an RKHS.

###### Proof.

The key point is to show that there exists a reproducing kernel for the inner product $\left\langle\cdot\,,\,\cdot\right\rangle_{bis}$ associated with $\left\|\cdot\right\|_{bis}$.
Our argument follows the simple reasoning in Paulsen & Raghupathi (2016, Definitions 1–2).

Since, for every $x\in X$, the point evaluation functional

$$ $\delta_{x}:u\in\mathcal{H}_{0}\mapsto u(x)$ (69) $$

is continuous with respect to $\left\|\cdot\right\|_{\mathcal{H}_{0}}$ by the definition of an RKHS, it is also continuous with respect to the equivalent norm $\left\|\cdot\right\|_{bis}$
Therefore, by the Riesz representation theorem, for each $x\in X$, there exists a unique element $k^{bis}_{x}\in\mathcal{H}_{0}$ such that for all $u\in\mathcal{H}_{0}$

$$ $\left\langle k^{bis}_{x}\,,\,u\right\rangle_{bis}=u(x).$ (70) $$

In particular, this defines a reproducing kernel for the norm $\left\|\cdot\right\|_{bis}$, given by

$$ $k_{bis}(x,y)=\left\langle k^{bis}_{x}\,,\,k^{bis}_{y}\right\rangle_{bis}=k^{bis}_{x}(y),\qquad\forall x,y\in X.$ (71) $$

Hence $\big(\mathcal{H}_{0},\left\|\cdot\right\|_{bis}\big)$ is indeed an RKHS.
∎

###### Lemma 2 .

Let $\mathcal{H}_{A},\mathcal{H}_{B}$ be two Hilbert spaces.
If $U:\mathcal{H}_{\mathcal{A}}\to\mathcal{H}_{\mathcal{B}}$ is an isometry, then

$$ $\displaystyle U^{*}U$ $\displaystyle=\textrm{I}_{\mathcal{H}_{\mathcal{A}}},\quad$ $\displaystyle UU^{*}$ $\displaystyle=\Pi_{\operatorname{Im}U}.$ (72) $$

In particular $\operatorname{Im}U$ is closed in $\mathcal{H}_{\mathcal{B}}$.

###### Proof.

The first identity follows from the fact that for all $x,y\in\mathcal{H}_{\mathcal{A}}$,

$$ $\left\langle x\,,\,U^{*}U[y]\right\rangle_{\mathcal{H}_{\mathcal{A}}}=\left\langle U[x]\,,\,U[y]\right\rangle_{\mathcal{H}_{\mathcal{B}}}=\left\langle x\,,\,y\right\rangle_{\mathcal{H}_{\mathcal{A}}}.$ (73) $$

Thus $\left(U^{*}U(y)-y\right)\in\mathcal{H}_{\mathcal{A}}^{\bot}$, i.e. $U^{*}U=\textrm{I}_{\mathcal{H}_{\mathcal{A}}}$.
For the second, the key point is to show that $\operatorname{Im}U$ is closed, i.e. $\operatorname{Im}U=\overline{\operatorname{Im}U}$.

Let $y\in\overline{\operatorname{Im}U}$, and $(y_{n})\in\operatorname{Im}U^{\mathbb{N}}$ with $y_{n}\to y$.
Since $(y_{n})$ is Cauchy, and $y_{n}=U(x_{n})$ for some $(x_{n})\in\mathcal{H}_{\mathcal{A}}^{\mathbb{N}}$, we have

$$ $\left\|U(x_{n})-U(x_{m})\right\|_{\mathcal{H}_{\mathcal{B}}}=\left\|x_{n}-x_{m}\right\|_{\mathcal{H}_{\mathcal{A}}},$ (74) $$

so $(x_{n})$ is also Cauchy and converges to $x\in\mathcal{H}_{\mathcal{A}}$, since $\mathcal{H}_{\mathcal{A}}$ is complete.
Since $U$ is an isometry, we have for all $x\in\mathcal{H}_{\mathcal{A}}$

$$ $\left\|U(x)\right\|_{\mathcal{H}_{\mathcal{B}}}=\left\|x\right\|_{\mathcal{H}_{\mathcal{A}}}.$ (75) $$

In particular, $U$ is bounded with operator norm $\left\|U\right\|=1$, and hence continuous.
Thus $U(x)=y$, hence $y\in\operatorname{Im}U$. We conclude that $\operatorname{Im}U$ is closed in $\mathcal{H}_{\mathcal{B}}$.
Finally:

- •
For $y\in\operatorname{Im}U$, say $y=U(x)$, we have
$UU^{*}(y)=U(U^{*}U)(x)=U(x)=y.$
(76)
- •
For $y\in(\operatorname{Im}U)^{\perp}$, we check that $UU^{*}(y)=0$. Indeed, for any $z\in\mathcal{H}_{\mathcal{B}}$,
$\left\langle z\,,\,UU^{*}(y)\right\rangle_{\mathcal{H}_{\mathcal{B}}}=\left\langle UU^{*}(z)\,,\,y\right\rangle_{\mathcal{H}_{\mathcal{B}}}=0,$
(77)
since $UU^{*}(z)\in\operatorname{Im}U$. Thus $UU^{*}(y)\in\mathcal{H}_{\mathcal{B}}^{\bot}$, i.e. $UU^{*}(y)=0$.

∎

We are interested in the restriction of $D$ to the domain $\mathcal{H}_{0}$.
Since the restriction $D^{*}D:\mathcal{H}_{D}\to\mathcal{H}_{D}$ does not, a priori , map $\mathcal{H}_{0}$ into itself, we first need to adapt the setting in order to apply the spectral theorem of [˜2](#ThmTh2).

Because $\mathcal{H}_{0}\subset\mathcal{H}_{D}$ isometrically, we have for all $u,v\in\mathcal{H}_{0}$:

$$ $\displaystyle\begin{split}\left\langle D[u]\,,\,D[v]\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}&=\left\langle D\big[\Pi_{\mathcal{H}_{0}}u\big]\,,\,D\big[\Pi_{\mathcal{H}_{0}}v\big]\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}\\ &=\left\langle\Pi_{\mathcal{H}_{0}}u\,,\,D^{*}D\big[\Pi_{\mathcal{H}_{0}}v\big]\right\rangle_{\mathcal{H}_{D}}\\ &=\left\langle u\,,\,\big(\Pi_{\mathcal{H}_{0}}D^{*}D\Pi_{\mathcal{H}_{0}}\big)[v]\right\rangle_{\mathcal{H}_{D}},\end{split}$ (78) $$

where we used in the last step that $\Pi_{\mathcal{H}_{0}}$ is self-adjoint.

We can therefore apply the spectral theorem [˜2](#ThmTh2) to the bounded self-adjoint operator $\Pi_{\mathcal{H}_{0}}D^{*}D\Pi_{\mathcal{H}_{0}}:\mathcal{H}_{0}\to\mathcal{H}_{0}$, obtaining the analogue of the decomposition in [Equation˜58](#A7.E58):

$$ $\displaystyle\Pi_{\mathcal{H}_{0}}D^{*}D\Pi_{\mathcal{H}_{0}}$ $\displaystyle=\int_{\mathbb{R}_{+}}\lambda\,\pi^{\mathcal{H}_{0}}_{D}(\textrm{d}\lambda),$ $\displaystyle I_{\mathcal{H}_{0}}$ $\displaystyle=\int_{\mathbb{R}_{+}}\pi^{\mathcal{H}_{0}}_{D}(\textrm{d}\lambda).$ (79) $$

###### Lemma 1 .

###### Proof.

###### Lemma 2 .

###### Proof.

##### Regularized spaces.

Fixing $\alpha>0$, and analogously to [Equations˜59](#A7.E59) and [60](#A7.E60), we define the regularized projection and subspace:

$$ $\displaystyle\Pi_{D,\mathcal{H}_{0}}^{\alpha}$ $\displaystyle:=\int_{\alpha^{2}}^{+\infty}\pi^{\mathcal{H}_{0}}_{D}(\textrm{d}\lambda),$ $\displaystyle\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$ $\displaystyle:=\operatorname{Im}\Pi_{D,\mathcal{H}_{0}}^{\alpha}\subset\mathcal{H}_{0}\subset\mathcal{H}_{D}.$ (80) $$

Let $k:\Omega\times\Omega\to\mathbb{R}$ be the reproducing kernel of $\mathcal{H}_{0}$. Then, by Paulsen & Raghupathi (2016, Theorem 2.5), $\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$ is an RKHS with reproducing kernel

$$ $k_{\alpha}(x,y):=\Pi_{D,\mathcal{H}_{0}}^{\alpha}[k(x,\cdot)](y),\qquad\forall x,y\in\Omega.$ (81) $$

##### Norm equivalence.

Since $\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}\subset\mathcal{H}_{0}\subset\mathcal{H}_{D}$, inequality in [Equation˜30](#A7.E30) remains valid, i.e. for all $u\in\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$:

$$ $\left\|D[u\right\|]_{\textrm{L}^{2}(\Omega,\mu)}\leq\left\|u\right\|_{\mathcal{H}_{D}}.$ (82) $$

Furthermore, by an argument entirely analogous to [Equation˜61](#A7.E61), we also have

$$ $\left\|D[u\right\|]_{\textrm{L}^{2}(\Omega,\mu)}\geq\alpha\left\|u\right\|_{\mathcal{H}_{D}},\qquad\forall u\in\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}.$ (83) $$

In particular, the functional

$$ $\left\|\cdot\right\|_{D}:\left\{\begin{array}[]{lll}\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}&\to&\mathbb{R}\\ u&\mapsto&\left\|D[u\right\|]_{\textrm{L}^{2}(\Omega,\mu)}\end{array}\right.$ (84) $$

defines a norm equivalent to $\left\|\cdot\right\|_{\mathcal{H}_{D}}$ on $\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$.
By [˜1](#ThmLM1), the pair $\big(\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha},\left\|\cdot\right\|_{D}\big)$ is itself an RKHS with a reproducing kernel $k_{D}$.

##### Isometry property.

The crucial observation is that $D$ is an isometry with respect to this norm. Indeed, for all $u,v\in\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}$,

$$ $\left\langle u\,,\,v\right\rangle_{D}=\left\langle D[u]\,,\,D[v]\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}.$ (85) $$

This allows us to characterize the associated Green’s function.

See [1](#ThmTh1)

###### Proof.

For all $f\in\textrm{L}^{2}(\Omega,\mu)$ and $x\in\Omega$,

$$ $\displaystyle\begin{split}\int_{\Omega}g_{D}(x,s)f(s)\mu(\textrm{d}s)&=\left\langle g_{D}(x,\cdot)\,,\,f\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}\\ &=\left\langle D[k_{D}(x,\cdot)]\,,\,f\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}\\ &=\left\langle k_{D}(x,\cdot)\,,\,D^{*}f\right\rangle_{D}\\ &=\big(D^{*}f\big)(x).\end{split}$ (86) $$

Since $D$ is an isometry, [˜2](#ThmLM2) gives $DD^{*}=\Pi_{D[\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}]}$. Therefore,

$$ $D\Big[x\mapsto\int_{\Omega}g_{D}(x,s)f(s)\mu(\textrm{d}s)\Big]=D\big[D^{*}f\big]=\Pi_{D[\mathcal{H}_{D,\mathcal{H}_{0}}^{\alpha}]}f,$ (87) $$

which precisely shows that $g_{D}$ is a generalized Green’s function.
∎

The key insight of [˜1](#ThmTh1) is that, in the PINNs setting—and most notably in our algorithm—we implicitly construct the reproducing kernel $k_{D}$ associated with the norm $\left\|\cdot\right\|_{D}$ on the regularized tangent space $T^{\alpha}_{\theta}\mathcal{M}$ of the neural network manifold $\mathcal{M}$, at cutoff level $\alpha$. This kernel is precisely the PINNs NNTK introduced by Schwencke & Furtlehner (2025).

A crucial consequence is that the regularization of the Gram matrix is not merely a “numerical trick” to guarantee stability: it is the very mechanism that ensures the Green’s function is well defined.

###### Proof.

##### Conceptual interpretation.

This perspective also offers a profound interpretation of the procedure: rather than attempting to invert the operator $D$ directly, we build a kernel $k_{D}$ whose associated metric makes $D$ an isometry, and thus ensures that $D^{*}$ acts as the generalized left-inverse of $D$.
The magic of the kernel lies in the following facts:

- (i)
We never need to compute $D^{*}$ explicitly, since it is implicitly encoded in the relation
$\left\langle D[k_{D}(x,\cdot)]\,,\,f\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}=\left\langle k_{D}(x,\cdot)\,,\,D^{*}f\right\rangle_{D}.$
(88)
- (ii)
The same formula allows us to directly evaluate the generalized solution $D^{*}f$: indeed, for all $x\in\Omega$, the reproducing property gives
$D^{*}f(x)=\left\langle k_{D}(x,\cdot)\,,\,D^{*}f\right\rangle_{D}.$
(89)

##### Comparison with Ridge regression.

An analogous analysis holds for Ridge regression. However, instead of inverting $D$ “via isometry,” we invert the augmented operator $\big(D,\alpha\textrm{I}_{\mathcal{H}_{D}}\big)$.

## Appendix H Proofs

### H.1 Statement and proof of ˜ 3

We start by recalling the following statements from Schwencke & Furtlehner (2025).

###### Definition ( Schwencke & Furtlehner, 2025 , Definition 4 ) .

A linear operator $A:\mathcal{H}\to\mathcal{H}$ is an integral operator given that there is $k:\Omega\times\Omega\to\mathbb{K}$, $\mathbb{K}\in\{\mathbb{R},\mathbb{C}\}$, such that: $\textrm{for all }f\in\mathcal{H}$, $\textrm{for all }x\in\Omega$

$$ $A\big(f\big)(x)=\left\langle k(x,\cdot)\,,\,f\right\rangle_{\mathcal{H}}.$ (90) $$

###### Lemma ( Schwencke & Furtlehner, 2025 , Lemma 1 ) .

Let us be $\mathcal{H}_{0}:=\operatorname*{Span}(u_{p}:1\leq p\leq P)\subset\mathcal{H}$ and consider the Gram matrix $G_{pq}:=\left\langle u_{p}\,,\,u_{q}\right\rangle_{\mathcal{H}}$ of $(u_{p})$ and its eigen-decomposition $G=U\Delta^{2}U^{t}$. Then:

$$ $L_{p}:=\sum_{1\leq q\leq P}u_{q}U_{q,p}\Delta^{\dagger}_{p},$ (91) $$

is an orthonormal basis of $\mathcal{H}_{0}$.
In particular, $\Pi_{\mathcal{H}_{0}}$ is an integral operator whose kernel is:

$$ $k(x,y)=\sum_{1\leq p,q\leq P}u_{p}(x)G^{\dagger}_{p,q}u_{q}(y).$ (92) $$

Furthermore $L_{p}$ are the left-singular vector of the so-called synthesis operator(^3^33Name and notation are taken from Adcock & Huybrechs (2019).):

$$ $\mathcal{T}:\left\{\begin{array}[]{lll}\mathbb{R}^{P}&\to&\mathcal{H}_{0}\\ \alpha&\mapsto&\sum\limits_{1\leq p\leq P}\alpha_{p}u_{p}\end{array}\right..$ (93) $$

###### Proposition 3 .

Given the scalar loss

$$ $\ell({\bm{\theta}}):=\mathcal{L}\!\left(u_{{\bm{\theta}}}\right)\stackrel{{\scriptstyle(\ref{eqn:functional-loss})}}{{=}}\tfrac{1}{2}\,\left\|u_{{\bm{\theta}}}-f\right\|_{\textrm{L}^{2}(\Omega,\mu)}^{2},$ (94) $$

the Natural Gradient update of [Equation˜12](#S2.E12)

$$ $\displaystyle u_{{\bm{\theta}}_{t+1}}$ $\displaystyle\leftarrow u_{{\bm{\theta}}_{t}}-\eta\,\Pi_{T_{{\bm{\theta}}_{t}}\mathcal{M}}\!\left(\nabla\mathcal{L}_{u_{{\bm{\theta}}_{t}}}\right)\,;$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,\mathrm{d}u_{{\bm{\theta}}_{t}}^{\dagger}\!\left(\Pi_{T_{{\bm{\theta}}_{t}}\mathcal{M}}\!\left(\nabla\mathcal{L}_{u_{{\bm{\theta}}_{t}}}\right)\right)$ ( 12 ) $$

can be equivalently written as

$$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,G_{{\bm{\theta}}_{t}}^{\dagger}\nabla\ell({\bm{\theta}}_{t})\,;$ $\displaystyle{G_{{\bm{\theta}}_{t}}}_{p,q}$ $\displaystyle:=\left\langle\partial_{p}u_{{\bm{\theta}}_{t}}\,,\,\partial_{q}u_{{\bm{\theta}}_{t}}\right\rangle_{\textrm{L}^{2}(\Omega,\mu)}.$ ( 13 ) $$

###### Proof.

Since the tangent space $T_{{\bm{\theta}}}\mathcal{M}$ of [Equation˜11](#S2.E11):

$$ $T_{{\bm{\theta}}}\mathcal{M}:=\mathrm{Im}(\mathrm{d}u_{{\bm{\theta}}})=\operatorname*{Span}\left(\partial_{p}u_{{\bm{\theta}}}:1\leq p\leq P\right)\subset\mathcal{H},$ ( 11 ) $$

is finite-dimensional, we may invoke Schwencke & Furtlehner (2025, Lemma 1). This result shows that the Natural Neural Tangent Kernel (NNTK), given by

$$ $\displaystyle{NNTK}_{\bm{\theta}}(x,y)$ $\displaystyle:=\sum_{1\leq p,q\leq P}\big(\partial_{p}u_{{\bm{\theta}}}(x)\big)\,{G_{\bm{\theta}}^{\dagger}}_{pq}\,\big(\partial_{q}u_{{\bm{\theta}}}(y)\big)^{t},$ $\displaystyle{G_{{\bm{\theta}}}}_{p,q}$ $\displaystyle:=\left\langle\partial_{p}u_{{\bm{\theta}}}\,,\,\partial_{q}u_{{\bm{\theta}}}\right\rangle_{\mathcal{H}},$ (95) $$

is the kernel of the orthogonal projection $\Pi^{\bot}_{T_{{\bm{\theta}}}\mathcal{M}}$ onto $T_{{\bm{\theta}}}\mathcal{M}$. Therefore, by [Equation˜90](#A8.E90), for all $x\in\Omega$,

$$ $\displaystyle\begin{split}\Pi^{\bot}_{T_{{\bm{\theta}}}\mathcal{M}}\!\left(\nabla\mathcal{L}_{|u_{{\bm{\theta}}}}\right)(x)&=\left\langle NNTK_{\bm{\theta}}(x,\cdot)\,,\,\nabla\mathcal{L}_{|u_{{\bm{\theta}}}}\right\rangle_{\mathcal{H}}\\ &\stackrel{{\scriptstyle(\ref{eqn:NNTK-def-gram})}}{{=}}\sum_{1\leq p,q\leq P}\partial_{p}u_{{\bm{\theta}}}(x)\,{G_{\bm{\theta}}^{\dagger}}_{pq}\,\left\langle\partial_{q}u_{{\bm{\theta}}}\,,\,\nabla\mathcal{L}_{|u_{{\bm{\theta}}}}\right\rangle_{\mathcal{H}}.\end{split}$ (96) $$

Next, note that

$$ $\left\langle\partial_{q}u_{{\bm{\theta}}}\,,\,\nabla\mathcal{L}_{|u_{{\bm{\theta}}}}\right\rangle_{\mathcal{H}}=\textrm{d}\mathcal{L}_{|u_{{\bm{\theta}}}}\!\big(\partial_{q}u_{{\bm{\theta}}}\big)\stackrel{{\scriptstyle\text{chain rule}}}{{=}}\partial_{q}\mathcal{L}(u_{{\bm{\theta}}})\stackrel{{\scriptstyle(\ref{eqn:scalar-loss-def})}}{{=}}\partial_{q}\ell({\bm{\theta}}).$ (97) $$

Therefore, by linearity of $\textrm{d}u_{{\bm{\theta}}}^{\dagger}$,

$$ $\textrm{d}u_{{\bm{\theta}}}^{\dagger}\!\left(\Pi^{\bot}_{T_{{\bm{\theta}}}\mathcal{M}}\!\left(\nabla\mathcal{L}_{|u_{{\bm{\theta}}}}\right)\right)\stackrel{{\scriptstyle(\ref{eqn:projection-formula}),(\ref{eqn:partial-ell-formula})}}{{=}}\sum_{1\leq p,q\leq P}\textrm{d}u_{{\bm{\theta}}}^{\dagger}\!\big(\partial_{p}u_{{\bm{\theta}}}\big)\,{G_{\bm{\theta}}^{\dagger}}_{pq}\,\partial_{q}\ell({\bm{\theta}}).$ (98) $$

Finally, observe that $\partial_{p}u_{{\bm{\theta}}}=\textrm{d}u_{{\bm{\theta}}}({\bm{e}}^{(p)})$, where ${\bm{e}}^{(p)}$ is the $p$-th canonical basis vector of $\mathbb{R}^{P}$.
If $\textrm{d}u_{{\bm{\theta}}}$ were invertible, we would directly obtain

$$ $\textrm{d}u_{{\bm{\theta}}}^{\dagger}\!\big(\partial_{p}u_{{\bm{\theta}}}\big)={\bm{e}}^{(p)},$ (99) $$

which would complete the argument. However, this invertibility does not hold in general.

To address this, recall that $\textrm{d}u_{{\bm{\theta}}}$ can be identified with the synthesis operator $\mathcal{T}$ introduced in [Equation˜93](#A8.E93) of Schwencke & Furtlehner (2025, Lemma 1). From the final part of that lemma, we know that $\operatorname{Im}\textrm{d}u_{{\bm{\theta}}}^{\dagger}=\operatorname{Im}G_{\bm{\theta}}^{\dagger}$. Consequently,

$$ $G_{\bm{\theta}}^{\dagger}{\bm{e}}^{(p)}=G_{\bm{\theta}}^{\dagger}\textrm{d}u_{{\bm{\theta}}}^{\dagger}\!\big(\partial_{p}u_{{\bm{\theta}}}\big).$ (100) $$

Putting all pieces together yields the desired update rule, thereby completing the proof.
∎

###### Definition ( Schwencke & Furtlehner, 2025 , Definition 4 ) .

###### Lemma ( Schwencke & Furtlehner, 2025 , Lemma 1 ) .

###### Proposition 3 .

###### Proof.

### H.2 Ridge-regression implementation ANaGRAM

In the following, we show that a Ridge-regression can be implemented in ANaGRAM’s update rule given by [Equation˜14](#S2.E14).

###### Proposition 4 .

A Ridge-regression can be implemented in the SVD-based update [Equation˜14](#S2.E14) by replacing the pseudo-inverse $\widehat{\Delta}^{\dagger}$ with

$$ $\left(\frac{\widehat{\Delta}_{t,i}}{\widehat{\Delta}_{t,i}^{2}+S\alpha}\right)_{1\leq i\leq\textrm{r}_{\text{svd}}}.$ (101) $$

###### Proof.

As shown in (Schwencke & Furtlehner, 2025, Section E), the ANaGRAM’s update of [Equation˜14](#S2.E14):

$$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,\widehat{\phi}^{\dagger}\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}_{t}};$ $\displaystyle\widehat{\phi}_{i,p}$ $\displaystyle:=\partial_{p}u_{{\bm{\theta}}}(x_{i});$ $\displaystyle\left(\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}}\right)_{i}$ $\displaystyle:=u_{{\bm{\theta}}}(x_{i})-f(x_{i}),$ ( 14 ) $$

is equivalent to the update with the empirical matrix $\widehat{\mathcal{G}}_{{\bm{\theta}}}$:

$$ $\displaystyle{\bm{\theta}}_{t+1}$ $\displaystyle\leftarrow{\bm{\theta}}_{t}-\eta\,\widehat{\mathcal{G}}_{{\bm{\theta}}_{t}}^{\dagger}\nabla\ell({\bm{\theta}}_{t})\,;$ $\displaystyle\widehat{\mathcal{G}}_{{\bm{\theta}}_{t}}$ $\displaystyle:=\frac{1}{S}\widehat{\phi}_{{\bm{\theta}}_{t}}^{t}\widehat{\phi}_{{\bm{\theta}}_{t}},$ (102) $$

where $\ell$ is defined in [Equation˜9](#S2.E9):

$$ $\ell({\bm{\theta}}):=\frac{1}{2S}\sum_{i=1}^{S}\left(u_{{\bm{\theta}}}(x_{i})-f(x_{i})\right)^{2}.$ ( 9 ) $$

Thus, we get immediately

$$ $\nabla\ell({\bm{\theta}}_{t})=\frac{1}{S}\widehat{\phi}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}=\frac{1}{S}\widehat{U}\widehat{\Delta}\widehat{V}^{t}_{\bm{\theta}}\widehat{\nabla\mathcal{L}}_{\bm{\theta}},$ (103) $$

where we used the SVD decomposition of $\widehat{\phi}$:

$$ . (104) $$

Using [Equation˜104](#A8.E104) again, we have

$$ $\widehat{\mathcal{G}}_{{\bm{\theta}}}=\frac{1}{S}\widehat{U}\widehat{\Delta}^{2}_{\bm{\theta}}\widehat{U}^{t}_{\bm{\theta}},$ (105) $$

thus for all $\alpha>0$

$$ $\widehat{\mathcal{G}}_{{\bm{\theta}}}+\alpha I_{d}=\frac{1}{S}\widehat{U}\widehat{\Delta}^{2}_{\bm{\theta}}\widehat{U}^{t}_{\bm{\theta}}+\alpha\widehat{U}\widehat{U}^{t}=\widehat{U}\left(\operatorname{diag}\left(\frac{\widehat{\Delta}^{2}_{{\bm{\theta}}_{i}}}{S}+\alpha\right)_{1\leq i\leq\textrm{r}_{\text{svd}}}\right)\widehat{U}^{t}_{\bm{\theta}},$ (106) $$

which implies

$$ $\left(\widehat{\mathcal{G}}_{{\bm{\theta}}}+\alpha I_{d}\right)^{-1}=\widehat{U}\left(\operatorname{diag}\left(\frac{S}{\widehat{\Delta}^{2}_{{\bm{\theta}}_{i}}+S\alpha}\right)_{1\leq i\leq\textrm{r}_{\text{svd}}}\right)\widehat{U}^{t}_{\bm{\theta}}.$ (107) $$

This finally yields

$$ $\displaystyle\left(\widehat{\mathcal{G}}_{{\bm{\theta}}}+\alpha I_{d}\right)^{-1}\nabla\ell({\bm{\theta}}_{t})$ $\displaystyle\stackrel{{\scriptstyle(\ref{eqn:nable-ell-expression})}}{{=}}\widehat{U}\left(\operatorname{diag}\left(\frac{S}{\widehat{\Delta}^{2}_{{\bm{\theta}}_{i}}+S\alpha}\right)_{1\leq i\leq\textrm{r}_{\text{svd}}}\right)\widehat{U}^{t}_{\bm{\theta}}\frac{1}{S}\widehat{U}\widehat{\Delta}\widehat{V}^{t}_{\bm{\theta}}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}$ (108) $\displaystyle=\widehat{U}\left(\operatorname{diag}\left(\frac{\widehat{\Delta}_{t,i}}{\widehat{\Delta}^{2}_{{\bm{\theta}}_{i}}+S\alpha}\right)_{1\leq i\leq\textrm{r}_{\text{svd}}}\right)\widehat{V}^{t}_{\bm{\theta}}\widehat{\nabla\mathcal{L}}_{\bm{\theta}},$ (109) $$

which conludes the proof.
∎

###### Proposition 4 .

###### Proof.

### H.3 Proof of ˜ 1

To prove [˜1](#ThmProp1), we need the following lemma:

###### Lemma 3 .

For $1\leq M\leq N\leq\textrm{r}_{\text{svd}}$:

$$ $\left(\textrm{RCE}^{S}_{M}\right)^{2}-\left(\textrm{RCE}^{S}_{N}\right)^{2}=\frac{1}{S}\left\|\Pi^{M}_{N}\widehat{V}^{T}\widehat{\nabla\mathcal{L}}\right\|^{2}_{\mathbb{R}^{S}}.$ (110) $$

###### Proof.

Let us first recall the definition of the $\textrm{RCE}_{N}^{S}$ in [Equation˜17](#S3.E17), namely

$$ $\textrm{RCE}^{S}_{N}:=\frac{1}{\sqrt{S}}\left\|\widehat{V}\Pi^{0}_{N}\widehat{V}^{T}\widehat{\nabla\mathcal{L}}-\widehat{\nabla\mathcal{L}}\right\|_{\mathbb{R}^{S}}.$ ( 17 ) $$

Fixing $1\leq N\leq M\leq\textrm{r}_{\text{svd}}$ and
applying the same reasoning as in [Equation˜120](#A8.E120) to $\textrm{RCE}_{M}^{S}$ and $\textrm{RCE}_{N}^{S}$ (see the proof of [˜1](#ThmProp1) in [Section˜H.3](#A8.SS3)), we get

$$ $\displaystyle S\left(\textrm{RCE}^{S}_{M}\right)^{2}$ $\displaystyle=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\sum_{p=1}^{M}\left(\widehat{V}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2};$ $\displaystyle S\left(\textrm{RCE}^{S}_{N}\right)^{2}$ $\displaystyle=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\sum_{p=1}^{N}\left(\widehat{V}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2},$ (111) $$

and therefore

$$ $\displaystyle\begin{split}S\left(\left(\textrm{RCE}^{S}_{N}\right)^{2}-\left(\textrm{RCE}^{S}_{M}\right)^{2}\right)&=\sum_{p=1}^{N}\left(\widehat{V}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2}-\sum_{p=1}^{M}\left(\widehat{V}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2}\\ &\stackrel{{\scriptstyle M\leq N}}{{=}}\sum_{p=M+1}^{N}\left(\widehat{V}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2}=\sum_{p=M+1}^{N}\left({\bm{e}}^{(p)^{t}}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2}\\ &=\sum_{p=M+1}^{N}\left(\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{V}{\bm{e}}^{(p)}\right)\left({\bm{e}}^{(p)^{t}}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)\\ &=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{V}\underbracket{\left(\sum_{p=M+1}^{N}{\bm{e}}^{(p)}{\bm{e}}^{(p)^{t}}\right)}_{=\Pi_{N}^{M}~by~\lx@cref{creftype~refnum}{eqn:PiMN-def}}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\\ &=\left\langle\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\,,\,\Pi_{N}^{M}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right\rangle_{\mathbb{R}^{S}}\\ &\stackrel{{\scriptstyle\begin{subarray}{c}\Pi_{N}^{M^{2}}=\Pi_{N}^{M}\\ \Pi_{N}^{M^{t}}=\Pi_{N}^{M}\end{subarray}}}{{=}}\left\langle\Pi_{N}^{M}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\,,\,\Pi_{N}^{M}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right\rangle_{\mathbb{R}^{S}}\\ &=\left\|\Pi_{N}^{M}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right\|_{\mathbb{R}^{S}}^{2},\end{split}$ (112) $$

where we use in the penultimate equality, the fact that $\Pi^{M}_{N}$ is an orthogonal projection. ∎

###### Remark 5 .

The above lemma provides an interesting property that gives a further understanding of what is happening during the flattening, i.e. $\textrm{RCE}^{S}_{M}-\textrm{RCE}^{S}_{N}\approx 0$. In particular, as $\left(\textrm{RCE}^{S}_{M}\right)^{2}-\left(\textrm{RCE}^{S}_{N}\right)^{2}=\left(\textrm{RCE}^{S}_{M}-\textrm{RCE}^{S}_{N}\right)\left(\textrm{RCE}^{S}_{M}+\textrm{RCE}^{S}_{N}\right)$, therefore flattening for the components in the range $[N_{\text{flat}},\textrm{r}_{\text{cutoff}}]$ means that $\frac{1}{S}\left\|\Pi^{M}_{N}\widehat{V}^{T}\widehat{\nabla\mathcal{L}}\right\|^{2}_{\mathbb{R}^{S}}\approx 0$. In other words, the problem is "learned" for those components, as the projection of the functional gradient (which is propotional to the error) on their corresponding span is null. The proof of this lemma is provided in [Section˜H.3](#A8.SS3).

See [1](#ThmProp1)

###### Proof.

The first statement is a direct consequence of [˜3](#ThmLM3) proven above.

Let us now show that the second statement takes place. Since $\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\in\textrm{L}^{2}(\Omega,\mu)$ and $\operatorname{Im}\textrm{d}u_{{\bm{\theta}}}\subset\textrm{L}^{2}(\Omega,\mu)$, the law of large numbers yields that for all $1\leq p,q\leq P$

$$ $\displaystyle\lim\limits_{S\to\infty}\frac{1}{S}\sum_{i=1}^{S}\left[\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x_{i})\right]^{2}$ $\displaystyle=\lim\limits_{S\to\infty}\frac{1}{S}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}=\int_{\Omega}\left[\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\right]^{2}\mu(\textrm{d}x)\quad\textit{a.s},$ (113) $\displaystyle\lim\limits_{S\to\infty}\frac{1}{S}\sum_{i=1}^{S}\partial_{p}u_{{\bm{\theta}}}(x_{i})\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x_{i})$ $\displaystyle=\lim\limits_{S\to\infty}\frac{1}{S}\widehat{\phi}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}=\int_{\Omega}\partial_{p}u_{{\bm{\theta}}}(x)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x)\quad\textit{a.s},$ (114) $\displaystyle\lim\limits_{S\to\infty}\frac{1}{S}\sum_{i=1}^{S}\partial_{p}u_{{\bm{\theta}}}(x_{i})\partial_{q}u_{{\bm{\theta}}}(x_{i})$ $\displaystyle=\lim\limits_{S\to\infty}\frac{1}{S}\widehat{\phi}_{{\bm{\theta}}_{p}}^{t}\widehat{\phi}_{{\bm{\theta}}_{q}}=\int_{\Omega}\partial_{p}u_{{\bm{\theta}}}(x)\partial_{q}u_{{\bm{\theta}}}(x)\mu(\textrm{d}x)\quad\textit{a.s}.$ (115) $$

In particular, this implies

$$ $\displaystyle\lim\limits_{S\to\infty}\frac{1}{S}\widehat{\phi}^{t}\widehat{\phi}=G_{\bm{\theta}}=U_{\bm{\theta}}\Delta_{\bm{\theta}}^{2}U_{\bm{\theta}}^{t}\quad\textit{a.s}.$ (116) $$

Since the eigenvectors (and eigenvalues) are continuous functions of the matrix coefficients (by polynomial dependence through the characteristic polynomial) and taking into account that $\frac{1}{S}\widehat{\phi}^{t}\widehat{\phi}=\frac{1}{S}\widehat{U}\Delta_{\bm{\theta}}^{2}\widehat{U}^{t}$, this yields

$$ $\displaystyle\lim\limits_{S\to\infty}\widehat{U}$ $\displaystyle=U_{\bm{\theta}}\quad\textit{a.s};$ $\displaystyle\lim\limits_{S\to\infty}\frac{1}{S}\widehat{\Delta}^{2}$ $\displaystyle=\Delta_{\bm{\theta}}^{2}\quad\textit{a.s}.$ (117) $$

By continuity of the square root and the inverse on $\mathbb{R}^{*}_{+}$, we get that for all $1\leq p\leq P$ such that $\Delta_{{\bm{\theta}}_{p}}>0$

$$ $\displaystyle\lim\limits_{S\to\infty}\sqrt{S}\widehat{\Delta}_{{\bm{\theta}}_{p}}^{-1}=\Delta_{{\bm{\theta}}_{p}}^{-1}\quad\textit{a.s},$ (118) $$

and thus for all $1\leq p\leq P$ such that $\Delta_{{\bm{\theta}}_{p}}>0$, we have almost surely

$$ $\displaystyle\begin{split}\lim\limits_{S\to\infty}\frac{1}{\sqrt{S}}\widehat{V}_{{\bm{\theta}}_{p}}^{T}\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}}&=\lim\limits_{S\to\infty}\sqrt{S}\widehat{\Delta}_{{\bm{\theta}}_{p}}^{-1}\widehat{U}_{{\bm{\theta}}_{p}}^{t}\left(\sum_{q=1}^{P}{\bm{e}}^{(q)}{\bm{e}}^{(q)^{t}}\right)\frac{1}{S}\widehat{\phi}^{t}\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}}\\ &=\sum_{q=1}^{P}\left(\lim\limits_{S\to\infty}\sqrt{S}\widehat{\Delta}_{{\bm{\theta}}_{p}}^{-1}\right)\left(\lim\limits_{S\to\infty}\widehat{U}_{{\bm{\theta}}_{p}}^{t}{\bm{e}}^{(q)}\right)\left(\lim\limits_{S\to\infty}\frac{1}{S}\widehat{\phi}_{{\bm{\theta}}_{q}}^{t}\widehat{\nabla\mathcal{L}}_{{\bm{\theta}}}\right)\\ &=\sum_{q=1}^{P}\Delta_{{\bm{\theta}}_{p}}^{-1}U_{{\bm{\theta}}_{p}}^{t}{\bm{e}}^{(q)}\int_{\Omega}\partial_{q}u_{{\bm{\theta}}}(x)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x)\\ &=\int_{\Omega}\textrm{d}u_{{\bm{\theta}}}\left(U_{{\bm{\theta}}_{p}}\Delta_{{\bm{\theta}}_{p}}^{-1}\right)(x)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x)\\ &=\int_{\Omega}V_{{\bm{\theta}}_{p}}(x)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x),\end{split}$ (119) $$

where we used in the last equality,
the identification of the singular vectors of $\textrm{d}u_{{\bm{\theta}}}$
in (Schwencke & Furtlehner, 2025, Lemma 1 p. 24, section C.2).
Returning to the definition of the $\textrm{RCE}_{N}^{S}$ in [Equation˜17](#S3.E17), namely

$$ $\textrm{RCE}^{S}_{N}:=\frac{1}{\sqrt{S}}\left\|\widehat{V}\Pi^{0}_{N}\widehat{V}^{T}\widehat{\nabla\mathcal{L}}-\widehat{\nabla\mathcal{L}}\right\|_{\mathbb{R}^{S}},$ ( 17 ) $$

we get

$$ $\displaystyle\begin{split}S\left(\textrm{RCE}^{S}_{N}\right)^{2}&=\left\langle\widehat{V}\Pi^{0}_{N}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\,,\,\widehat{V}\Pi^{0}_{N}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right\rangle_{\mathbb{R}^{S}}\\ &=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}+\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{V}\overbracket{\Pi^{0}_{N}\underbracket{\widehat{V}^{t}\widehat{V}}_{=I_{d}}\Pi^{0}_{N}}^{=\Pi^{0}_{N}}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-2\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{V}\Pi^{0}_{N}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\\ &=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{V}\Pi^{0}_{N}\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\\ &=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{V}\left(\sum\limits_{p=1}^{N}{\bm{e}}^{(p)}{\bm{e}}^{(p)^{t}}\right)\widehat{V}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\\ &=\widehat{\nabla\mathcal{L}}_{\bm{\theta}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}-\sum_{p=1}^{N}\left(\widehat{V}_{{\bm{\theta}}_{p}}^{t}\widehat{\nabla\mathcal{L}}_{\bm{\theta}}\right)^{2},\end{split}$ (120) $$

where in the second equality, we use the fact that $\widehat{V}$ is orthogonal and $\Pi_{N}^{0}$ is a projection. Combining [Equations˜113](#A8.E113) and [119](#A8.E119), this yields

$$ $\displaystyle\lim\limits_{S\to\infty}\left(\textrm{RCE}^{S}_{N}\right)^{2}$ $\displaystyle=\int_{\Omega}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)^{2}\mu(\textrm{d}x)-\sum_{p=1}^{N}\left(\int_{\Omega}V_{{\bm{\theta}}_{p}}(x)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x)\right)^{2}\quad\textit{a.s}.$ (121) $$

By Fubini’s theorem, we have almost surely

$$ $\displaystyle\sum_{p=1}^{N}\left(\int_{\Omega}V_{{\bm{\theta}}_{p}}(x)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x)\right)^{2}$ $\displaystyle=\int_{\Omega^{2}}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\left(\sum_{p=1}^{N}V_{{\bm{\theta}}_{p}}(x)V_{{\bm{\theta}}_{p}}(y)\right)\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(y)\mu(\textrm{d}x)\mu(\textrm{d}y)$ $\displaystyle=\int_{\Omega}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}(x)\mu(\textrm{d}x)$ (122) $\displaystyle=\left\|\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)},$ $$

where in the second equality, we used (Schwencke & Furtlehner, 2025, Theorem 4 p. 23, section C.2) and the fact that $\left(\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\right)^{2}=\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}$ in the third. Therefore,
from [Equation˜121](#A8.E121) and [Equation˜122](#A8.E122)

$$ $\displaystyle\lim\limits_{S\to\infty}\left(\textrm{RCE}^{S}_{N}\right)^{2}$ $\displaystyle=\left\|\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}-\left\|\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}\quad\textit{a.s},$ (123) $\displaystyle=\left\|\nabla\mathcal{L}_{u_{{\bm{\theta}}}}-\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}\quad\textit{a.s},$ (124) $$

where in the second equality, we use in the reverse order a reasoning similar to [Equation˜120](#A8.E120).
Finally, we obtain

$$ $\displaystyle\left\|\nabla\mathcal{L}_{u_{{\bm{\theta}}}}-\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)}=\left\|\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)^{\bot}}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|^{2}_{\textrm{L}^{2}(\Omega,\mu)},$ (125) $$

which comes from the canonical decomposition in Hilbert spaces, i.e. using that $\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)$ is a closed subspace and

$$ $\nabla\mathcal{L}_{u_{{\bm{\theta}}}}=\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}+\Pi^{\bot}_{\operatorname*{Span}(V_{{\bm{\theta}}_{i}}\,:\,1\leq i\leq N)^{\bot}}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}.$ (126) $$

This completes the proof.
∎

###### Corollary 1 .

For $1\leq M\leq N\leq\textrm{r}_{\text{svd}}$:

$$ $\lim\limits_{S\to\infty}\left(\textrm{RCE}^{S}_{M}\right)^{2}-\left(\textrm{RCE}^{S}_{N}\right)^{2}=\left\|\Pi^{\bot}_{T^{M}_{N}\mathcal{M}}\nabla\mathcal{L}_{u_{{\bm{\theta}}}}\right\|_{\textrm{L}^{2}(\Omega)}^{2}$ (127) $$

###### Proof.

Apply [˜1](#ThmProp1) to [Equation˜110](#A8.E110) of [˜3](#ThmLM3).
∎

###### Lemma 3 .

###### Proof.

###### Remark 5 .

###### Proof.

###### Corollary 1 .

###### Proof.

## Peer Reviews

confidence: 5
contribution: 4
initial_rating: 8
presentation: 4
questions: No questions, great paper.
soundness: 4
strengths: See summary. The paper provides a notable and timely contribution for what has been an outstanding issue with PINNs. The technique is well-rooted in analysis and provides a principled solution to the problem.
summary: The authors analyze and propose an extension of ANaGRAM, a natural gradient method technique for solving PDEs in the physics-informed neural network setting. I have extensively studied the original ANaGRAM work, and the current work represents a non-trivial extension. Notably, the authors are able to obtain convergence "to machine-precision". This is a timely paper, as many are attempting to achieve machine-precision accuracy from PINNs. While this is a trivial achievement for conventional finite element solutions to PDEs, PINNs have struggled to achieve this and this is viewed by many as an important stepping stone to achieving AI/ML models with comparable quality of science/engineering-relevance. 

For example, deepmind just launched an effort to develop a numerical counterexample for the millinium prize of finite time blowup for navier-stokes. In contrast, my impression is that this has the desirable feature of being an "optimization only" solution that achieves desired accuracy with no ad hoc modification of the loss or architecture. 

I have a few suggestions about how the text could be further improved, but overall this is a strong paper and should be accepted.
weaknesses: As noted, the paper could benefit from some minor modifications.

(1) The work should be contextualized in the context of other methods. The only one to my knowledge that has been able to achieve this accuracy is from Ching-Yao Lai's group at Stanford (Wang, Yongji, and Ching-Yao Lai. "Multi-stage neural networks: Function approximator of machine precision." Journal of Computational Physics 504 (2024): 112865.), which is the technique deepmind used for their work (Wang, Yongji, et al. "Discovery of Unstable Singularities." arXiv preprint arXiv:2509.14185 (2025). These are very recent works but would help frame the importance of the fundamental issue and highlights the way in which PINNs can answer ill-posed questions that are challenging for conventional FEM. The authors may also be interested in (Trask, Nathaniel, Amelia Henriksen, Carianne Martinez, and Eric Cyr. "Hierarchical partition of unity networks: fast multilevel training." In Mathematical and Scientific Machine Learning, pp. 271-286. PMLR, 2022.), where authors achieved machine precision but only for regression problems. All of these works achieve accuracy through **hierarchy**, while the current paper achieves it purely through the optimizer, which is interesting.

(2) One more study should be added to the end to disentangle the accuracy of the solution from the sampling error. Formally I would say this achieves a machine precision **residual**, not **error**. The final plot showing that the error is of comparable wavelength to the grid spacing is an interesting one. The authors could easily add another plot showing the residual (log x axis) vs the L2 error (log y axis) for verifying collocation densities (h = sqrt(N_collocationpoints)). Alternatively, the authors could sample residual points from a uniform distribution on the square at each step of training to avoid grid dependence altogether. If the authors include this study I will further increase my score.

---

confidence: 4
contribution: 3
initial_rating: 6
presentation: 4
questions: (1) I found that ANaGRAM exactly resembles classical natural gradient descent. I think the main technique is applying SVD to the matrix and adding the cutoff in the singular values. Besides this, it is exactly the natural gradient descent. Do I understand correctly?
soundness: 3
strengths: (1) The work is well motivated, as the adaptive processing for singular values should benefit the natural gradient method.

(2) The paper is well written.

(3) The empirical results are strong that the proposed method significantly improves the accuracy of training PINNs over other optimizers.
summary: The paper proposes AMStraMGRAM, an adaptive multi-cutoff variant of ANaGRAM, which is a natural-gradient-based optimizer for training Physics-Informed Neural Networks (PINNs). ANaGRAM regularizes the pseudo-inverse of the feature matrix by a fixed cutoff in singular values. Based on this, AMStraMGRAM generalizes it by introducing a dynamic, multi-cutoff strategy that tracks the evolution of the reconstruction error (RCE) and singular values during training.
weaknesses: (1) While the intuition behind adaptive cutoff dynamics is persuasive, the paper lacks a formal convergence analysis. Is the algorithm convergent to the stationary points?

(2) AMStraMGRAM requires SVD at every iteration, similar to ANaGRAM. Although it reuses intermediate quantities, a discussion of computational overhead is missing. I think the computational cost for natural gradient is high if the sample size is large.

(3) The algorithm is deterministic. However, due to the limited GPU memory as well as the large sample size, stochastic optimizers are favourable in training PINNs. Is the algorithm valid with stochastic gradients/mini-batch?

---

confidence: 4
contribution: 2
initial_rating: 4
presentation: 2
questions: 1. AMStraMGRAM introduces a new hyperparameter—target accuracy $\varepsilon$. How sensitive is the algorithm to the choice of $\varepsilon$? How does the setting of $\varepsilon$ relate to accuracy?

2. Section 6 of the paper presents an excellent discussion on overfitting in the Allen-Cahn equation, highlighting the need for “co-designing the optimizer and sampler.” A natural question arises: Can AMStraMGRAM seamlessly collaborate with existing adaptive sampling methods, such as Residual-Adaptive Distribution (RAD)?

3. All experiments in this paper are based on standard multilayer perceptrons. However, many modern PINN studies employ specialized network architectures to mitigate issues such as spectral bias, including Fourier feature networks or modulated networks. Is AMStraMGRAM's adaptive strategy compatible with these architectures? Does its performance improvement depend on specific properties of the MLP? How will testing on networks with different spectral characteristics affect the emergence of the “flattening” phenomenon and the algorithm's final performance? This concerns the method's general applicability.

4. Beyond the overfitting observed in the Allen-Cahn equation, have the authors encountered instances of poor algorithmic performance or failure under other settings—such as varying network width/depth, activation functions, or equation parameters? Could the “decision boundary” of AMStraMGRAP be described more systematically—specifically, under what conditions (problem characteristics, hyperparameter selection) does the algorithm operate stably, and under which conditions might it fail? This is crucial for potential users applying the method in practice.
soundness: 3
strengths: 1.	**Identified and defined the “flattening phenomenon” in the training dynamics**: by introducing the novel metric of “reconstruction error,” it clearly reveals the underlying mechanism behind the sudden drop in the loss function during the late stages of ANaGRAM training.

2.	**Proposed AMStraMGRAM Algorithm**: By dynamically tracking the intersection point between the RCE and singular value curve, the algorithm intelligently selects the truncation rank, thereby avoiding issues of “premature flattening” (training stagnation) or “delayed convergence” (insufficient accuracy) caused by fixed truncation.
summary: This paper focuses on enhancing the optimization process of Physical Information Neural Networks (PINNs). Specifically, it improves ANaGRAM, an advanced natural gradient-based optimizer. Through an in-depth analysis of ANaGRAM's training dynamics, this paper identifies and defines a critical phenomenon termed “Flattening.” Building upon this theoretical insight, the authors propose AMStraMGRAM, an adaptive multi-cutoff strategy that dynamically and automatically adjusts the cutoff value, thereby eliminating the need for manual tuning of this hyperparameter.

The contributions of this paper can be divided into three closely interrelated aspects:
1) **Discovering and Explaining the “Flattening Phenomenon”**: The paper introduces a novel metric called “Reconstruction Error (RCE).” Through RCE, it reveals the training dynamics of ANaGRAM: Fixed truncation hinders or delays this process, whereas adaptive strategies can actively guide and leverage it.
2) **Proposing the AMStraMGRAM Adaptive Optimizer**: The algorithm can dynamically determine the optimal truncation value for each iteration by tracking the intersection point between the RCE and singular value curve in real time.
3) **Empirical Performance**: On multiple classical PDE benchmark problems, AMStraMGRAM reduces errors by several orders of magnitude compared to ANaGRAM, achieving machine accuracy on some problems and outperforming other types of advanced optimizers.
weaknesses: The main issues with the paper are as follows:

1.	The entire effort of this paper focuses on the adaptation of a single hyperparameter (the cutoff value) within the optimizer. While executed with great skill, this represents essentially a refinement in algorithmic engineering rather than a breakthrough in scientific concepts.

2.	While standard benchmark PDEs such as Laplace, Burgers, and Allen-Cahn equations are commonly used for testing, they lack more challenging cases. PINNs typically perform poorly for equations with solutions exhibiting large variations or strong multiscale characteristics. Without such tests, it remains unclear whether the method enhances the model's intrinsic fitting capability or merely finds smoother solutions more efficiently.

3.	All experiments were conducted within regular domains (e.g., rectangular domains). However, the vast majority of problems in science and engineering are defined over complex geometries (e.g., domains with internal cavities, non-convex regions, or manifolds). On complex domains, collocation sampling and boundary condition treatment themselves pose significant challenges. The paper urgently needs to incorporate experiments on irregular domains to test the robustness of its methods under more realistic settings.

4.	The paper demonstrates AMStraMGRAM's significant advantage in final accuracy but overlooks the computational time required to achieve this precision.

---

confidence: 4
contribution: 2
initial_rating: 2
presentation: 2
questions: 1. See weakness above.

2.How does the float precision affect the results of AMStraMGRAM? The L2 error is very small. The author uses 64p in Jax, while many works conduct experiments with default 32p in Pytorch.

3.In line 392, why the author compare ANaGRAM with the same equation but with ‘modified dataset’?

4.Is the feature matrix the same as the jacobian matrix J in arXiv:2408.00573? What is the difference between your (10) with existing natural gradient methods (arXiv: 1905.11675, arXiv:2408.00573) besides the SVD decomposition instead of the pesodu-inverse of the jacobian matrix?

5.Can it be combined with approximated natural gradient descent (NGD) methods like F-KAC?  In general deep learning, K-FAC is often used instead of 'exact' NGD.In particular, theoretical analysis of K-FAC has already been developed for shallow neural networks.
soundness: 3
strengths: 1. The analysis of the "flattening phenomenon" via the Reconstruction Error (RCE) metric is novel.

2. The paper is theoretically depth as it provides a rigorous functional-analytic framework, explaining why regularization is necessary and drawing a clear distinction between ridge and cutoff regularization.

3. The proposed AMStraMGRAM algorithm is a direct simple but effective modification of the ANaGRAM. The cutoff of the rank is dynamically adjusted, based on the theoretical analysis.
summary: This paper introduces AMStraMGRAM, an adaptive optimization algorithm designed to improve the training of Physics-Informed Neural Networks (PINNs). The work builds directly upon ANaGRAM, a natural-gradient descent method that uses a fixed cutoff to regularize the singular value decomposition (SVD) of the feature matrix during training. AMStraMGRAM addresses the “flattening phenomenon” by dynamically adjusting the cutoff at each iteration (compared to the fixed cutoff in AnaGRAM). The algorithm tracks the intersection between the RCE curve and the singular values, switching between an "intersection rank" and a "precision rank" to ensure the training progresses efficiently towards a target precision. This adaptive, multi-cutoff strategy aims to avoid early stagnation and exploit the flattening phenomenon to force a final drop in the training loss. The method is validated on several benchmark PDEs (Heat, Laplace, Allen-Cahn), showing significant performance improvements over ANaGRAM.
weaknesses: 1. While the theoretical framework is impressive, it primarily provides an explanation for the necessity of regularization and an interpretation of the method. It does not provide rigorous convergence analysis of the training loss for the adaptive cutoff scheme itself. It seems that convergence analysis has been established for natural gradient method to train PINN problems (arXiv:2408.00573) at the same time. So the main incrementation of AMStraMGRAM to ANaGRAM is the adaptive cutoff strategy of the rank, which is too limited for a top conference.
 
2. The experiments are not comprehensive in today's PINNs’ community. Although the AMStraMGRAM shows ‘L2 error improvements of up to 8 orders of magnitude’ compared to ANaGRAM, it is on very simple Heat and Laplace equation with simple smooth solutions, which is less meaningful. No experiments conducted on complex PINN problems (such as the NS equation with turbulence, KS equation with chaotic property, as in PINNacle NIPS24’) to show the practical improvements over other optimizers on practical problems.


Decision: Reject