Spectral bias in physics-informed and operator learning: Analysis and mitigation guidelines
Authors: Siavash Khodakarami, Vivek Oommen, Nazanin Ahmadi Daryakenari, Maxim Beekenkamp, George Em Karniadakis
Organizations: Division of Applied Mathematics, Brown University, Providence, RI 02912, USA · School of Engineering, Brown University, Providence, RI 02912, USA
Solving partial differential equations (PDEs) by neural networks as well as Kolmogorov-Arnold Networks (KANs), including physics-informed neural networks (PINNs), physics-informed KANs (PIKANs), and neural operators, are known to exhibit spectral bias, whereby low-frequency components of the solution are learned significantly faster than high-frequency modes. While spectral bias is often treated as an intrinsic representational limitation of neural architectures, its interaction with optimization dynamics and physics-based loss formulations remains poorly understood. In this work, we provide a systematic investigation of spectral bias in physics-informed and operator learning frameworks, with emphasis on the coupled roles of network architecture, activation functions, loss design, and optimization strategy. We quantify spectral bias through frequency-resolved error metrics, Barron-norm diagnostics, and higher-order statistical moments, enabling a unified analysis across elliptic, hyperbolic, and dispersive PDEs. Through diverse benchmark problems, including the Korteweg-de Vries, wave and steady-state diffusion-reaction equations, turbulent flow reconstruction, and earthquake dynamics, we demonstrate that spectral bias is not simply representational but fundamentally dynamical. In particular, second-order optimization methods substantially alter the spectral learning order, enabling earlier and more accurate recovery of high-frequency modes for all PDE types. For neural operators, we further show that spectral bias is dependent on the neural operator architecture and can also be effectively mitigated through spectral-aware loss formulations without increasing the inference cost.
Figures & tables
Architecture
Rel. L2 Error
Barron Norm Rel. Error
log(eF,p=0)
log(eF,p=2)
log(eF,p=4)
Case 1: Discontinuous Function
MLP–Tanh
4.96×10−4
2.35×10−4
−3.72
1.31
6.38
MLP–SIREN
4.68×10−4
6.48×10−4
−3.78
0.01
4.78
cKAN
3.67×10−2
1.62×10−2
0.01
5.07
10.17
Tanh-cKAN
5.19×10−4
3.28×10−4
−3.68
−0.04
4.86
Case 2: Smooth Multi-Frequency Function
Table 1: Data-driven spectral benchmark errors. Comparison of neural network architectures for Case 1 (discontinuous function) and Case 2 (smooth multi-frequency function). Performance is evaluated using physical- and frequency-domain error metrics.
Figure 1: Case 1: Training-time spectral evolution for the discontinuous benchmark. Comparison of MLP–Tanh, MLP–SIREN, cKAN, and Tanh-cKAN on the piecewise-discontinuous target function. Columns correspond to different training epochs, illustrating the evolution of the predicted signal in physical space (top row of each block) and the corresponding Fourier amplitude spectrum (bottom row of each block). The ground-truth solution is shown in black, while network predictions are shown in red. The discontinuity induces broad-band high-frequency content in the frequency domain, highlighting differences in how each architecture and optimization strategy capture sharp transitions and recover high-frequency modes over training.
Figure 2: Case 2: Training-time spectral evolution for the smooth multi-frequency benchmark. Comparison of MLP–Tanh, MLP–SIREN, cKAN, and Tanh-cKAN on the multi-scale multi-frequency target function. Columns correspond to different training epochs to illustrate the evolution of the predicted signal in physical space (top row of each block) and the corresponding Fourier amplitude spectrum (bottom row of each block). The ground-truth solution is shown in black, while network predictions are shown in red. Although the target function is smooth, the presence of multiple well-separated frequency components makes this benchmark sensitive to spectral bias to reveal differences in how each architecture captures and balances low- and high-frequency modes during training.
Figure 3: KdV equation: PINN training history with different optimizers (Adam, L-BFGS, SOAP, and SS-Broyden) and activation functions (Tanh, and SIREN) for KdV equation. The sub-figure shows the zoomed-in training history for the SS-Broyden optimizer.
Figure 4: KdV equation: Analysis of the predictions. (a, b) First four moments (mean, variance, skewness, and kurtosis) of the PINN predictions with (a) Adam optimizer and Tanh activation, and (b) Soap or SS-Broyden optimizer with Tanh activations. (c) Time-averaged absolute error at each moment for PINN predictions with different optimizers and activations. (d) Barron norm at each time-step of PINN predictions with different optimizers.
Figure 5: KdV equation: Effect of optimizer and activation function. First row: Ground truth (GT) solution, gradient of the GT solution, and Laplacian of the GT solution. The other rows in first, second and third columns show the errors in PINN solution, errors in gradient magnitude of the solution, and errors in Laplacian of the solution, respectively. The rows are organized from the most to the least accurate results (top to bottom). The axis ranges ( x∈[−1,1] , t∈[0,1] ) shown in top left subplot is applicable to all.
Figure 6: KdV equation: Spectral bias dynamics during the training. (a-d) Errors of first four statistical moments during the training with different optimizers and activation functions. (e) Error in the frequency domain during the training. The small red arrows show the approximate slope of learning the high-frequencies (shown in third and fourth moments). The dotted black lines in the four moments shows the iteration at which the Tanh and SIREN networks trained with SOAP intersect. Note how the intersection happens earlier at higher moments. The legends in (a) are applicable to all.
Rel. L2 Error
Barron Norm Rel. L2 Error
Log(eF,p=0)
Log(eF,p=2)
Log(eF,p=4)
Convergence Time [hr.]
Time to 1% Error
PINNs
Adam (Tanh)
2.47×10−1
1.58×10−1
-1.52
0.77
5.34
0.056
NA
Adam+LBFGS (Tanh)
8.12×10−2
4.14×10−2
-2.48
-0.15
4.35
0.056+0.21
NA
Adam (SIREN)
6.35×10−2
1.06×10−2
-2.69
-0.19
4.37
0.60
NA
Adam+LBFGS (SIREN)
1.76×10−3
1.28×10−3
-5.48
-2.94
1.37
0.60+0.045
0.29
SOAP (Tanh)
7.05×10−4
5.62×10−4
-6.59
-3.91
0.39
0.40
0.033
Table 2: KdV equation prediction errors. Comparison of the performance of PINNs with different activation functions and PIKANs with different univariate functions trained with different optimizers. All models except those trained with SS-Broyden have ∼ 40,000 parameters. The models trained with SS-Broyden have ∼ 12,900 parameters.
Figure 7: Wave equation results with different optimizers and activation functions (a-b). Relative L2 error and fourth moment error of each of the wave problems trained with different optimizers and activation functions. (c) Skewness and Kurtosis of PINN prediction with SOAP optimizer and SIREN activation. (d) Skewness and Kurtosis of PINN prediction with SS-Broyden optimizer and SIREN activation.
Method
Rel. L2 Error
Barron Norm Rel. L2 Error
log(eF,p=0)
log(eF,p=2)
log(eF,p=4)
PINN
Adam (Tanh)
1.00
0.99
-1.65
2.23
6.80
Adam (SIREN)
1.00
0.95
-1.65
2.23
6.80
SOAP (Tanh)
4.03×10−2
2.12×10−2
-4.44
-0.69
3.92
SOAP (SIREN)
4.34×10−3
2.54×10−3
-6.37
-2.27
2.28
SS-Broyden (Tanh)
3.05×10−5
2.37×10−4
-10.67
-5.22
0.56
Table 3: Wave equation case (iv) prediction errors. Comparison of PINN and PIKAN performance for the high-frequency wave problem. Results for simpler cases (i)–(iii) are reported in the appendix. All models except those trained with SS-Broyden have ∼ 40,000 parameters. The models trained with SS-Broyden have ∼ 12,900 parameters.
Method
Rel. L2 Error
Barron Norm Rel. L2 Error
log(eF,p=0)
log(eF,p=2)
log(eF,p=4)
#params K
Convergence Time (hr)
PINN
Adam (Tanh)
1.41×10−1
1.21×10−4
1.69
-2.39
1.99
14.6
0.18
Adam (SIREN)
3.14×10−3
2.66×10−3
-1.60
0.23
3.44
14.6
0.11
SOAP (Tanh)
6.31×10−2
3.17×10−5
0.99
-3.26
-0.56
14.6
0.38
SOAP (SIREN)
1.34×10−2
3.93×10−6
-0.35
-3.99
-1.39
14.6
0.29
SS-Broyden (Tanh)
1.27×10−6
7.08×10−10
-8.39
-12.63
-10.11
14.6
0.06
Table 4: Steady state diffusion-reaction equation prediction errors. Comparison of PINN and PIKAN performance for the diffusion-reaction problem.
Figure 8: Steady-state diffusion–reaction equation: Effect of optimizer and activation function. Comparison of ground-truth field, gradient magnitude, and Laplacian (top row) with corresponding prediction errors (bottom row) for representative PINN and PIKAN configurations. The rows are organized from the most to the least accurate results (top to bottom). The error maps illustrate qualitative differences in spatial and derivative accuracy that are not fully captured by scalar norms alone. These visualizations complement the quantitative metrics reported in Table 4 .
Figure 9: Turbulent jet: Spectral bias in neural operators. Columns show (from left to right) schlieren images of an impinging jet, the corresponding energy spectrum E(k), spatial gradients, and the Laplacian. Rows compare the ground truth against predictions from different neural operator architectures (with and without BSP training). We observe that using Binned Spectral Power (BSP) loss helps mitigate spectral bias, improving agreement with the true energy spectrum across wavenumbers. All the results are based on SOAP optimizer.
Model
Field Error
Energy- Spectrum Error
Barron- Norm Error
# params (M)
Inference time (s)
DeepONet
0.0529
0.1566
0.4171
3.8
0.00086
DeepONet (BSP)
0.0535
0.1675
0.4083
3.8
0.00086
DeepOKAN
0.0538
0.1289
0.3659
3.0
0.00991
DeepOKAN (BSP)
0.0547
0.0416
0.1662
3.0
0.00991
FNO
0.0517
0.0865
0.2641
3.4
0.00433
FNO (BSP)
0.0556
0.0256
0.0828
3.4
0.00433
Table 5: Comparison of neural operator errors and costs for the impinging jet problem.
Model
Field Error
Log Spectral Error
Barron Norm Error
Baseline
BSP
Baseline
BSP
Baseline
BSP
DeepONet (SIREN)
0.0008
0.0009
0.1257
0.1366
0.0019
0.0027
DeepONet (Tanh)
−
0.0006
−
0.1103
−
0.0012
DeepOKAN
0.0022
0.0044
0.1230
0.1750
0.0029
0.0077
FNO
0.0038
0.0041
0.1485
0.0819
0.0045
0.0030
CNO
0.0133
0.0145
0.1323
0.0758
0.0040
0.0043
Table 6: Architecture and loss comparison for earthquake structural response prediction. SOAP optimizer, causal training for DeepONet/DeepOKAN.
Figure 10: Representative earthquake response predictions illustrating key design choices. Each row shows the time-domain response (left) and energy spectrum (right) for the same test sample. (a) DeepONet (SIREN) achieves the lowest field error (NRMSE =0.0008 ; Table 6 ), followed by DeepOKAN ( 0.0022 ), FNO ( 0.0038 ), and CNO ( 0.0133 ). (b) BSP improves FNO spectral fidelity, reducing log spectral error from 0.149 to 0.082 ( 1.8× ), visible in the high-frequency alignment. (c) SOAP enables DeepOKAN convergence (NRMSE =0.0022 ), whereas Adam fails across all seeds. (d) SIREN activation enables DeepONet convergence without BSP (NRMSE =0.0008 ), while Tanh fails to converge under the same conditions. Black dashed lines denote ground truth; extended comparisons are in Figures E.4 – E.6 .
Appendix figures & tables11 assets
Supplementary material from the paper’s appendix.
Appendix
Figure B.1: PINN with adaptive Tanh activation function training history for KdV equation. Relative L2 error during training with adaptive activation functions including slope recovery loss term with different weights ( wsr ), trained with Adam. Note that by increasing the weight for the slope recovery term, the errors starts dropping down earlier.
PINN
Rel. L2 Error
Barron Norm Rel. L2 Error
log(eF,p=0)
log(eF,p=2)
log(eF,p=4)
Case (i)
Adam (Tanh)
5.08×10−3
2.94×10−3
-5.19
-1.83
2.27
Adam (SIREN)
6.11×10−3
4.45×10−3
-5.03
-1.58
2.54
SOAP (Tanh)
1.19×10−3
2.86×10−4
-6.44
-4.28
0.80
SOAP (SIREN, w0=30 )
1.59×10−3
2.19×10−4
-5.11
-4.24
0.49
SOAP (SIREN, w0=5 )
8.45×10−4
1.81×10−4
-6.59
-4.20
0.50
Appendix
Table C.1: Wave equation prediction errors. Comparison of PINN performance with different activation functions and optimizers for cases (i)-(iii). The results for case (iv) are shown in Table 3 of the main text. Note that for case (i) with no multi-scale features or high frequencies, using SIREN with w0=5 is better than using SIREN with w0=30 .
Figure D.2: Optimizer effects on spectral bias in neural operators. Columns show (from left to right) Schlieren fields, the isotropic energy spectrum E(k) , spatial gradients, and the Laplacian. Rows compare the ground truth against DeepONet and CNO trained with Adam and SOAP. SOAP yields improved agreement with the reference spectrum at high wavenumbers and recovers sharper derivative-based features, indicating reduced attenuation of fine scales relative to Adam.
Figure D.3: SOAP vs Adam for earthquake response prediction. Time response and energy spectrum for a representative test sample comparing Adam (left pair) and SOAP (right pair). Adam DeepOKAN fails to converge (near-flat prediction and flat spectral tail), while SOAP enables convergence for all architectures with improved spectral accuracy at high frequencies.
Model
Field Error
Log Spectral Error
Barron Norm Error
Accel. Error
Adam
SOAP
Adam
SOAP
Adam
SOAP
Adam
SOAP
DeepONet (SIREN)
0.0015
0.0008
0.0974
0.1257
0.0025
0.0019
0.0381
0.0506
DeepONet (SIREN, BSP)
0.0013
0.0009
0.1026
0.1366
0.0020
0.0027
0.0314
0.0619
DeepONet (Tanh)
−
−
−
−
−
−
−
−
DeepONet (Tanh, BSP)
0.0013
0.0006
0.1090
0.1103
0.0021
0.0012
0.0373
0.0303
DeepOKAN
−
0.0022
−
0.1230
−
0.0029
−
0.0584
Appendix
Table E.2: Optimizer, activation, and loss comparison for earthquake response prediction (causal training). − Failed to converge.
Figure E.4: Effect of BSP on spectral fidelity for earthquake response prediction. Time response (columns 1, 3) and energy spectrum (columns 2, 4) for a representative test sample comparing baseline (left pair) and BSP (right pair) training across all four architectures. Note the difference in scale of energy values (y-axis) when comparing these results to the impinging-jet problem.
Figure E.5: Activation function effects on DeepONet. Time response and energy spectrum comparing Tanh (left pair) and SIREN (right pair) activations for both baseline (top row) and Log-BSP (bottom row) training. Tanh baseline fails to capture the structural response, but BSP enables convergence to the best overall accuracy. SIREN converges reliably without BSP but shows diminishing returns when BSP is added.
Figure E.6: Causal vs non-causal training for branch-trunk architectures. Time response and energy spectrum comparing non-causal training (left pair) and causal training (right pair) for DeepONet and DeepOKAN. Non-causal DeepOKAN fails to converge entirely, while non-causal DeepONet achieves only marginal convergence with substantially degraded spectral fidelity compared to its causal counterpart.
Model
Field Error
Log Spectral Error
Barron Norm Error
Accel. Error
C
NC
C
NC
C
NC
C
NC
DeepONet (SIREN)
0.0008
0.0211
0.1257
0.0961
0.0019
0.0054
0.0506
0.0277
DeepONet (SIREN, BSP)
0.0009
−
0.1366
−
0.0027
−
0.0619
−
DeepONet (Tanh)
−
−
−
−
−
−
−
−
DeepONet (Tanh, BSP)
0.0006
−
0.1103
−
0.0012
−
0.0303
−
DeepOKAN
0.0022
−
0.1230
−
0.0029
−
0.0584
−
Appendix
Table E.3: Causal (C) vs non-causal (NC) training for branch-trunk architectures (SOAP optimizer). − Failed to converge (field error ≥0.023 ).
Figure E.7: Gradient cosine similarity over training. Cosine similarity between L2 field-loss and BSP spectral-loss gradients as a function of training completion. Dashed lines indicate causal training; solid lines indicate non-causal training. DeepONet (SIREN) under causal training diverges toward strongly negative values, confirming increasing gradient conflict over the course of training.
Partial differential equations (PDEs) often have high-frequency and multi-scale features that neural networks struggle to approximate. Physics-Informed Neural Networks (PINNs) build the governing equations directly into training, but suffer from spectral bias: they learn low-frequency components faster than high-frequency ones. Techniques such as Fourier feature embeddings and sinusoidal activations address this, but most studies assume they help across the board without checking which spectral regimes actually benefit. We introduce a dual-branch, spectrally-gated architecture (DBSG-PINN) that splits low- and high-frequency components into separate subnetworks joined by an adaptive gate, and use it to run a partially controlled ablation of frequency decomposition and spectral routing. We test this on five one-dimensional benchmark PDEs, ranging from smooth, single-scale problems to oscillatory, multi-scale ones. Frequency decomposition helps most on the spectrally complex benchmarks, cutting relative L2 error by up to 59.2% on a multimodal wave problem, but gives little benefit on smoother PDEs. On one benchmark (1D Wave), it performs substantially worse than a simpler fixed-combination variant. The gate's benefit scales with how spectrally rich the target solution is: the full model's advantage over the ablations is largest on multi-scale benchmarks and smallest (or negative) on single-scale ones, consistent with the gate exploiting frequency structure rather than acting as noise,though we do not directly visualize or quantify its spatial activations in this study. All results come from a single training seed across five 1D benchmarks, so we present this as an exploratory study meant to raise questions rather than answer them, and outline the additional seeds and benchmarks needed to test whether the pattern holds.
Physics-Informed Neural Networks (PINNs) typically exhibit spectral bias, where some frequencies of the target function converge more slowly than others. In this work, we analyze the training dynamics of Fourier Feature PINNs in the Neural Tangent Kernel regime to address this limitation. We derive an explicit evolution equation to estimate the residual error in the frequency domain, demonstrating that the convergence rate of specific frequencies is primarily governed by the product of the differential operator's symbol and the spectral density of the initialization weights. Leveraging this theoretical insight, we propose an informative initialization strategy that tailors the initial weight distribution to the specific PDE being solved. With this method, we can diminish the operator-induced spectral bias, balancing the convergence rates across the frequency spectrum and achieving better prediction accuracy. Numerical experiments on linear and nonlinear partial differential equations confirm that this initialization strategy improves learning dynamics and approximation accuracy across frequencies compared to standard initialization methods, with no additional training cost.
Juan Molina, Paris Perdikaris, Mircea Petrache +2
National Center for Artificial Intelligence, Santiago, Chile · Department of Mechanical Engineering and Applied Mechanics, University of Pennsylvania, Philadelphia, PA, USA · Department of Mathematics, Pontificia Universidad Católica de Chile, Santiago, Chile +3
We present pseudo-differential enhanced physics-informed neural networks (PINNs), an extension of gradient enhancement but in Fourier space. Gradient enhancement of PINNs dictates that the PDE residual is taken to a higher differential order than prescribed by the PDE, added to the objective as an augmented term in order to improve training and overall learning fidelity. We propose the same procedure after application via Fourier transforms, since differentiating in Fourier space is multiplication with the Fourier wavenumber under suitable decay. Our methods are fast and efficient. Our methods oftentimes achieve superior PINN versus numerical error in fewer training iterations, potentially pair well with few/fixed samples in collocation, and can achieve breakthrough rather than gradual effects on the error. Moreover, our methods are suitable for fractional derivatives. We establish that our methods, due to the dynamical effects, improve spectral eigenvalue decay of the neural tangent kernel (NTK), and so our methods contribute towards the learning of high frequencies in early training, mitigating the effects of frequency bias up to the polynomial order and possibly greater with smooth activations. In particular, our primary contribution can be characterized as gradient enhancement affects spectral bias, and we specialize to Fourier space, although our theory does not prove the ambient case exactly due to a result with the Plancherel theorem on the descent paths. Our methods accommodate advanced techniques in PINNs, such as Fourier feature embeddings. A pitfall of discrete Fourier transforms via the Fast Fourier Transform (FFT) is mesh subjugation, and so we demonstrate compatibility of our methods for greater mesh flexibility and invariance on alternative Euclidean and non-Euclidean domains via Monte Carlo methods and otherwise, although possessing trade-offs of their own.
Andrew Gracyk
Department of Mathematics, Purdue University, West Lafayette, IN 47907, United States