Multi-stage Neural Networks: Function Approximator of Machine Precision

1. INTRODUCTION

• In scientific and physical applications, neural-network training can plateau and few existing methods reduce training error below O(10^-5).The introduction contrasts this limitation with classical numerical methods that can improve accuracy by reducing grid size.Source: §1 · p. 1
• Each training stage adds a separate network optimized to fit the residue from the preceding stage.The paper says this stage-by-stage process maintains a high convergence rate throughout iterations.Source: §1 · pp. 1–2
• The combined networks can approach double-precision machine error, O(10^-16).The paper presents multi-stage training as a way to sustain convergence across iterations.Source: §1 · p. 2
• The paper extends the multi-stage scheme from regression to PINNs and combined forward-and-inverse problems.It also outlines stage-setting algorithms and further discussion of the method's challenges.Source: §1 · p. 2

2.1. Limitation of regular neural network training

• After 5,000 iterations, the single tanh network’s loss plateaus around O(10^-7), leaving error around 10^-4 despite larger-network and additional-data tests.The passage reports that the plateau was not optimizer-specific.Source: §2.1 · p. 2

2.2.1. Magnitude of the second neural network

• Each higher-stage network fits the preceding combined model’s residual after normalization by that residual’s root mean square magnitude.The final model combines the stage outputs weighted by their respective residual magnitudes.Source: §2.2.1 · p. 4, §2.3 · p. 8

2.2.2. Frequency of the second neural network

• For high-frequency residuals, sine activation in the first hidden layer and a scaled input weight improve fitting; the modified scale factor should satisfy πf_d < κhat < N_d/6 to limit overfitting.When κhat ≥ 60, the example networks overfit, despite faster training-loss reduction.Source: Figure 4 · p. 6, §2.2.2 · p. 7

2.3. Algorithm of multi-stage training scheme for regression problems

• The second-stage network fits the first network’s residual; the two-stage loss reaches O(10^-14), and the combined fit reaches single-float machine precision with fewer weights than the single-stage network.Further loss reduction requires double precision.Source: Figure 1 · p. 3, §2.3 · p. 7
• After four stages, the combined regression model approaches double-float machine precision, with multi-stage error decay close to exponential rather than the single-stage 1/n_iters decay.The residue magnitude and dominant frequency follow a power law whose exponent is reported as independent of target functions, network size, and data-point count.Source: §2.3 · p. 8, Figure 5 · p. 9
• For 2D regression, multi-stage training is slower than the 1D case but remains faster than single-stage linear error decay, and the residue relation follows the same power law with α ≈ 1/6.The required number of training points also increases with the residual’s dominant frequency to support higher stages.Source: §2.3 · pp. 8, 10

3. MULTI-STAGE TRAINING FOR PHYSICS-INFORMED NEURAL NETWORK

• PINNs pose two higher-stage challenges: estimating error magnitude and frequency without the exact solution, and choosing the equation weight that balances data and equation losses.The equation residue is the available basis for inferring the lower-stage prediction error.Source: §3.1 · p. 11, §3 · p. 11

3.1.1. A simple example

• In the ODE example, the equation residue and prediction error share a dominant frequency near 1.5, while the residue magnitude is about ten times larger.The estimated error magnitude follows from dividing the residue magnitude by 2πfd.Source: §3.1.1 · pp. 12–13

3.1.2. Magnitude and frequency estimation for general differential equations

• For general equations, higher-stage error frequency is tied to the equation residue, and small errors make higher-stage governing equations essentially linear.The magnitude and frequency estimates are sufficiently accurate for most equations, but may not hold exactly for two types of nonlinear equations.Source: §3.1.2 · pp. 14–15

3.1.3. Importance of magnitude and frequency for higher-stage PINN training

• Only PINNs with both magnitude and frequency set according to the derived relations rapidly converge to the exact solution in the tested Poisson problem.The comparison assumes the correct equation weight γ.Source: Figure 8 · p. 15, §3.1.3 · p. 16

3.1.4. Algorithm for determining the solution magnitude for higher-stage PINN training

• Algorithm 2 estimates solution magnitude by iteratively adjusting a scale coefficient until the differential-operator magnitude is comparable to the source magnitude.It uses the source function as a guess solution, based on the premise that solution and source share a dominant frequency.Source: §3.1.4 · p. 16
• Algorithm 2 estimates the Poisson solution magnitude as 1.41 × 10^-3, close to the exact 1.43 × 10^-3.The algorithm mainly applies to linear differential equations with one dependent variable; nonlinear equations or systems with multiple dependent variables may require a more advanced algorithm.Source: §3.1.4 · p. 16

3.2. Second challenge: equation weight γ for higher-stage network

• For high-frequency solutions, equation loss can be six orders of magnitude larger than data loss, so γ = 0.5 leaves data loss decaying slowly and prediction error around 5%.The section contrasts this with low-frequency equations, where 0.1 < γ < 0.5 can make the data-loss contribution larger.Source: §3.2 · pp. 17–18

3.2.1. Theoretical approach of determining γ

• For equation (3.20), the theoretical expression gives γ = 2 × 10^-6 and reduces relative error by more than one hundred times versus γ = 0.5, but requires loss estimates that are generally unavailable without prior solution knowledge.Both losses decrease by more than five orders of magnitude with this γ.Source: §3.2.1 · p. 18

3.2.2. Algorithm for determining γ for general equations

• Algorithm 3 estimates the two loss-convergence rates after brief pre-training, adjusts γ until their ratio Rc meets a criterion, and requires retraining from the beginning after each update.The method is presented as more general and adaptable than the theoretical expression.Source: §3.2.2 · pp. 19–20
• γ ≈ 10^-4 from Algorithm 3 gives the smallest prediction errors among γ = 0.5, γ = 2 × 10^-6, and the algorithm-selected value, with 0.1 < Rc < 5 marking the low-error range.The algorithm-selected γ maximizes the convergence rates of both losses in the reported comparison.Source: §3.2.2 · p. 20

3.3.1. Optimization method and re-sampling collocation points

• For high-frequency solutions, Adam has a better overall convergence rate than L-BFGS, and re-sampling collocation points every few iterations further increases convergence rate.The cited optimizer comparison concerns solving equation (3.20).Source: §3.3.1 · pp. 21–22
• Adam with stochastic gradient descent (SGD) can attain optimal precision even when collocation points Nc fall below Ncrit ≈ 3100, unlike fixed-point methods that are more sensitive to point count.For the studied high-frequency equation, fixed-point methods improve sharply near Ncrit, with only marginal gains from additional points.Source: Figure 10 · p. 21, §3.3.1 · p. 22

3.3.2. Advanced methods from the literature: RAR and gPINNs

• Residual-based adaptive refinement (RAR) adds collocation points in regions of high equation residue, addressing local regions where gradients exceed the domain-averaged gradient.The text presents RAR as a way to reduce equation residue across the domain.Source: §3.3.2 · p. 22
• gPINNs add the gradient of the equation residue to the loss, constraining higher-derivative information and potentially improving convergence when the gradient-loss weight is chosen appropriately.The section relates residue magnitude to its gradient as ||r||/||∇r|| ∼ O(2πfd)^-1 for high-frequency solutions.Source: §3.3.2 · pp. 22–23

3.4. Algorithm of multi-stage training for PINNs

• The scheme estimates each stage’s magnitude prefactor and frequency scale from the preceding equation residue, then repeats training until the residue is sufficiently small and combines all stages.The equation is normalized so its largest term is around O(1).Source: §3.4 · pp. 23, 25
• Because PINNs lack solution data, the next network is trained through a combined ansatz with the prior trained network fixed, using the original governing equation rather than directly fitting the unknown solution error.At each higher stage, the new network and the combined network from the previous stage form the ansatz.Source: §3.4 · pp. 23, 25
• Multi-stage training reduces both data loss and equation loss by more than 20 orders of magnitude versus single-stage training in the same iterations, and three-stage gPINN training reaches double-floating-point machine precision.The single-stage network has approximately the same total number of weights as the three-stage networks; the final error oscillation is attributed primarily to round-off error.Source: §3.4 · pp. 23–24, Figure 11 · p. 24, Figure 12 · p. 26
• For multi-stage PINNs, dominant frequency fd and error ϵ follow fd ∼ ϵ^-α, with α around 1/7 for regular PINNs and 1/8 for gPINNs, compared with 1/6 for regression.Each stage’s prediction error can be reduced by 3–5 orders of magnitude within 10^5 iterations in the reported three-stage ODE example.Source: §3.4 · pp. 23–25

3.5. Application to 2D partial differential equations

• After three stages using RAR and gPINNs, the 2D prediction error is around O(10^-11), seven orders of magnitude lower than single-stage training.The convergence rate is slightly slower than for the 1D problem.Source: §3.5 · p. 26, Figure 13 · p. 27
• For 2D equations with gPINNs, the dominant-frequency and error relation follows the same power law as 1D problems, with exponent α ≈ 1/8.The inset reports that 2D problems require more iterations than 1D problems to reach the same accuracy.Source: §3.5 · p. 26, Figure 13 · p. 27
• A classical finite-difference estimate for O(10^-11) accuracy requires O(10^9) grid points across the domain, while the reported multi-stage network uses around 3 × 4 × 30^2 ≈ 10^4 weights and biases.The paper describes the neural-network parameter count as five orders of magnitude smaller than the discretized grid-point count.Source: §3.5 · pp. 26–27

4. GENERALIZATION OF MULTISTAGE PINN TO A COMBINED FORWARD-AND-INVERSE PROBLEM

• The tested task seeks a smooth solution of the inviscid Burgers’ equation while simultaneously inferring λ; smooth solutions occur only at λ = 1/(2+2i) for i = 1, 2, 3....The PINN loss includes a smoothness constraint on the third derivative of the equation residue, evaluated near the origin, which acts as an additional boundary condition for determining λ.Source: §4 · pp. 27–28
• Second-stage training improves prediction error and inferred λ by four orders of magnitude, but a dominant low-frequency error profile prevents further reduction with the existing scheme.The figure compares the first-stage network at inferred λ1 with the exact smooth solution at λg = 0.5.Source: §4 · p. 28, Figure 14 · p. 29
• The third stage’s two-network approach reduces both high-frequency equation-residue error and low-frequency inferred-parameter error, bringing prediction toward double-precision machine precision.The combined networks reduce prediction error by another three orders of magnitude.Source: §4 · pp. 29–30
• The low-frequency source associated with inferred-λ error contributes prediction error of approximately O(10^-12), exceeding the approximately O(10^-13) error associated with the high-frequency equation residue.This explains why the low-frequency profile dominates the second-stage prediction error.Source: §4 · p. 30
• The method estimates inferred-parameter error using √Ls when the exact parameter is unknown, since the reported inference error is proportional to that quantity.The estimate also supplies the magnitude prefactor for the low-frequency network.Source: §4 · p. 30
• A single low-frequency correction is insufficient because the smoothness constraint depends on higher derivatives and the high-frequency equation residue must also be reduced.The proposed third stage therefore uses one network for each source in the higher-stage equation.Source: §4 · p. 30

5. DISCUSSION

• The error magnitude follows an inverse power law with dominant error frequency, with α ≈ 1/6 for regression, α ≈ 1/7 for regular PINNs, and α ≈ 1/8 for gPINNs.The authors compare this behavior to adding higher-frequency terms in a Fourier series.Source: §5 · p. 31
• Fourier series have α ≈ 0.5, larger than MSNN exponents, so MSNN errors can have frequencies several orders lower at the same error magnitude.The authors interpret this as evidence that MSNNs can accurately approximate target functions and their high-derivative information.Source: §5 · p. 31
• MSNNs remain a proof of concept rather than a substitute for classical numerical methods, and their convergence is slower for 2D than 1D problems, with greater challenges expected at higher dimensions.The authors describe MSNNs as complementary to classical methods.Source: §5 · p. 31
• Steep-gradient regions can produce local error or equation-residue peaks that hinder successive-stage error reduction and must be removed before training proceeds to the next stage.PINN solutions for stiff, nonlinear, or singular-perturbation equations are beyond the paper’s scope.Source: §5 · p. 31
• Choosing when to move to the next stage trades faster convergence against maximizing error reduction within the current stage; selecting stage-wise network sizes to reduce total parameters also remains open.The authors identify both as questions for further investigation.Source: §5 · p. 31

6. CONCLUSION

• MSNNs add stages sequentially, with each new network trained to capture the residue from the previous stage.The authors frame the approach as inspired by perturbation theory and report that three stages can reach machine precision.Source: §6 · p. 31
• Staged training reduces previously plateaued error at each successive stage, enabling progressively finer target-function details to be captured despite spectral bias in single-stage networks.The authors attribute spectral bias to difficulty capturing the target’s full spectrum in one training stage.Source: §6 · p. 32
• Each new-stage network is designed from prior-stage error, using an optimal magnitude prefactor and scale factor with sine activation to represent small-magnitude, high-frequency predictions.For regression, the prefactor is the prior-stage error’s RMS magnitude and the scale factor is proportional to its dominant frequency.Source: §6 · p. 32
• For PINNs, prior-stage equation residues are used to infer prediction-error magnitude and frequency because the prediction error is not directly available.The conclusion also discusses tuning equation weight, collocation-point count, optimizer, and techniques such as RAR and gPINNs.Source: §6 · p. 32
• MSNNs substantially reduce prediction error for regression and PINNs, approaching machine precision, and are reported to solve combined forward-and-inverse problems to that accuracy.The authors describe the latter task as typically challenging for classical numerical methods.Source: §6 · p. 32

Appendix A. NEURAL NETWORK ERROR UNDER DIFFERENT SETTINGS

• Increasing layers or hidden units leaves RMS error unchanged, while added training data yields a reduction smaller than the standard deviation across eight random initializations.The comparison uses RMS error between the trained network and target-function data; the figure reports standard deviations across eight experiments.Source: Appendix A · p. 33, Figure 15 · p. 33
• L-BFGS has a higher convergence rate than Adam but quickly reaches the same loss plateau, indicating that the plateau is not optimizer-specific.L-BFGS is described as a quasi-Newton method and Adam as a first-order gradient-descent method.Source: Appendix A · p. 33, Figure 15 · p. 33

• Without normalization, network fitting is difficult when data magnitudes are far from 1.Source: Figure 16 · p. 33
• Data with magnitudes much larger than 1 are hard for the network to fit.Source: Figure 16 · p. 33
• Data with magnitudes much smaller than 1 are also hard for the network to fit.Source: Figure 16 · p. 33

• For high-frequency u1, its nonlinear term can become dominant, so its magnitude and frequency must be reassessed by balancing that term against the preceding equation residue.Source: Appendix C · p. 34
• For the high-frequency case, the scale factor can be set using the larger dominant frequency of the equation residue.The passage relates this choice to a criterion for fitting high-frequency functions.Source: Appendix C · p. 34
• A singular perturbation has a coefficient before the highest-order derivative that is much smaller than the other coefficients; boundary-layer problems are an example.Source: Appendix C · p. 34
• When α < (2πf_d)^-2, the second-derivative term, rather than the highest-order derivative, dominates the right-hand side.Source: Appendix C · p. 34
• For singular perturbations, narrow inner regions with very large local gradients make both first-stage and higher-stage network training difficult.The paper states that solving this challenge is beyond its scope.Source: Appendix C · p. 35

• As the blow-up time approaches, the self-similar profile is assumed to reach a steady state governed by the transformed equation.Source: Appendix D · p. 35
• The multi-stage neural networks are used to infer λ, the rate at which the Burgers’ equation’s singularity forms.Source: Appendix D · p. 35
• Figure 17 contrasts the non-smooth λ=0.4 solution, whose third derivative has cusps and fourth derivative is discontinuous at the origin, with the all-order-continuous λ=0.5 solution.Source: Figure 17 · p. 36
• λ=1/(2+2i), for i=0,1,2,..., gives solutions smooth at every derivative order; other λ values can produce non-smooth solutions at the origin.The example λ=0.4 has a discontinuous fourth derivative there.Source: Appendix D · p. 36
• The PINN implementation imposes odd symmetry and a smoothness constraint based on the equation residue to identify the smooth solution and determine λ.The smoothness constraint is treated as a boundary condition.Source: Appendix D · p. 36