Modeling the Dynamics of PDE Systems with Physics-Constrained Deep Auto-Regressive Networks

• Surrogate models are motivated by the high cost of accurate PDE simulations, especially when repeated runs are needed for optimization or inverse problems.They aim to imitate standard numerical solvers at substantially lower computational cost.Source: §1 · p. 2
• Physics-constrained loss functions enable surrogate models to use known governing equations without training data, addressing the expense or scarcity of simulation data.Prior physics-constrained approaches mainly targeted individual PDE solutions or time-independent surrogate models.Source: §1 · pp. 2–3
• The paper proposes a deep auto-regressive dense encoder-decoder and physics-constrained training algorithm for dynamical PDE surrogate modeling.This generalizes earlier physics-constrained surrogate models from time-independent to transient PDE systems.Source: §1 · p. 3
• The proposed model is extended with stochastic weight averaging Gaussian to quantify both epistemic and aleatoric uncertainty.The framework evaluates prediction accuracy and associated uncertainty across the studied PDE systems.Source: §1 · p. 4
• The framework is tested on chaotic or turbulent dynamics, multiple shock-wave interactions, and a two-dimensional coupled nonlinear PDE system, with comparisons against state-of-the-art numerical methods.The paper also compares the proposed surrogate’s computational efficiency with traditional numerical approaches.Source: §1 · p. 4

2. Problem Definition

• The problem is to learn surrogate models for arbitrary transient PDE systems from governing equations, without requiring output training data and while allowing unknown, variable, or stochastic initial states.The PDE state evolves over a spatial domain with boundary conditions, while the initial state is a real-valued random field whose density may be unknown.Source: §2 · pp. 4–5
• The proposed surrogate performs discrete time integration with a convolutional encoder-decoder, mapping k + 1 previous states to the next state at a specified time-step Δt.The learned dynamics are intended to be invariant to the current reference time and avoid matrix inversions, iterative relaxations, and multi-step numerical processes.Source: §2 · pp. 5–6
• The model predicts an entire response autoregressively by repeatedly applying the same network, requiring only the discretized initial state to generate an arbitrary number of time-steps.The first input contains only the initial state to kick-start the time series; later inputs are formed from the model’s own predictions.Source: §2 · p. 6
• Physics-constrained training minimizes the discrepancy between the network’s one-step prediction and a discrete numerical time integrator for the governing equation.When the integrator is consistent with the governing equation, this objective can be interpreted as minimizing the residual of the discretized PDE.Source: §2 · pp. 6–7

2.1. Surrogate Modeling of Dynamical Systems

• The deterministic surrogate is trained on a finite set of sampled initial conditions and should predict the governing-PDE solution for any new initial condition drawn from the same underlying distribution.Because labeled outputs are unnecessary, this physics-constrained surrogate formulation is interpreted as unsupervised learning.Source: §2.1 · p. 7
• The probabilistic surrogate predicts a response density rather than a point estimate, with samples from that distribution required to satisfy the governing PDEs for their respective initial states.This extension addresses the limited coverage of possible initial conditions and responses by finite training data and model expressibility.Source: §2.1 · pp. 7–8

3. Auto-regressive Dense Encoder-Decoder

• AR-DenseED sequentially predicts future states from k + 1 previous states, reusing each prediction as the next input without latent recurrent connections or LSTM-specific gates.The model uses convolutional encoders and decoders interleaved with dense blocks to map previous states through a lower-dimensional latent representation and predict the next state.Source: §3 · pp. 8, 10, Figure 1 · p. 9
• Multiple previous time-steps improve training stability, but using too many slows training without considerable predictive gains.The number of previous states is a tunable hyper-parameter, and periodic boundaries are enforced with circular padding in the model convolutions and loss evaluation.Source: §3 · p. 8, §3.1 · p. 13

3.1. Physics-Constrained Loss Function

• Minimizing the L2 difference between network predictions and a consistent numerical time integrator is equivalent to minimizing the discretized PDE residual across the domain.The loss can use any time-integration algorithm; for large time-steps with CFL > 1, predictions evaluate spatial gradients to provide implicit-like integration without matrix inversions.Source: §3.1 · pp. 11–13

3.2. AR-DenseED Training

• Training progressively increases the number of autoregressively predicted time-steps, enabling multi-step back-propagation while limiting the interval to avoid vanishing gradients.Because time is not an explicit input, the model has no fundamental limitation on predictive range or initial conditions and can learn dynamics absent from the provided initial states.Source: §3.2 · p. 13, Figure 3 · p. 14

4. Bayesian AR-DenseED

• BAR-DenseED interprets epistemic uncertainty as model confidence and aleatoric uncertainty as truncation error from the physics-constrained numerical time-integrator.Epistemic uncertainty reflects limited scenarios or model expressibility, while aleatoric uncertainty represents temporal and spatial discretization errors.Source: §4 · p. 15

4.1. Posterior Formulation

• The posterior formulation treats numerical-integrator predictions as target data, adds Gaussian output noise with learnable precision, and remains data-less because targets are generated on-the-fly.Under a Markov assumption, the sequence likelihood factors across time-steps, and maximum likelihood recovers the standard L2 discretized-PDE residual.Source: §4.1 · p. 16

4.2. Posterior Approximation

• SWAG approximates the BAR-DenseED posterior by collecting stochastic-gradient parameter samples after MAP training, providing simple, non-invasive, low-overhead Bayesian inference.Its sampling learning rate must explore the local minimum without jumping between modes because the approximation uses a single-mode Gaussian.Source: §4.2 · pp. 19–20

4.3. Predictive Statistics

• Predictive statistics are obtained by Monte Carlo model averaging, with independently autoregressed posterior samples producing uncertainty that is expected to increase at later time-steps.The predictive covariance combines conditional output uncertainty with parameter uncertainty, while zero-mean additive output noise does not affect the predictive expectation.Source: §4.3 · pp. 21–22

5. Kuramoto-Sivashinsky Equation

• The K-S experiment uses a chaotic periodic domain [0, 22π], 96 spatial cells, ∆t = 0.1, and fully developed turbulent initial states with T ≥100.Training used approximated Fourier-series scenarios, whereas testing used 200 true turbulent initial states from a numerical simulator.Source: §5 · pp. 22–23

5.1. AR-DenseED Deterministic Predictions

• AR-DenseED remains qualitatively physical and stable for at least 5000 time-steps, although it matches the numerical solver only through approximately t = 30 before chaotic divergence.The model’s energy spectrum agrees most closely at larger wavelengths and places peak energy near Hz∈[0.1, 0.15].Source: §5.1 · pp. 24–25
• The AR-DenseED simulation is computationally less effective than the spectral ETDRK4 scheme for 5000 time-steps, despite demonstrating complex chaotic-system modeling.The comparison uses averaged wall-clock times from 10 independent simulation runs.Source: §5.1 · p. 25, Table 1 · p. 27

5.2. BAR-DenseED Probabilistic Predictions

• BAR-DenseED posterior samples begin similarly for t < 40 but diverge into unique responses afterward, reflecting the K-S system’s chaotic butterfly effect.SWAG used 100 parameter samples, with very small learning rates required to avoid unstable autoregressive responses.Source: §5.2 · p. 27
• BAR-DenseED and AR-DenseED have essentially equivalent predictive performance because later-time fields average toward zero, despite apparently lower Bayesian mean squared error.The Bayesian model’s apparent improvement at later times is therefore attributed to field averaging rather than superior prediction.Source: Figure 6 · p. 26, §5.2 · p. 27
• BAR-DenseED uncertainty produces reasonable error bars capturing the true spectral solution at largest wavelengths, Hz < 0.3, while smaller-wavelength energy remains sample-consistent.The spectral density was estimated from 30 posterior model samples.Source: Figure 7 · p. 26, §5.2 · p. 27

6. 1D Viscous Burgers’ Equation

• The Burgers’ test problem uses random multi-wave initial conditions at viscosity ν = 0.0025, producing multiple shocks that intersect, combine, and follow new trajectories.The model uses five previous time-steps, 512 spatial points, ∆t = 0.005, 2560 training scenarios, and 200 test scenarios.Source: §6 · pp. 28–30

6.1. AR-DenseED Deterministic Predictions

• AR-DenseED accurately predicts sharp shock formations, intersections, and long-horizon behavior without training data, but small shock-intersection errors can compound over time.Testing predicts 400 time-steps, including extrapolation beyond the training range, against high-fidelity FEM targets using a surrogate time-step five times larger than FEM’s.Source: §6.1 · p. 30, Figure 10 · p. 31
• For the 200-case test set, MSE and ESE initially spike during shock formation and intersection; ESE remains stable before growing in extrapolation, while mean errors exceed medians because of poor outlier cases.AR-DenseED is significantly cheaper computationally than FEM and FDM for simulating 400 time-steps.Source: §6.1 · pp. 31–32, Table 2 · p. 33

6.2. BAR-DenseED Probabilistic Predictions

• BAR-DenseED generally outperforms the deterministic model by reducing poor outlier predictions through model averaging, with uncertainty concentrated at shock formation and intersection locations.The predictive expectation remains structurally accurate beyond the training range, while widening error bars indicate reduced confidence during extrapolation.Source: §6.2 · pp. 33, 36, Figure 13 · p. 35
• BAR-DenseED’s uncertainty increases as predictions leave the training range, yet its error bars capture the true solution’s structure during extrapolated shock dynamics.At earlier times, BAR-DenseED compares well with FEM and FDM, whose solutions differ because of their numerical discretizations.Source: §6.2 · pp. 35–36, Figure 14 · p. 36, Figure 15 · p. 37

7.1. AR-DenseED Deterministic Predictions

• AR-DenseED accurately predicts complex 2D coupled Burgers’ structures and discontinuities through t > 0.5, with errors concentrated at shock interfaces and wave fronts.The target is a high-fidelity FEM solution interpolated from 128 × 128 to 64 × 64; predictions cover 200 time-steps from t = 0.Source: §7.1 · p. 40, Figure 18 · p. 41, Figure 19 · p. 42
• AR-DenseED’s mean and median errors remain stable during training-time prediction, but mean error grows faster than median during extrapolation, indicating poor outlier cases.MSE and ESE are evaluated over 200 test scenarios at each time-step.Source: §7.1 · p. 42, Figure 20 · p. 43
• AR-DenseED simulates 200 time-steps substantially faster than FEM across spatial discretizations, with GPU acceleration and batching offering further efficiency gains.The comparison uses averaged wall-clock times from 10 independent simulation runs.Source: §7.1 · p. 42, Table 3 · p. 43

7.2. BAR-DenseED Probabilistic Predictions

• BAR-DenseED’s predictive standard deviation captures the true solution for almost all times, while uncertainty concentrates on the leading faces of shocks and waves.Posterior predictions use SWAG with 90 parameter samples, and predictive expectation and variance use 30 model samples.Source: §7.2 · pp. 43–44, Figure 23 · p. 45, Figure 24 · p. 46
• BAR-DenseED achieves smaller and more stable MSE and ESE during extrapolation beyond the training-time range.These errors are computed using the Bayesian model’s predictive expectation.Source: Figure 20 · p. 43, §7.2 · p. 44

Conclusion

• The model accurately reproduced turbulent statistics, multi-shock formation and intersections, and complex wave dynamics across three PDE systems, while extrapolating beyond the training time range.The systems were the Kuramoto-Sivashinsky equation, low-viscosity 1D Burgers’ equation, and 2D coupled Burgers’ equations.Source: §8 · p. 47
• Physics-constrained learning trains the deep auto-regressive convolutional surrogate with zero output training data and Bayesian uncertainty quantification covering epistemic and aleatoric uncertainty.The governing equations formulate the training loss, while the Bayesian framework is built on top of the deep learning model.Source: §8 · p. 47
• Convolutional neural networks can build PDE surrogate models that are order of magnitudes faster than state-of-the-art numerical solvers.Their use requires discretizing spatial and temporal derivatives, introducing challenges including truncation error, oscillations, and convergence criteria.Source: §8 · p. 48
• Future work should extend the model to larger systems while reducing exponentially increasing time-series training costs, potentially through self-attention mechanisms.Suggested applications include Navier-Stokes equations, coupled transport through porous media, and combustion.Source: §8 · p. 48
• Future extensions could combine data with physics-constrained learning to predict fine-resolution dynamics when system states are available only at sparse time intervals.This hybrid framework is proposed to bridge missing temporal observations.Source: §8 · p. 48

Model architecture

• The model encodes the 1D input {u_n, u_n−1} into latent variables of dimensionality d/2, then decodes them to predict u_n+1.Its three components are an encoding convolution, dense block, and decoding block.Source: Appendix A · p. 49

Training

• Smaller models were more stable and faster to train, while circular padding supported learning periodic boundaries.Training used ADAM with an initial learning rate of 1e − 3; exponential learning-rate decay produced the most stable, consistent, and accurate results.Source: Appendix A · p. 49

Model architecture

• The AR-DenseED model has 4,821 learnable parameters and uses the two previous time-steps as inputs, followed by encoding, a dense block, and decoding.The dense block has growth rate 4 and length 4.Source: Figure A.26 · p. 49

Network blocks

• The dense block uses residual connections that stack convolution outputs and increase feature channels, while the decoding block uses nearest-neighbor up-sampling to reach target dimensions.Both blocks include batch normalization, ReLU activations, and convolutions.Source: Figure A.27 · p. 50, Figure A.28 · p. 50

Training Initial States

• Training uses truncated Fourier series with random coefficients instead of simulator data, preserving the model’s zero-training-data philosophy.These initial states serve as the training scenarios from which the model starts.Source: Appendix A.1 · p. 51
• The initialization samples amplitude as a = N(0, 0.5) + a0, with a0 set to the prior mean amplitude estimate 2.5.The construction also specifies λn = [1, N(0, 2), 1] and c = 2π U[0, 1].Source: Appendix A.1 · p. 51
• The resulting function is designed to provide a physically realizable initial condition for AR-DenseED to explore the Kuramoto–Sivashinsky system.Its construction uses the domain size L and the estimated number of unstable modes k0.Source: Appendix A.1 · p. 51

Training difficulty

• Sharp gradients made the 1D viscous Burgers’ system the most difficult to train, likely causing exploding gradients during optimization.The Burgers model was approximately three times larger than the K-S model, reflecting the system’s greater training difficulty.Source: Appendix B · p. 51

Training stability

• Using more input time-steps than the K-S system improved learning stability, but making the model too large caused unstable predictions because SWAG could not approximate a good posterior.This instability occurred regardless of the learning rate for sampled models.Source: Appendix B · p. 51

Model architecture

• The 1D Burgers’ AR-DenseED model uses five previous time-steps, an encoding block, one dense block, and a decoding block with 13,442 learnable parameters.The encoder maps the 1D input to latent variables of dimensionality d/2 before decoding the next prediction.Source: Appendix B · p. 51, Figure B.29 · p. 52

Training parameters

• Table B.5 lists the AR-DenseED and BAR-DenseED training parameters used for the 1D Burgers’ system.The table provides the appendix’s training-parameter specification for both model variants.Source: Table B.5 · p. 52

Appendix C. 2D Coupled Burgers’ System

• The model predicts x- and y-velocity components from both components at previous time steps, giving six input channels when using three histories.It encodes {u_n, v_n, u_{n−1}, v_{n−1}, u_{n−2}, v_{n−2}} into latent variables of size d/2 × d/2 and decodes {u_{n+1}, v_{n+1}}.Source: Appendix C · p. 52
• Although the 2D system has the most complex dynamics, it was easier to train and substantially more stable and consistent than the other 1D systems.The authors largely attribute this behavior to Sobel filters that approximate gradients and damp higher frequencies, preventing oscillations.Source: Appendix C · p. 52
• The AR-DenseED architecture has 71,953 learnable parameters and uses an encoding block, one dense block, and a decoding block with five previous time steps as inputs.The dense block has growth rate 4 and length 4.Source: Figure C.30 · p. 53
• Additional AR-DenseED and BAR-DenseED training parameters for the 2D coupled Burgers’ system are provided in Table C.6.The supplied passage identifies the table but does not report its parameter values.Source: Appendix C · p. 52, Table C.6 · p. 53