Published as a conference paper at ICLR 2027
L EARNING P HYSICS FROM AN I MPERFECT A NCESTOR S. Mohammad Mousavi1,2 , Teeratorn Kadeethum2 , Nikolaos Bouklas1 , Somdatta Goswami3 1
Sibley School of Mechanical and Aerospace Engineering, Cornell University AI Lab, Siemens Energy 3 Department of Civil and Systems Engineering, Johns Hopkins University 2
arXiv:2609.24947v1 [cs.LG] 21 Sep 2026
A BSTRACT Neural operators evaluate parametric partial differential equations (PDEs) cheaply but degrade sharply outside their training distribution. Physics-informed neural networks (PINNs) avoid dependence on labeled data, yet their optimization can be basin-fragile: when the governing residual admits multiple solutions, a PINN trained from scratch may converge to a physically incorrect state despite achieving a small residual. We show that these failure modes can be addressed jointly: an imperfect NO provides the structural prior needed to place a PINN in the correct solution basin, while the PDE residual refines the solution beyond the operator’s accuracy. We introduce a three-stage framework that freezes the spatial basis of a physics-informed NO, extrapolates its solution branch to an out-of-distribution parameter using a polynomial continuation prior, and distills the resulting field into a fresh PINN. The NO (teacher) need not be accurate at the target; it transfers solution-branch information, while PDE residual minimization in the PINN (student) governs convergence. We evaluate the framework on three nonlinear PDEs: 1D viscous Burgers, 2D steady Allen–Cahn near a pitchfork bifurcation, and 2D steady lid-driven cavity flow. For Allen–Cahn, where the trivial solution satisfies the PDE residual exactly, a standard PINN consistently collapses to the trivial zero branch, whereas distillation from the crude extrapolated operator recovers the non-trivial branch that matches the finite-difference reference. For the lid-driven cavity, extrapolating to a Reynolds number of Re = 3200 accelerates convergence to the correct physical state, achieving competitive accuracy using fewer parameters and optimization steps than recent literature baselines. These results establish a simple principle: an NO need not accurately predict the solution to be useful; it only needs to identify the correct basin from which PINN optimization can recover it.
1
I NTRODUCTION
Parametric partial differential equations (PDEs) underpin many scientific and engineering workflows, where the same governing equations must be solved repeatedly across varying parameters. Neural operators (NOs) amortize this cost by learning parameter-to-solution maps, enabling inexpensive inference after training. However, their accuracy can deteriorate sharply outside the training distribution, while adapting them to new parameters typically requires additional simulations or labeled data. Physics-informed neural networks (PINNs) provide a complementary alternative by enforcing the governing equations directly, but their optimization is sensitive to initialization and can converge to physically incorrect solutions when the PDE residual admits multiple minima. In such cases, a small residual does not guarantee selection of the physically relevant solution branch. These failure modes suggest a complementary division of roles: an NO can provide structural information about the solution manifold, while a PINN can enforce the governing physics at the target parameter. We therefore ask whether an imperfect, out-of-distribution NO can guide a PINN without constraining its final accuracy. We show that it can. The key is to use the NO not as a predictor of the final solution, but as a prior for solution-basin selection. To this end, we introduce a three-stage framework: a physics-informed NO is trained and its spatial basis frozen; the branch is extrapolated to an out-of-distribution target using a polynomial continuation prior, without any reference solution; and the extrapolated field is distilled into a freshly 1
Published as a conference paper at ICLR 2027
initialized PINN, whose subsequent residual minimization governs final accuracy. The resulting NO acts as a teacher, but unlike conventional knowledge distillation, it is not required to provide an accurate prediction. Its role is to transfer solution-branch information and place the PINN in a favorable basin of attraction; subsequent minimization of the PDE residual allows the PINN to refine the field and, in principle, surpass the accuracy of the teacher. This separation of roles is central to the proposed framework. We evaluate the framework on three nonlinear PDEs. First, the 1D viscous Burgers equation establishes the baseline on a problem with a single solution branch. Second, the 2D steady Allen–Cahn equation probes basin selection near a pitchfork bifurcation, a scenario where the trivial solution u ≡ 0 satisfies the residual exactly. Third, to test substantial out-of-distribution capabilities, we apply the framework to a 2D steady lid-driven cavity, extending an operator trained on a Reynolds number range of Re ∈ [100, 800] to Re = 3200. Finally, an ablation study isolates the underlying mechanism by removing the extrapolation stage, designed to determine whether the critical challenge in these failure modes is strictly residual optimization or physical basin selection.
2
R ELATED WORK
Operator learning. NOs learn maps between function spaces and are, in principle, discretizationinvariant (Kovachki et al., 2023). DeepONet (Lu et al., 2021a; Haghighat et al., 2025) realizes this through a trunk–branch decomposition; the Fourier neural operator (Li et al., 2020a) and the earlier graph-kernel network (Li et al., 2020b) operate in spectral and graph domains respectively, with Wavelet (Tripura & Chakraborty, 2023),Laplace (Cao et al., 2024a) and other variants targeting localized features (Lu et al., 2022; Cao et al., 2024b); no single architecture dominates across problem classes. Physics-informed DeepONets (Wang et al., 2021c; Mandl et al., 2025) replace data losses with residual losses but retain the fixed-parameter-box limitation, and transfer-learning approaches under conditional shift (Weiss et al., 2016; Goswami et al., 2022) extend a trained operator to new domains only with labeled target data at the shifted parameters. Our Part 2 reuses a frozen basis at an extrapolated parameter without any labeled data and without retraining the trunk. PINNs and their failure modes. Physics-informed neural networks (Raissi et al., 2019) (alongside related residual-minimization schemes such as the deep Galerkin method (Sirignano & Spiliopoulos, 2018), and surveyed broadly in (Karniadakis et al., 2021)) solve PDEs without labeled data but are notoriously hard to optimize. A substantial literature documents training pathologies in such coordinate-based and physics-informed networks: spectral bias (Tancik et al., 2020; Wang et al., 2021b), unbalanced loss gradients (Wang et al., 2021a), a neural-tangent-kernel account of why training stalls (Wang et al., 2022), causality violations in time (Wang et al., 2024b), and instability at depth (Wang et al., 2024a). Remedies span loss re-weighting and self-adaptive weights (Wang et al., 2021a; 2022; McClenny & Braga-Neto, 2023), residual-based adaptive sampling (Wang et al., 2026; Wu et al., 2023), adaptive activations (Jagtap et al., 2020), gradient-enhanced objectives (Yu et al., 2022), architectural changes (Wang et al., 2024a), and consolidated training recipes (Wang et al., 2023a). Less studied are failures of basin selection: on problems whose residuals admit multiple physically distinct minima (bifurcating branches, symmetry-related states, trivial constant solutions) a from-scratch PINN can converge to a spurious solution branch that does not match the reference field, even when the residual is near zero. While temporal and adaptive sampling remedies address erroneous local minima in specific PDEs like the Allen–Cahn equation (Wang et al., 2024b; 2026), our results address broader failures of basin selection that the other remedies above do not target. Solving the lid-driven cavity. The lid-driven cavity is the canonical incompressible-flow benchmark: multigrid finite-difference tables (Ghia et al., 1982), spectral benchmarks (Botella & Peyret, 1998), high-Reynolds fine-grid studies (Erturk et al., 2005), and an 8192 × 8192 Richardsonextrapolated solution (Marchi et al., 2021) remain the standard yardstick. PINNs struggle far more: rising Reynolds number sharpens the corner singularities and destabilizes residual minimization. Wang et al. (2023b) show vanilla PINNs at Re = 2000–5000 admit two residual minimizers (one matching DNS, one unphysical) restored to uniqueness via an added entropy-viscosity term. Stateof-the-art single-case accuracy comes from the residual-adaptive PirateNet architecture (Wang et al., 2
Published as a conference paper at ICLR 2027
2024a), our point of comparison at Re = 3200. These works solve one cavity at a time; we instead extrapolate a frozen operator basis, letting a crude teacher supply the physical basin. Knowledge distillation. Traditionally, knowledge distillation originates in model compression (Hinton et al., 2015), where a smaller student network is trained to match a larger, highly accurate teacher. In this classical setup, the student’s performance is inherently capped by a nearperfect teacher. Our approach differs fundamentally: the teacher is a neural operator providing only a crude out-of-distribution estimate, and the student is explicitly expected to surpass it through subsequent physics-based optimization. Accordingly, the distillation weight decays and ultimately vanishes rather than remaining active throughout training. This shares intuition with scientific computing approaches that use learned models as initial guesses, including warm-start strategies for PINNs (Wang et al., 2024a) and iterative solvers (Eshaghi et al., 2026). However, in contrast to methods that either apply the prior only at initialization (Eshaghi et al., 2026) or rigidly freeze the prior representations (Desai et al., 2021; Goswami et al., 2020; Chakraborty, 2021), our framework retains the teacher as a persistent but decaying constraint. This provides necessary guidance while the student explores the solution landscape, progressively relaxing the imperfect teacher’s influence before transferring full control to the governing physics.
3
M ETHODOLOGY
3.1
S HARED ARCHITECTURE
Let s(x; p) denote the unknown solution field on domain Ω(p) for parameters p ∈ Rdp . We employ strong boundary enforcement by defining the ansatz as: s(x; p) = g(x; p) + M(x)
K X
bk (p) φk (x),
(1)
k=1
where {φk }K k=1 is a spatial basis produced by a trunk network, K is the number of trunk basis functions (equivalently, the number of branch outputs), and bk (p) are the branch coefficients, collected into the vector b(p) = (b1 (p), . . . , bK (p)). Here, g(x; p) is a known function satisfying the boundary conditions (the boundary lift), and M(x) is a distance multiplier that evaluates to zero on the boundary ∂Ω. Enforcing boundary data exactly via this construction eliminates the need for boundary loss terms and stabilizes training (Lu et al., 2021b). For the Allen–Cahn problem, the lift g also incorporates a non-zero initial profile. Because the trivial solution (s = 0) satisfies the Allen–Cahn PDE exactly, this addition is required to prevent the network from collapsing into the non-physical zero branch. Consequently, no boundary loss term enters training, and the objective relies strictly on the mean-squared interior PDE residual, R (the problem-specific residual operator, given explicitly for each of the three PDEs in Appendix A): LPDE (p) =
2 R s(·; p) 2
,
L (Ω(p))
(2)
averaged over parameter batches drawn from a bounded base box B ⊂ Rdp . Throughout, we report a relative-L2 accuracy against a reference field sref (analytical or finite-difference, per problem): Acc = 1 − εL2 ,
εL2 =
∥ŝ − sref ∥L2 , ∥sref ∥L2
(3)
where ŝ is the predicted field. Table 2 reports εL2 (labeled error) for the operator and both PINNs, alongside the absolute performance gain (labeled Improvement), defined as εL2 ,baseline − εL2 ,distilled . Everywhere else in the paper, performance is reported in terms of Acc (labeled accuracy). All three problems in §4 instantiate the ansatz equation 1 and train against the residual objective equation 2; they differ in three architectural aspects summarized in Table 1: (i) the procedure used to construct the anchor coefficients after freezing the spatial basis; (ii) the order to which the multiplier M vanishes on the boundary; and (iii) the spatial bandwidth of the trunk features. Additional implementation-level differences, including residual discretization and a problem-specific regularizer, are specified in Appendix A.6. 3
Published as a conference paper at ICLR 2027
Table 1: Problem-specific design choices across the evaluated PDEs.
Field / operator Parameters p Order / dim. Design choice (i): Part 1 anchors Design choice (ii): boundary lift Design choice (iii): trunk bandwidth Part 2 and 3 target
3.2
Burgers
Allen–Cahn
Cavity
scalar u; uux − νuxx = 0 (ν, L) 2nd-order, 1D self-supplied from branch
scalar u; ∇2 u + λ(u−u3 )=0 (λ, H) 2nd-order, 2D data-free Gauss–Newton on frozen basis seed lift ℓ0 + 1st-order mult. mild, isotropic (σ=[3, 3])
stream ψ; vorticity transport (Re, H) 4th-order, 2D FD reference reprojected onto basis (8 solves) 2nd-order multiplier + lid lift moderate, anisotropic σ = [σξ , ση ] larger Re (a geometry variant is reported in Appendix E)
1st-order multiplier moderate, isotropic σ = 6 larger L
λ below base box (near onset)
PART 1: BASE OPERATOR
Trunk and branch are trained jointly by minimizing equation 2 over p ∼ Unif(B). Upon convergence, the trunk is frozen, rendering {φk } a fixed function space. An anchor set {ci } is then constructed at a small collection of base points {pi } ⊂ B. These anchors act as trusted reference coefficients within the training domain, which will later be interpolated to form a structural prior for out-of-distribution extrapolation in Part 2. The three problems differ only in how these anchors are formed, exhibiting a spectrum of increasing reliance on external information (Table 1, Design choice i). Burgers stores the trained branch output directly. Allen–Cahn refines it through a damped Gauss–Newton solve of the residual on the frozen basis, utilizing only the residual and its exact Jacobian (with no reference field). The cavity computes eight finite-difference reference fields at Re ∈ {100, 200, . . . , 800} with fixed H=1, reprojecting each onto the frozen basis via ridge least squares in the velocity metric. These eight solves represent the only external data the main cavity pipeline ever consumes.1 3.3
PART 2: EXTRAPOLATION BEYOND THE BASE BOX
While parameter continuation is common in PINN curricula to avoid failure modes (Krishnapriyan et al., 2021), unguided extrapolations often drift into non-physical configurations that spuriously minimize the residual. To prevent this, the base anchors are interpolated by a low-degree polynomial ccont (p) fit independently for each of the K modes, acting as a continuation prior. To reach a target p⋆ ∈ Rdp \ B, the branch is retrained at p⋆ under two regularizing mechanisms that require no reference solution, keeping the extrapolation entirely data-free. Bounded correction. At the target, the branch network is re-optimized to output a raw coefficient correction ∆b(p) (analogous to b(p) in equation 1, but produced by a fresh Part-2 branch head rather than the Part-1 branch). Rather than using ∆b(p) directly as the coefficients, it is added to the continuation prior as a bounded perturbation, bcomp (p) = ccont (p) + α tanh ∆b(p) , (4) with α representing a hard cap on the departure from the prior. If the prior is correct, ∆b = 0 is optimal. Gram-weighted continuation penalty.
Training minimizes 2
Ltotal = LPDE (bcomp ) + β bcomp − ccont G(H) ,
(5)
1 The aspect-ratio extrapolation reported separately in Appendix E trains an independent base operator over the full (Re, H) box with 12 anchors on {100, 200, 300, 400} × {1.0, 1.5, 2.0}; see Table 7 for both configurations.
4
Published as a conference paper at ICLR 2027
where G(H) is the field-metric Gram matrix of the frozen basis, assembled on the physical geometry at aspect ratio H via the isoparametric map of Appendix A (for the 1D and fixed-geometry problems, G(H) reduces to a constant matrix G). Penalizing coefficient deviations in the field metric (rather than the raw coefficient space) ensures the objective strictly penalizes modifications that alter the physical field. The initial scalar weight β is calibrated once per target so that the penalty term contributes a fixed, predetermined fraction of the total initial objective loss. From that calibrated value, β is either held constant for the remainder of the step or annealed exponentially toward a small non-zero floor over the target’s training budget, letting the PDE residual increasingly govern convergence as training proceeds; the problem-specific choices are given in Appendix C. The fundamental extrapolation operation projects the base anchors to a new target parameter. If the target p⋆ is sufficiently close to the training domain, this is executed in a single shot: the prior is fit once on the Part 1 anchors, and the branch is retrained directly at p⋆ using equation 4 and equation 5. However, when a single step spans too large a parameter distance, the extrapolation employs an incremental self-supply strategy. In this multi-stage approach, the model extrapolates to intermediate targets; upon convergence, these intermediate solutions are appended to the anchor pool, and the prior is refitted before advancing to the next stage. For example, the principal cavity Re-extrapolation utilizes a four-stage continuation (with a 1D prior over 1/Re), the geometry variant fits a 2D prior over (1/Re, H) (Appendix A), and Allen–Cahn requires a two-stage sequence (Appendix A.6). The goal of this stage is not high target accuracy, but to generate a coefficient field sufficient to initialize Part 3, which then refines the solution. The Allen–Cahn problem illustrates this: near √ the pitchfork bifurcation, the true branch amplitude scales as λ − λ1 (with λ1 the critical rate), a form the analytic polynomial prior cannot represent. As a result, the extrapolated field degrades predictably at the near-onset target. 3.4
PART 3: G UARDED D ISTILLATION INTO A F RESH S OLVER
At the target p⋆ , the composited field from equation 4 is passed to a fresh MLP-PINN uϕ trained from scratch on the physics residual. Crucially, this student network abandons the trunk-branch decomposition and the specialized boundary lift utilized in Part 1. The training loss over optimization step t is defined as: L(ϕ; t) =
2
2
R[uϕ ] L2 (Ω) + λBC LBC + w(t) uϕ − uteacher Γ(x) ,
(6)
where LBC penalizes boundary-condition violations, uteacher is the Part 2 prediction, w(t) is the distillation schedule, and Γ(x) is an optional spatial reliability gate (Appendix B). To demonstrate modularity, we employ two strategies. Allen–Cahn and the primary Re-extrapolated cavity use an ungated (Γ ≡ 1) piecewise-linear w(t) that vanishes during the final training phase. Conversely, the Burgers and geometry-extrapolated cavity configurations utilize a physics-adaptive gate: w(t) tracks the student’s residual, while Γ(x) locally attenuates the teacher’s influence at high-residual points to tolerate localized teacher inaccuracies. Both schedules strictly vanish before training concludes, leaving the PDE and boundary constraints to govern final convergence. Figure 1 summarizes the complete three-stage framework.
4
E XPERIMENTS
We evaluate the framework through a progressive sequence of experiments designed to address four core objectives. First, we use the viscous Burgers equation to establish whether an approximate, outof-distribution teacher can accelerate a PINN (§4.1). Second, we analyze the Allen–Cahn equation to determine why this transfer succeeds in basin selection where pure residual minimization fails (§4.2). Third, we evaluate the lid-driven cavity to test if the mechanism scales to complex, highdimensional systems far outside their training distribution (§4.3). Finally, an ablation study isolates precisely which component of the framework drives this success (§4.4). 4.1
E XPERIMENT 1 (BASIC TRANSFER ): V ISCOUS B URGERS (1D)
The 1D viscous Burgers equation establishes the fundamental transfer mechanism in a controlled setting featuring a unique monotone solution manifold, one internal transition layer, and a Cole–Hopf 5
Published as a conference paper at ICLR 2027
Figure 1: Three-stage framework schematic. Table 2: Performance summary across the three evaluated PDEs. Burgers
Allen–Cahn
Cavity
(ν, L) ∈ [0.10, 0.20] × [1, 2] (ν, L) = (0.05, 4) 2×
(λ, H) ∈ [35, 80] × [1, 2] (λ, H) = (12.5, 3) 0.35× lower box edge
(Re, H) ∈ [100, 800] × {1} (Re, H) = (3200, 1) 4×
Operator: in-distribution error Operator: OOD error Baseline PINN error Distilled PINN error (Ours)
0.003 ± 0.002 0.12 0.4123 ± 0.1601 0.0209 ± 0.0043
0.0015 ± 0.0004 0.78 1.00 ± 0.00 0.0037 ± 0.001
0.2586 ± 0.1294 0.53 0.9800 ± 0.0248 0.0663 ± 0.0099
Improvement
0.3914 ± 0.1635
0.9963 ± 0.001
0.9137 ± 0.0341
Training range Target OOD∗ distance
∗
OOD stands for out of distribution.
closed form for exact evaluation. The base operator is trained on L ∈ [1, 2] in Part 1 and extrapolated to L = 4 at ν = 0.05 in Part 2, extending twice beyond the maximum training length. Figure 2(a) illustrates the frozen-trunk operator at this out-of-distribution target: it captures the correct shock location and outer branches but exhibits an overshoot at the transition layer. This represents the precise approximate-yet-structurally-accurate prior the framework is designed to exploit. Across n = 5 independent seeds, the distilled student removes the overshoot and matches the analytical solution to graphical accuracy in every run (Figure 2(f)). It reaches 97.91% ± 0.43% mean accuracy versus 58.77% ± 16.01% for the from-scratch baseline, an improvement of 39.14% ± 16.35% that holds in all 5/5 seeds. Beyond the higher mean, distillation markedly reduces run-to-run variance: the baseline’s final accuracy is highly seed-dependent and occasionally collapses to a poor local solution, while the distilled student consistently converges to a tight neighborhood of the true solution regardless of initialization. Distillation also converges faster and to a lower final residual, with the reliability gate closing as the student’s residual falls (Figure 2 d). 4.2
E XPERIMENT 2 (BASIN SELECTION ): S TEADY A LLEN –C AHN (2D)
This experiment demonstrates how the teacher resolves a solution-branch ambiguity that pure resid2 2 ual minimization cannot. At λ close to the pitchfork √ bifurcation λ1 (H) = π (1 + 1/H ), the ⋆ principal positive branch u > 0 exhibits an O( λ − λ1 ) amplitude, while the trivial solution (u ≡ 0) satisfies the residual exactly. We evaluate the framework at λ = 12.5 and H = 3 (where λ1 (3) ≈ 10.97). This target resides approximately 14% above the onset, with a true peak amplitude of ≈ 0.45. Figure 3 (a) illustrates this dynamic. Baseline PINN collapses to the flat u ≡ 0 field, identifying a mathematically correct but physically trivial global minimizer; across n = 5 seeds, this collapse is total and deterministic, with 0.00% ± 0.00% accuracy and a peak amplitude of only 0.0002 ± 0.0002 against a true peak of 0.4530. The Part 2 operator (second panel) correctly identifies the single-hump √ structure but over-predicts the amplitude, as the polynomial prior cannot precisely capture the λ − λ1 onset behavior, giving it only 32% accuracy on its own. Distilling this operator provides the necessary prior to pull the student into the correct basin of attraction despite 6
Published as a conference paper at ICLR 2027
Figure 2: Burgers equation training metrics aggregated over n = 5 seeds. Panel (a) illustrates the Part 2 operator extrapolation versus the analytical reference. The remaining five panels detail the Part 3 distilled PINN performance.
the teacher’s own imprecision. Subsequent residual minimization sharpens the field: Distilled PINN matches the finite-difference reference in both shape and amplitude, reaching 99.63%±0.10% mean accuracy across seeds. Figure 3 (b) confirms this quantitatively: accuracy climbs as the distillation weight anneals to zero, and the centerline profile aligns with the reference. This confirms the central hypothesis of the framework: an imperfect operator is sufficient because its role is solely to identify the correct physical basin, not to provide the exact solution. Here, a teacher accurate to only 32% still reliably rescues every seed from the trivial minimizer.
Figure 3: Allen–Cahn training metrics (λ = 12.5, H = 3), aggregated over n = 5 seeds, versus finite difference reference.
4.3
E XPERIMENT 3 (L ARGE OOD EXTRAPOLATION ): L ID - DRIVEN CAVITY (2D)
To demonstrate scalability, the lid-driven cavity introduces a fourth-order, two-dimensional system extrapolated substantially outside its base box. The primary target is set to Re = 3200 (far outside 7
Published as a conference paper at ICLR 2027
the Re ∈ [100, 800] training domain) to match the state-of-the-art single-case PINN benchmark established by Wang et al. (2024a). A secondary experiment, extrapolating aspect ratio rather than Reynolds number, is reported in Appendix E.
Figure 4: Lid-driven cavity at Re = 3200. (a) Predicted u and v fields. (b) Velocity magnitude and streamlines for the distilled and baseline students. (c) Centerline profiles at x = 0.5 and y = 0.5 against Ghia et al. (1982). Fields and streamlines show the ensemble mean across n = 5 seeds. Figure 4 illustrates the distilled student at Re = 3200. Aggregated across n = 5 independent seeds, the u and v fields (a) successfully capture the primary vortex and thin wall jets. The streamlines (b) reproduce the primary recirculation and the two lower corner vortices in the distilled model, features that the baseline completely fails to capture. Consequently, the centerlines (c) tightly align with the benchmark points established by Ghia et al. (1982) and Marchi et al. (2021), whereas the baseline deviates severely.
Figure 5: Lid-driven cavity (Re = 3200) training metrics. Solid lines and shaded regions represent µ ± σ across n = 5 independent seeds. Figure 5 details the optimization trajectory: the baseline stalls at near-zero accuracy with a noisy residual, while the distilled student steadily climbs to high accuracy at a lower final residual, utilizing the identical architecture and step budget. Quantitatively, the baseline PINN entirely fails to converge to the correct physical state, yielding a mean accuracy of just 2.00% ± 2.48%. In contrast, the distilled student consistently resolves the solution-branch ambiguity, achieving a mean accuracy of 93.37% ± 0.99%. This represents a robust accuracy improvement of 91.4% ± 3.4%, an acceleration that holds consistently across all 5/5 initialized seeds. Beyond the higher mean performance, distillation markedly reduces run-to-run variance, confirming that the operator initialization effectively guards against the inherent sensitivity of PINNs to random parameter initialization. For computational context, Table 3 compares the Part 3 student with Wang et al. (2024a): our student uses roughly 20% of the parameters and 40% of the optimization steps (250K vs. 630K), but more point evaluations. End-to-end, Parts 1–3 require 546K optimization steps, but note that Part 1 is a one-time base-box cost reused across subsequent targets. 4.4
E XPERIMENT 4 (M ECHANISTIC ABLATION ): T HE ROLE OF PART 2
Section 4.3 presents the full framework, where Part 2 extrapolates the frozen-trunk operator over a four-stage Re sequence (§3.3) before distillation. This raises a key ablation: does incremental 8
Published as a conference paper at ICLR 2027
Table 3: Part 3 computational comparison at Re = 3200 with Wang et al. (2024a).
2
Rel. L error against Ghia et al. (1982) Depth × width Trainable parameters Total training steps Total point-evaluations
Distilled PINN (Ours)
PirateNet
0.0663 ± 0.0099 5 × 256 297 K 250 K 6.5 × 109
0.0421 18 × 256 ∼1.32 M 630 K 2.6 × 109
continuation uniquely determine basin selection, or is a direct evaluation of the Part 1 operator at Re = 3200 an adequate prior? To test this, we keep the Part 3 architecture, the piecewise-linear distillation schedule, and the 250K-step training budget fixed, but initialize the teacher coefficients by directly evaluating the Part 1 branch at the target: c = bPart 1 (Re=3200, H=1). This extrapolates the neural operator fourfold beyond its training range (Re ∈ [100, 800]) without the regularizing polynomial continuation prior equation 4. We then compare this ablated model with the full framework and an unguided baseline PINN. Table 4 and Figure 8 summarize the ensemble results across n = 5 independent seeds. Relying exclusively on the unguided Part 1 teacher consistently fails to identify the correct physical basin. As depicted in the mean velocity fields, the student’s primary vortex is mislocated and the lower corner vortices are entirely absent, resulting in severe deviations from the benchmark (Ghia et al., 1982). Quantitatively, the ablated model yields a highly variable mean accuracy of 27.09% ± 13.49%. While this slightly outperforms the catastrophic collapse of the baseline PINN (2.00% ± 2.48%), it falls drastically short of the full framework’s 93.37% ± 0.99%. Table 4: Ablation of the Part 2 continuation stage for the cavity at Re = 3200, H = 1.
Teacher Relative L2 accuracy PDE residual
Full framework
No Part 2
Baseline PINN
Continued operator 93.37% ± 0.99% 6.3 × 10−6
Direct Part-1 operator 27.09% ± 13.49% 1.3 × 10−5
None 2.00% ± 2.48% 2.9 × 10−5
Crucially, this represents a basin-selection failure rather than an optimization failure. The mean PDE residual for the ablated student converges cleanly to 1.3 × 10−5 , which is strictly comparable to both the full framework’s 6.3 × 10−6 and the baseline PINN’s 2.9 × 10−5 (Table 4). This demonstrates that all three networks effectively minimize the governing physical equations, yet both the unguided baseline and the ablated model converge to steady states that deviate substantially from the Ghia et al. (1982) reference. This behavior perfectly characterizes the basin-selection vulnerability the proposed framework is designed to circumvent (§4.2).
5
C ONCLUSION
We introduced a three-stage framework that decouples the distinct failure modes of operator learning and physics-informed neural networks. The central finding of this work, demonstrated across three nonlinear PDEs, is that an approximate neural operator need not be highly accurate to be effective; it needs only to initialize a student network within the correct basin of attraction, after which residual minimization governs final optimization. Our analyses confirm this mechanism, demonstrating that distilled students robustly recover the true physics in scenarios where standard PINNs deterministically collapse or stall at non-physical steady states. Future directions include extending the continuation prior to better capture bifurcation onsets, evaluating the framework in time-dependent multi-physics environments, and integrating this distillation strategy with state-of-the-art architectures (Wang et al., 2024a) to yield compounding improvements in optimization stability. 9
Published as a conference paper at ICLR 2027
R EPRODUCIBILITY STATEMENT All experimental hyperparameters are given in Appendix C. The finite-difference reference solvers used for scoring and, for the cavity problem, anchor formation are standard, and their configurations are specified in Appendix D. An anonymized implementation and pretrained models are available at https://anonymous.4open.science/r/ PINN-Distillation-prototype-D3AF/. AI USE STATEMENT During the preparation of this work, the authors used Claude (Anthropic) for assistance in writing the code used for analysis and evaluation, and to draft and revise manuscript text; Gemini (Google) was used to help improve manuscript text. The authors did not use generative AI tools for research ideation or experimental design; the three-stage framework and the extrapolation experiments described in this paper reflect the authors’ own conceptual contributions. The authors also did not use generative AI tools for data analysis or interpretation of results. All AI-assisted work was reviewed before inclusion: LLM-generated code was verified and tested for correctness by the authors, and AI-drafted and AI-revised manuscript text was reviewed and edited by all the authors for accuracy and clarity. The authors take full responsibility for the content of the published article.
R EFERENCES O. Botella and R. Peyret. Benchmark spectral results on the lid-driven cavity flow. Computers & Fluids, 27(4):421–433, 1998. doi: 10.1016/S0045-7930(98)00002-4. Qianying Cao, Somdatta Goswami, and George Em Karniadakis. Laplace neural operator for solving differential equations. Nature Machine Intelligence, 6(6):631–640, 2024a. Qianying Cao, Somdatta Goswami, Tapas Tripura, Souvik Chakraborty, and George Em Karniadakis. Deep neural operators can predict the real-time response of floating offshore structures under irregular waves. Computers & Structures, 291:107228, 2024b. Souvik Chakraborty. Transfer learning based multi-fidelity physics informed deep neural network. Journal of Computational Physics, 426:109942, 2021. Shaan Desai, Marios Mattheakis, Hayden Joy, Pavlos Protopapas, and Stephen Roberts. One-shot transfer learning of physics-informed neural networks. arXiv preprint arXiv:2110.11286, 2021. Ercan Erturk, Thomas C Corke, and Cihan Gökçöl. Numerical solutions of 2-d steady incompressible driven cavity flow at high reynolds numbers. International journal for Numerical Methods in fluids, 48(7):747–774, 2005. Mohammad Sadegh Eshaghi, Cosmin Anitescu, Navid Valizadeh, Yizheng Wang, Xiaoying Zhuang, and Timon Rabczuk. Nows: Neural operator warm starts for accelerating iterative solvers. Computer Methods in Applied Mechanics and Engineering, 458:118989, 2026. UKNG Ghia, Kirti N Ghia, and CT Shin. High-re solutions for incompressible flow using the navierstokes equations and a multigrid method. Journal of computational physics, 48(3):387–411, 1982. Somdatta Goswami, Cosmin Anitescu, Souvik Chakraborty, and Timon Rabczuk. Transfer learning enhanced physics informed neural network for phase-field modeling of fracture. Theoretical and Applied Fracture Mechanics, 106:102447, 2020. Somdatta Goswami, Katiana Kontolati, Michael D Shields, and George Em Karniadakis. Deep transfer operator learning for partial differential equations under conditional shift. Nature Machine Intelligence, 4(12):1155–1164, 2022. Ehsan Haghighat, Mohammad Hesan Adeli, S Mohammad Mousavi, and Ruben Juanes. Stonet: A neural operator for modeling solute transport in micro-cracked reservoirs. Advances in Water Resources, pp. 105046, 2025. Geoffrey Hinton, Oriol Vinyals, and Jeff Dean. Distilling the knowledge in a neural network. arXiv preprint arXiv:1503.02531, 2015. 10
Published as a conference paper at ICLR 2027
Ameya D Jagtap, Kenji Kawaguchi, and George Em Karniadakis. Adaptive activation functions accelerate convergence in deep and physics-informed neural networks. Journal of Computational Physics, 404:109136, 2020. George Em Karniadakis, Ioannis G Kevrekidis, Lu Lu, Paris Perdikaris, Sifan Wang, and Liu Yang. Physics-informed machine learning. Nature Reviews Physics, 3(6):422–440, 2021. Nikola Kovachki, Zongyi Li, Burigede Liu, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Neural operator: Learning maps between function spaces with applications to pdes. Journal of Machine Learning Research, 24(89):1–97, 2023. Aditi Krishnapriyan, Amir Gholami, Shandian Zhe, Robert Kirby, and Michael Mahoney. Characterizing possible failure modes in physics-informed neural networks. Advances in neural information processing systems, 34:26548–26560, 2021. Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Fourier neural operator for parametric partial differential equations. arXiv preprint arXiv:2010.08895, 2020a. Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Neural operator: Graph kernel network for partial differential equations. arXiv preprint arXiv:2003.03485, 2020b. Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis. Learning nonlinear operators via deeponet based on the universal approximation theorem of operators. Nature machine intelligence, 3(3):218–229, 2021a. Lu Lu, Raphael Pestourie, Wenjie Yao, Zhicheng Wang, Francesc Verdugo, and Steven G Johnson. Physics-informed neural networks with hard constraints for inverse design. SIAM Journal on Scientific Computing, 43(6):B1105–B1132, 2021b. Lu Lu, Xuhui Meng, Shengze Cai, Zhiping Mao, Somdatta Goswami, Zhongqiang Zhang, and George Em Karniadakis. A comprehensive and fair comparison of two neural operators (with practical extensions) based on fair data. Computer Methods in Applied Mechanics and Engineering, 393:114778, 2022. Luis Mandl, Somdatta Goswami, Lena Lambers, and Tim Ricken. Separable physics-informed deeponet: Breaking the curse of dimensionality in physics-informed machine learning. Computer Methods in Applied Mechanics and Engineering, 434:117586, 2025. Carlos Henrique Marchi, Cosmo Damião Santiago, and Carlos Alberto Rezende de Carvalho, Jr. Lid-driven square cavity flow: A benchmark solution with an 8192× 8192 grid. Journal of Verification, Validation and Uncertainty Quantification, 6(4):041004, 2021. Levi D McClenny and Ulisses M Braga-Neto. Self-adaptive physics-informed neural networks. Journal of Computational Physics, 474:111722, 2023. Ken Perlin. Improving noise. In Proceedings of the 29th annual conference on Computer graphics and interactive techniques, pp. 681–682, 2002. Maziar Raissi, Paris Perdikaris, and George 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, 2019. Justin Sirignano and Konstantinos Spiliopoulos. Dgm: A deep learning algorithm for solving partial differential equations. Journal of computational physics, 375:1339–1364, 2018. Matthew Tancik, Pratul Srinivasan, Ben Mildenhall, Sara Fridovich-Keil, Nithin Raghavan, Utkarsh Singhal, Ravi Ramamoorthi, Jonathan Barron, and Ren Ng. Fourier features let networks learn high frequency functions in low dimensional domains. Advances in neural information processing systems, 33:7537–7547, 2020. 11
Published as a conference paper at ICLR 2027
Tapas Tripura and Souvik Chakraborty. Wavelet neural operator for solving parametric partial differential equations in computational mechanics problems. Computer Methods in Applied Mechanics and Engineering, 404:115783, 2023. 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, 2021a. Sifan Wang, Hanwen Wang, and Paris Perdikaris. On the eigenvector bias of fourier feature networks: From regression to solving multi-scale pdes with physics-informed neural networks. Computer Methods in Applied Mechanics and Engineering, 384:113938, 2021b. Sifan Wang, Hanwen Wang, and Paris Perdikaris. Learning the solution operator of parametric partial differential equations with physics-informed deeponets. Science advances, 7(40):eabi8605, 2021c. Sifan Wang, Xinling Yu, and Paris Perdikaris. When and why pinns fail to train: A neural tangent kernel perspective. Journal of Computational Physics, 449:110768, 2022. Sifan Wang, Shyam Sankaran, Hanwen Wang, and Paris Perdikaris. An expert’s guide to training physics-informed neural networks. arXiv preprint arXiv:2308.08468, 2023a. 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, 2024a. Sifan Wang, Shyam Sankaran, and Paris Perdikaris. Respecting causality for training physicsinformed neural networks. Computer Methods in Applied Mechanics and Engineering, 421: 116813, 2024b. Wei Wang, Tang Paai Wong, Haihui Ruan, and Somdatta Goswami. Causality-respecting adaptive refinement for pinns: enabling precise interface evolution in phase field modeling. Machine Learning for Computational Science and Engineering, 2(1):10, 2026. Zhicheng Wang, Xuhui Meng, Xiaomo Jiang, Hui Xiang, and George Em Karniadakis. Solution multiplicity and effects of data and eddy viscosity on navier-stokes solutions inferred by physicsinformed neural networks. arXiv preprint arXiv:2309.06010, 2023b. Karl Weiss, Taghi M Khoshgoftaar, and DingDing Wang. A survey of transfer learning. Journal of Big data, 3(1):9, 2016. Chenxi Wu, Min Zhu, Qinyang Tan, Yadhu Kartha, and Lu Lu. A comprehensive study of nonadaptive and residual-based adaptive sampling for physics-informed neural networks. Computer Methods in Applied Mechanics and Engineering, 403:115671, 2023. Jeremy Yu, Lu Lu, Xuhui Meng, and George Em Karniadakis. Gradient-enhanced physics-informed neural networks for forward and inverse pde problems. Computer Methods in Applied Mechanics and Engineering, 393:114823, 2022.
A
P ER - PROBLEM INSTANTIATIONS OF THE SHARED ARCHITECTURE
A.1
D OMAINS AND THE ISOPARAMETRIC MAP
The two-dimensional problems are posed on ΩH = (0, 1)×(0, H); the one-dimensional problem on (0, L). A reference domain (ξ, η) ∈ [0, 1]2 is mapped to physical coordinates by (x, y) = (ξ, Hη), yielding ∂x = ∂ξ and ∂y = H −1 ∂η . Analogously, the one-dimensional reference coordinate ξ ∈ [0, 1] is mapped to the physical domain by x = Lξ, yielding ∂x = L−1 ∂ξ . The trunk basis is defined exclusively on the reference domain, such that geometry enters solely through algebraic chainrule weights. Reference-space partial derivatives of the frozen basis remain geometry-independent; consequently, assembling the physical operator at a new geometry (i.e. a new L or H) reduces to algebraic scaling. 12
Published as a conference paper at ICLR 2027
A.2
B URGERS : PDE, ANSATZ , RESIDUAL
On (0, L) with u(0) = +A, u(L) = −A (we use A = 1 in all experiments), u ∂x u − ν ∂xx u = 0, (7) P −1 −2 with ansatz u(ξ; p) = A(1−2ξ)+4ξ(1−ξ) k bk (p)φk (ξ) and residual r = L u∂ξ u−νL ∂ξξ u. In the notation of equation 1, the boundary lift is g(ξ; p) = A(1 − 2ξ) and the multiplier is M(ξ) = 4ξ(1 − ξ), which vanishes to first order at ξ ∈ {0, 1} (Table 1, Design choice ii). A.3
A LLEN –C AHN : PDE, ANSATZ , RESIDUAL , BIFURCATION
On ΩH with u = 0 on ∂ΩH ,
∇2 u + λ(u − u3 ) = 0. (8) 2 2 2 The critical rate is λ1 (H) = π (1 + 1/H ), representing the first Dirichlet eigenvalue of −∇ on the domain. Below λ1 (H), the only solution is the trivial state u ≡ 0; above it, a pair ±u⋆ bifurcates √ ⋆ along the principal (nodeless, sign-definite) eigendirection of −∇2 , scaling as ∥u√ ∥∞ ∼ λ − λ 1 near the onset and saturating to ∥u⋆ ∥∞ → 1 with boundary layers of width ∼ 1/ λ further above the critical rate. Because the seed lift ℓ0 > 0 (defined below) biases the ansatz toward the positive member of this pair, we take u⋆ > 0, the principal positive branch, as the target solution throughout. The base training domain λ ∈ [35, 80], H ∈ [1, 2] is restricted to the comfortably supercritical regime. The ansatz employs a fixed non-trivial seed lift ℓ0 (ξ, η) = A0 sin(πξ) sin(πη), with A0 = 0.3 used in all experiments: X u(ξ, η; p) = A0 sin(πξ) sin(πη) + 16ξ(1 − ξ)η(1 − η) bk (p)φk (ξ, η). (9) k
In the notation of equation 1, g(ξ, η; p) = ℓ0 (ξ, η) (constant in p) and M(ξ, η) = 16ξ(1−ξ)η(1−η), which vanishes to first order on ∂ΩH (Table 1, Design choice ii). This formulation explicitly breaks the u ≡ 0 degeneracy at initialization; setting b ≡ 0 yields a positive interior profile rather than a trivial zero field. The seed function vanishes on ∂ΩH and satisfies ∇2 ℓ0 = −π 2 (1 + H −2 )ℓ0 = −λ1 (H) ℓ0 : it is an eigenfunction of the linear Laplacian, with eigenvalue exactly equal to the critical bifurcation rate, but it is not an eigenfunction of the nonlinear residual equation 8, so it is never itself a spurious fixed point of training. A.4
C AVITY: PDE, ANSATZ , RESIDUAL
The streamfunction ψ satisfies the vorticity-transport equation: −u ∂x (∇2 ψ) − v ∂y (∇2 ψ) + Re−1 ∇4 ψ = 0,
u = ∂y ψ, v = −∂x ψ,
(10)
enforcing velocity Dirichlet Pconditions on all walls alongside a prescribed lid velocity ulid . The ansatz ψ = g(x) + M(x) k bk φk (x) employs a lift function g(ξ, η) = Hulid (ξ)(η 3 − η 2 ) and a multiplier M = (16ξ(1 − ξ)η(1 − η))2 . Because M vanishes to the second order on all boundaries (Table 1, Design choice ii), both ψ and its normal derivative ∂n ψ exactly satisfy the boundary conditions independent of the network coefficients. Regularized lid profile. A uniformly driven lid (u = 1) introduces a velocity discontinuity at the two top corners where the moving boundary intersects the stationary side walls. This corner singularity generates unbounded vorticity, which frequently destabilizes high-Re residual training (Wang et al., 2023b). To circumvent this, we regularize the boundary condition. The lid speed is maintained at exactly 1 over the interior of the top wall and transitions smoothly to 0 at each corner across a margin of width δ = 0.12: S(t) = 0, ulid (ξ) = S δξ S 1−ξ , δ
t < 0, 5
4
3
S(t) = 6t − 15t + 10t , S(t) = 1,
0 ≤ t ≤ 1, t > 1.
(11)
where S is Perlin’s quintic smootherstep function (Perlin, 2002). This formulation ensures matching values, slopes, and curvatures (S(0)=0, S(1)=1, S ′ =S ′′ =0 at both domain ends), guaranteeing that 13
Published as a conference paper at ICLR 2027
the transition into the flat core introduces no gradient discontinuities up to the second derivative. With δ = 0.12, the lid remains exactly 1 over the central 76% of the span, modifying the physical driving force only within the narrow corner margins. This identical ulid profile is applied to the lift g, the finite-difference reference solver, and the Part 3 student’s boundary loss, ensuring the teacher and student solve the exact same regularized boundary-value problem. Part 3 primitive-variable formulation. Unlike Parts 1–2, which operate on the streamfunction ψ, the Part 3 cavity student is a primitive-variable PINN uϕ (x) = (uϕ , vϕ , pϕ ) trained directly on the steady incompressible Navier–Stokes equations: uϕ ∂x uϕ + vϕ ∂y uϕ = −∂x pϕ + Re−1 ∂xx uϕ + ∂yy uϕ , (12) −1 uϕ ∂x vϕ + vϕ ∂y vϕ = −∂y pϕ + Re ∂xx vϕ + ∂yy vϕ , (13) ∂x uϕ + ∂y vϕ = 0.
(14)
Dirichlet velocity conditions are enforced softly: uϕ = vϕ = 0 on the left, right, and bottom walls, and uϕ = ulid (x), vϕ = 0 on the top wall, using the identical regularized lid profile equation 11 applied throughout the pipeline. Because pϕ is determined only up to an additive constant, we impose a soft gauge penalty pinning its value at the domain center to zero, pϕ (0.5, H/2) ≈ 0. The complete Part 3 objective is L(ϕ; t) = Lmom,x + Lmom,y + Lcont + λBC LBC + λgauge pϕ (0.5, H/2)2 + w(t) Ldistill (t), (15) where Lmom,x , Lmom,y , and Lcont are the mean-squared residuals of equation 12–equation 14 over interior collocation points, LBC is the mean-squared boundary-velocity violation, and λBC , λgauge are the momentum-normalized BC and gauge weights reported in Table 7 (10.0 and 10.0 respectively for the main Re-extrapolation target; the momentum terms carry unit weight). The distillation term is applied only to velocity: Ldistill (t) =
2
uϕ − uteacher Γ(x) +
2
vϕ − vteacher Γ(x) ,
(16)
where (uteacher , vteacher ) = (∂y ψteacher , −∂x ψteacher ) are the induced velocities of the Part 2 streamfunction field, obtained by automatic differentiation of the frozen Part 2 network. No distillation target is imposed on pϕ : the streamfunction teacher carries no pressure information, so the pressure field is constrained solely by the momentum residuals and the gauge penalty. This primitive-variable formulation for Part 3 is used identically for both the Re-extrapolation (§4.3) and H-extrapolation (Appendix E) cavity experiments; the Burgers and Allen–Cahn students, by contrast, retain the scalar-field loss equation 6 directly on uϕ , since those problems have no pressure or continuity constraint. Interior collocation points for the momentum and continuity residuals are not sampled uniformly: the wall-layer thickness scales as δwall ∼ Re−1/2 , comparable to or smaller than the spacing of a uniform sample, which otherwise lets the student smooth away the near-wall velocity extrema. We instead draw 40% of interior points from thin near-wall bands (width 4Re−1/2 , biased toward each wall) and the remaining 60% uniformly, applied identically to the baseline and distilled arms so it does not bias their comparison. Distillation query points are sampled uniformly rather than wall-clustered, since the teacher’s velocity error is largest near the walls and clustering there would concentrate the distillation pull on its least reliable signal. A.5
S PATIAL BANDWIDTH AND SPECTRAL BIAS
Standard coordinate-based networks exhibit spectral bias, so a random Fourier embedding is used to shift the burden of representing high-frequency features onto the analytical derivatives of the basis sinusoids. The required bandwidth scales with the PDE’s order: a frequency-k component entering an n-th-order residual is amplified by a factor of ∼ k n , bounded by k 2 for the second-order Burgers and Allen–Cahn systems but reaching k 4 for the fourth-order cavity system, further exacerbated by a factor of H −4 at high aspect ratios. The cavity’s trunk bandwidth is accordingly the most constrained of the three, using an anisotropic σ = [σξ , ση ]; Allen–Cahn and Burgers use a milder, isotropic bandwidth. Exact values for all three are given in Tables 5–7. 14
Published as a conference paper at ICLR 2027
A.6
W HERE THE THREE INSTANTIATIONS DIFFER
Beyond the three primary design choices outlined in Table 1, the implementations differ across four specific algorithmic parameters (detailed in Appendices B through C): (a) The Burgers base operator employs a pure autograd-derived residual. Conversely, Allen–Cahn and the cavity utilize a weighted finite-difference stencil residual, which significantly reduces computational overhead for the biharmonic cavity operator and the batched Allen–Cahn Laplacian. (b) The Allen–Cahn framework incorporates a soft amplitude floor into the Part 1 loss, penalizing batch samples whose fields collapse toward u≡0. This serves as a complementary regularization to the seed lift ℓ0 , mitigating the risk of the network encountering the competing trivial minimizer during early training iterations. (c) The extrapolation sequence detailed in §3.3 requires one stage for the Burgers equation and the cavity’s aspect-ratio variant (Appendix E), two stages for the Allen–Cahn problem, and four stages for the principal cavity Re-extrapolation. (d) Allen–Cahn and the primary Re-extrapolated cavity deploy the ungated, piecewise-linear distillation schedule. Burgers and the cavity H-extrapolation use the physics-adaptive, reliability-gated schedule (§3.4).
B
D ISTILLATION SCHEDULE AND RELIABILITY GATE
The Part 3 loss equation 6 weights the distillation term using two distinct components: a temporal scalar schedule w(t) and an optional spatial per-point gate Γ(x). The Allen–Cahn system and the cavity’s Re-extrapolation (main text) utilize an ungated, piecewise-linear temporal schedule (Γ ≡ 1, the identity); the Burgers equation and the cavity’s H-extrapolation (Appendix E) utilize a physicsadaptive temporal schedule coupled with the spatial reliability gate defined below. Piecewise-linear weight. Over a fraction fmin of the total training budget T , w(t) decays from an initial value w0 to an intermediate value wmin , and subsequently decays to wend =0 over the remainder of training: t t ≤ fmin T, w0 + (wmin − w0 ) f T , min w(t) = (17) t − fmin T wmin + (wend − wmin ) , t > fmin T. (1 − fmin )T For the cavity’s Re-extrapolation (Re=3200), the parameters are set to (w0 , fmin , wend ) = (0.7, 0.7, 0) with wmin =0. This constitutes a linear decay to zero over the first 70% of training, followed by exclusive residual minimization. Truncating w(t) to exactly zero during the final training window allows the student network to refine the solution independent of the teacher’s structural errors. Physics-adaptive weight. While the piecewise-linear schedule is sufficient for Allen–Cahn and the cavity’s Re-extrapolation, we use a physics-adaptive schedule for the Burgers equation and the cavity’s H-extrapolation to demonstrate the framework’s flexibility in accommodating alternative distillation strategies. After an initial warm-up period of Twarm epochs, during which w(t) increases linearly to wmax αw , the weight dynamically tracks the convergence of the student network: w = wmax α · min(1, r̄ϕ /r̄ϕref ). Here, r̄ϕ represents an exponential moving average (with a decay factor of 0.5) of the student’s PDE loss, and r̄ϕref is its reference value at the conclusion of the warmup phase. The weight w(t) is strictly zeroed once training surpasses a specified cutoff fraction of T . For the Burgers equation, the parameters are configured as (αw , wmax , Twarm , cutoff frac.) = (0.8, 1.0, 200, 0.5); the cavity’s H-extrapolation uses the same (α, wmax , Twarm , cutoff frac.) values. This formulation dynamically attenuates the teacher’s influence in response to the student’s improving physics residual, rather than relying on a predetermined epoch schedule. Spatial reliability gate Γ(x). When enabled (as in the Burgers and cavity H-extrapolation configurations), the spatial weight Γ acts as a per-collocation-point factor that attenuates the distillation constraint in regions where the student’s PDE residual rϕ (x) is large: Γ(x) = 1 − tanh(τ rϕ (x)) ∨ Γfloor , Γfloor = 0.1. (18) This value is clamped to a minimum floor Γfloor to ensure the teacher’s prior is never entirely deactivated. In regions where the student adequately satisfies the governing physics (rϕ is small), the gate 15
Published as a conference paper at ICLR 2027
Table 5: Burgers hyperparameters (§4.1). Part 1 — base operator [0.10, 0.20] × [1.0, 2.0] {0.10, 0.15, 0.20} × {1.0, 1.5, 2.0} (9) 128 × 5 tanh / 128 × 4 tanh / 96 64 / 6.0 2 × 104 / 10−3 → 10−6 256 / 16 / 1.0
Base box (ν, L) Anchors Trunk / branch / basis K Fourier features m / σ Operator steps / LR (→ min) Collocation pts / param. batch / grad clip Part 2 — continuation
(ν, L) = (0.05, 4.0) 0.63 [10−3 , 108 ] 1.2 × 104 / 8 × 10−4 / 8
Target Bounded-correction cap α β range Steps / LR / batch Part 3 — guarded distillation
128 × 4 tanh / 32, 8.0 4000 / 5 × 10−4 → 10−6 1000 / 256 / 1000 10.0 / 1.0 physics-adaptive + gate 2.0/0.1; (0.8, 1.0, 200, 0.5) 5 / 1729
Student / Fourier m, σ Epochs / LR (→ min) Collocation / BC / distil pts BC weight / grad clip Distillation schedule Gate τ /floor; (αw , wmax , Twarm , cutoff) Seeds (n / start)
approaches 1 and the teacher’s field is strictly enforced. Conversely, where the residual is large, the gate decays toward the floor value. For the Burgers equation, we utilize a lowered temperature of τ = 2.0 as the residual magnitude |rϕ | peaks inside the internal transition layer. Because this layer contains the critical structural prior provided by the 1D teacher, a high-temperature gate would inappropriately deactivate distillation in this region and restrict knowledge transfer exclusively to the trivial outer domains. Lowering τ ensures the gate remains sufficiently open at moderate residual values, a behavior we empirically verified by monitoring the gate-open fraction during training. The cavity’s H-extrapolation instead uses τ = 20.0 (with the same floor Γfloor = 0.1), reflecting the different residual scale of the two-dimensional momentum system.
C
H YPERPARAMETERS AND TRAINING DETAILS
Table 5 gives the Burgers hyperparameters (§4.1), Table 6 the Allen–Cahn hyperparameters (§4.2), and Table 7 the cavity hyperparameters for both the main-text Re-extrapolation (§4.3) and the Appendix E H-extrapolation, which are independently trained base operators. All values are taken directly from the released configuration; code is released alongside the paper for further implementation details.
D
R EFERENCE - SOLVER CONFIGURATIONS
Reference fields are used only for post-hoc scoring and, for the cavity, for forming the eight maintext Part 1 anchors (twelve for the separate Appendix E base operator); they never enter any training loss. For Burgers the reference is the Cole–Hopf closed form, evaluated exactly. For Allen–Cahn the reference is a Newton solve of the residual on a five-point Laplacian stencil. For the cavity the reference is a regularized-lid, streamfunction–vorticity finite-difference solver on a uniform Nx =129 grid (isoparametrically stretched to Ny = ⌈(Nx − 1)H⌉ + 1 for H ̸= 1), iterated to a steady state; it uses the same regularized lid profile ulid as the operator (Appendix A.4). Each cavity anchor solve is gate-checked against the Marchi et al. (Marchi et al., 2021) high-resolution benchmark where tabulated; the few-percent offset expected from the regularized-versus-singular lid is accounted for in that check. 16
Published as a conference paper at ICLR 2027
Table 6: Allen–Cahn hyperparameters (§4.2). Part 1 — base operator [35, 80] × [1.0, 2.0] {35, 50, 65, 80} × {1.0, 1.5, 2.0} (12), GNpolished 128 × 4 tanh / 128 × 4 tanh / 128 32 / [3.0, 3.0] 3 × 104 / 10−3 → 10−6 48 / 8 rms(u) ≥ 0.25 / 100 60 / 10−10 / 10−4 / 1.0
Base box (λ, H) Anchors Trunk / branch / basis K Fourier features m / σ Operator steps / LR (→ min) FD residual grid / param. batch Anti-collapse amplitude floor / weight GN iters / tol / damping / grad clip Part 2 — continuation
two-rung ladder: λ=18 → λ=12.5 (both at H=3) 0.63 [1.0, 2000] 2 1.2 × 104 / 5 × 10−4 / 8
Targets Bounded-correction cap α β range Poly. degree Steps / LR / batch Part 3 — guarded distillation
64 × 2 tanh / 32, 2.0 104 / 5 × 10−4 → 10−7 2000 / 800 / 2000 10.0 / 1.0 piecewise-linear (1.0, 0.5, 0.1, 0.0) 5 / 1729
Student / Fourier m, σ Epochs / LR (→ min) Collocation / BC / distil pts BC weight / grad clip Distillation schedule (w0 , fmin , wmin , wend ) Seeds (n / start)
E
A DDITIONAL RESULTS : ASPECT- RATIO EXTRAPOLATION
Section 4.3 reports the cavity’s Part 3 result at the Re-extrapolated target (Re = 3200, H = 1), reached from a base operator trained on Re ∈ [100, 800] at fixed H=1 with 8 anchors. Here we report a separate experiment stressing the complementary axis: a base operator trained on the full (Re, H) box with 12 anchors (Table 7), single-shot extrapolated from H ∈ [1, 2] out to H = 3 at fixed Re = 400, followed by Part 3 distillation of a fresh student PINN at that target. Every other mechanism is exactly as in the main pipeline (§3.3–§3.4). This example shows that the framework’s benefit from distillation is not specific to the Re-extrapolation axis reported in the main text; the same guarded-distillation mechanism accelerates convergence when the extrapolation instead stresses geometry. Figure 6 compares the baseline and distilled students against the finite-difference reference across training progress (25/50/75/100% of the training budget) for both velocity components. Figure 7 reports the corresponding total loss, accuracy against the FD reference, and PDE residual histories. As in the main-text Re = 3200 case, the distilled student converges faster and reaches a lower final residual and higher final accuracy than the baseline trained under an identical architecture and step budget.
17
Published as a conference paper at ICLR 2027
Table 7: Cavity hyperparameters: main-text Re-extrapolation (Re=3200, H=1, §4.3) vs. the Appendix E H-extrapolation (H=3, Re=400). Main: Re=3200, H=1
App. E: Re=400, H=3
Part 1 — base operator Base box (Re, H) Anchors Trunk / branch / basis K
[100, 800] × {1.0} 8 128×8 tanh / 128×8 tanh / 192
Fourier features m / σ Operator steps / LR (→ min) FD residual grid / param. batch Ridge µ / poly degree Grad clip
64 / [4.0, 7.0] 2 × 105 / 10−3 → 10−6 160 / 4 10−2 / 2 1.0
[100, 400] × [1.0, 2.0] 12 256 × 6 tanh / 256 × 4 tanh / 192 48 / [4.0, 7.0] 3 × 105 / 6 × 10−4 → 10−6 64 / 8 10−2 / 3 1.0
Bounded-correction cap α β calibration clamp [βmin , βmax ] β within-step decay (start → end) Steps / LR / batch
4-rung Re ladder: 1400→2000→2500→3200 0.63 [20, 2000] β0 → β0 /20 2.4 × 104 / 5 × 10−4 / 8
single-shot: (Re, H)=(400, 3) 0.63 [20, 2000] β0 → β0 /20 2.4 × 104 / 5 × 10−4 / 8
Part 3 — guarded distillation Student / Fourier m, σ Epochs / LR (→ min) Collocation / BC / distil pts Momentum / BC / gauge weight Grad clip Distillation schedule Seeds (n / start)
256 × 5 tanh / 64, 6.0 2.5 × 105 / 10−3 → 10−8 6000 / 4000 / 4000 1.0 / 10 / 10.0 1.0 piecewise-linear 5 / 1729
128 × 3 tanh / 32, 4.0 5 × 104 / 5 × 10−4 → 10−6 4000 / 800 / 4000 1.0 / 10 / 1.0 1.0 physics-adaptive + gate single run
Part 2 — continuation Targets
Figure 6: Lid-driven cavity, H-extrapolation example at Re = 400, H = 3 (base box H ∈ [1, 2]). Predicted horizontal (u, top row) and vertical (v, bottom row) velocity fields at 25/50/75/100% of training for the baseline (left block) and distilled (right block) students, against the finite-difference reference (rightmost column).
18
Published as a conference paper at ICLR 2027
Figure 7: Lid-driven cavity, H-extrapolation example at Re = 400, H = 3.
Figure 8: Cavity Part-2 removal ablation at Re = 3200, H = 1. Solid lines and shaded regions represent µ ± σ over n = 5 seeds.
19