Non-Conservative Error in Numerical PDEs

• Non-conservative error is the domain-dependent discrepancy arising when conservation or balance laws are not preserved in numerical formulations.
• It affects hyperbolic PDEs, interface-capturing, and long-time integrations by causing wrong shock speeds, spurious interfaces, and secular drift of invariants.
• Researchers address these errors through corrective schemes such as global-flux reconstructions, a posteriori level-set shifts, and energy-residual adjustments to restore physical consistency.

Non-conservative error is a domain-dependent term for the discrepancy that appears when the relevant conservation law, balance law, or weak-solution structure is not preserved by the governing formulation, by the numerical discretization, or by the learning objective. In hyperbolic PDEs it commonly denotes wrong jump conditions, wrong shock speeds, spurious interface terms, or imbalance between fluxes, sources, and non-conservative products; in interface-capturing it denotes spurious mass or volume change; in long-time integration it denotes secular drift of invariants; in generalized dynamics it denotes the gap between a conserved time-symmetry quantity and the usual mechanical energy, or a destabilizing rotational contribution; and in machine learning it denotes residuals induced by training on non-conservative forms without path-consistent or entropy-consistent structure (Neelan et al., 27 Jun 2025, Long et al., 2022, Ushveridze, 2022, Ochoa et al., 22 Jan 2025).

1. Weak solutions, shocks, and non-conservative products

For hyperbolic balance laws, the canonical distinction is between conservative form, $\partial_t U + \partial_x F(U)=0$, and non-conservative or quasi-linear forms such as

$\partial_t V + A(V)\partial_x V = 0$

or $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$. These formulations are equivalent for smooth solutions, but they are not equivalent at shocks because the chain rule fails in the distributional sense at jump discontinuities. In that regime, conservative systems admit Rankine–Hugoniot jump conditions, $s[[U]]=[[F(U)]]$, whereas non-conservative products require an explicit path-dependent definition, typically in the Dal Maso–LeFloch–Murat framework,

$$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$

If that path-consistent structure is omitted, the resulting weak solution can have wrong shock speeds, smeared shocks, non-propagating discontinuities, or spurious interface terms (Neelan et al., 27 Jun 2025, Neelan et al., 2 Apr 2026).

This issue is not restricted to textbook conservation laws. In shallow-water moment equations, the term $B(U)\partial_x U$ makes analytical steady states difficult or impossible to characterize, and the numerical error appears as path dependence, spurious interface contributions, and flux–source imbalance unless conservative fluxes, source terms, and non-conservative products are discretized consistently. The global-flux WENO construction addresses this by writing $G=F+R$, integrating both sources and non-conservative products into the same reconstructed quantity, and using interface jumps computed by a linear path, which is described as the only choice preserving lake at rest in that framework (Ciallella et al., 1 Jul 2025).

A related but geometrically distinct setting occurs for balance laws with discontinuous pipe geometry. There, products such as $p\,\partial_x A$ are not defined classically when $A$ is only of bounded variation. The ambiguity is resolved by separating discrete Dirac contributions at geometric jumps from the continuous measure part, using a coupling function $\Xi$ at jumps and a directional derivative $\partial_t V + A(V)\partial_x V = 0$0 for the continuous part. This selects a physically consistent weak solution as the limit of piecewise-constant geometries and prevents the non-uniqueness, wrong wave structures, and incorrect shock speeds associated with an undefined non-conservative product (Colombo et al., 2021).

A common misconception is that any divergence-form rewriting is automatically physically correct. The cited CFD analyses explicitly distinguish mathematically conservative rewritings from physically conservative variables. For shallow water, rewriting the system in terms of $\partial_t V + A(V)\partial_x V = 0$1 instead of $\partial_t V + A(V)\partial_x V = 0$2 changes the shock speed because $\partial_t V + A(V)\partial_x V = 0$3 is not a conserved density, even though both forms can be written in divergence form (Neelan et al., 27 Jun 2025, Neelan et al., 2 Apr 2026).

2. Interface-capturing and volume-loss error in level-set methods

In incompressible multiphase flow simulations, the level-set method represents the interface $\partial_t V + A(V)\partial_x V = 0$4 as the zero isocontour of a signed-distance field $\partial_t V + A(V)\partial_x V = 0$5. Its central non-conservative error is mass or volume loss. The paper on conservative level-set correction identifies two sources for the enclosed volume $\partial_t V + A(V)\partial_x V = 0$6 of one phase: discretization error in advection and continuity, $\partial_t V + A(V)\partial_x V = 0$7, and reinitialization error, $\partial_t V + A(V)\partial_x V = 0$8, because Hamilton–Jacobi redistancing can slightly move the zero contour. Using a smoothed Heaviside $\partial_t V + A(V)\partial_x V = 0$9, the phase volume is

$$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$0

with absolute and relative errors

$$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$1

The rate identity

$$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$2

shows that exact incompressibility and exact level-set advection would imply $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$3; in practice the discrete residuals drive nonzero volume change (Long et al., 2022).

The proposed remedy is an a posteriori mass correction by a uniform perturbation $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$4 to the level-set field, $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$5, chosen so that the corrected volume matches the target volume. Because $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$6 is spatially uniform, $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$7, so the signed-distance property is preserved exactly. The scalar $$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$8 is computed by Newton iteration applied to

$$\partial_t U + \partial_x F(U)=S(U,x)+B(U)\partial_x U$$9

with

$s[[U]]=[[F(U)]]$0

The paper reports that $s[[U]]=[[F(U)]]$1–$s[[U]]=[[F(U)]]$2 Newton iterations are sufficient, the convergence criterion is $s[[U]]=[[F(U)]]$3, and correction every $s[[U]]=[[F(U)]]$4 time steps performs as well as correction every step. Reinitialization is likewise performed sparingly, with $s[[U]]=[[F(U)]]$5 pseudo-time iterations every $s[[U]]=[[F(U)]]$6 physical steps. Across linear advection, Zalesak’s disk, vortex deformation, falling droplet, Rayleigh–Taylor instability, bubble rising, and binary droplet collision, the conservation error is reduced to the order of machine accuracy with negligibly extra cost, often below $s[[U]]=[[F(U)]]$7 and in some cases about $s[[U]]=[[F(U)]]$8 (Long et al., 2022).

The broader significance is that this is a non-conservative error in an interface representation rather than in a hyperbolic weak solution. Geometry and force coupling remain unchanged because interface normals and curvature are computed from the reinitialized $s[[U]]=[[F(U)]]$9, and the uniform shift leaves $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$0 unchanged. This suggests a useful distinction: some non-conservative errors are defects in the transport of a geometric marker rather than in the balance law for the underlying momentum equations (Long et al., 2022).

3. Conservation defects, balance-law preservation, and long-time drift

Several numerical-analysis papers formalize non-conservative error as a local conservation defect. In arbitrary high-order ADER-DG schemes written in primitive variables, the cellwise defect is denoted $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$1 and measures the mismatch between the conservative update of cell averages and the update induced by the primitive-variable DG corrector. The modified Lax–Wendroff theorem proved in that setting states that if $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$2, together with weak consistency, TVB-type bounds, and mesh regularity, then any convergent limit is a weak solution of the conservation law. In smooth regions the paper shows $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$3; in shock-triggered troubled subcells the a posteriori conservative finite-volume correction makes the subcell defect vanish to machine precision, so $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$4 there after summation (Gaburro et al., 2024).

A related strategy appears in pressure-based nonconservative schemes for Euler equations with nonlinear equations of state. There the leading non-conservative error is the mismatch between the discrete energy balance implied by the nonconservative updates and the conservative energy residual. The residual-distribution method introduces local corrections $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$5 or $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$6 so that the internal-energy or pressure residual satisfies an exact discrete energy identity at the element level. The purpose is explicit: restore correct shock speeds, preserve pressure and velocity continuity across contacts, and ensure convergence to the relevant weak solution of the conservative Euler system even though the primary variables are nonconservative (Abgrall et al., 2017).

The same principle extends to steady-state preservation. For shallow-water moment equations with non-conservative products, the high-order global-flux WENO method reconstructs the global flux $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$7 rather than the conservative variables. The sufficient condition for arbitrary steady-state preservation is that the numerical dissipation depend only on jumps in $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$8 and vanish when $$s[[U]]=\int_0^1 A(\varphi(\sigma))\,\partial_\sigma \varphi(\sigma)\,d\sigma.$$9. This removes one major form of non-conservative error: residual imbalance between conservative fluxes, source terms, and non-conservative products at equilibrium (Ciallella et al., 1 Jul 2025).

Long-time integration of dispersive waves exhibits another manifestation. For solitary waves viewed as relative equilibria, conservative fully discrete schemes that preserve the invariant associated with translation symmetry show global error growing approximately linearly in time, whereas non-conservative methods typically exhibit quadratic growth until saturation. The paper documents this across Fornberg–Whitham, Camassa–Holm, Degasperis–Procesi, Holm–Hone, BBM–BBM, a variable-coefficient $B(U)\partial_x U$0-system, and a 2D shallow-water model. It also reports that conserving energy alone is not always sufficient: in the KdV comparison, projection methods that conserve energy but violate mass show quadratic growth, which identifies the missing linear invariant as essential in practice (Ranocha et al., 2021).

An important corrective nuance is that physical non-conservativity does not imply numerical non-conservative error. For the non-conservative nonlinear Schrödinger equation

$B(U)\partial_x U$1

mass and energy are not conserved, but they obey explicit balance laws,

$B(U)\partial_x U$2

and

$B(U)\partial_x U$3

The relevant numerical objective is therefore not discrete conservation but exact discrete analogues of these balance laws. Both the linearly implicit relaxation scheme and the modified Delfour–Fortin–Payre scheme satisfy discrete mass and energy balance identities and are second-order accurate in time; in this setting, non-conservative error is the deviation of discrete mass and energy from the exact physically evolving values, not a failure to conserve an invariant that the PDE itself does not possess (Athanassoulis et al., 2024).

Self-gravitating Eulerian hydrodynamics provides a final contrast. The traditional source-term treatment of gravity is physically conservative at the continuum level but non-conservative at the discrete level. Rewriting the change in total energy

$B(U)\partial_x U$4

in flux-divergence form via the gravitational energy flux

$B(U)\partial_x U$5

eliminates the secular total-energy drift of the source-term scheme and conserves total energy to round-off error. The paper shows that the non-conservative algorithm can materially change dynamics in collapse problems, including symmetry breaking and spurious rotation, especially when results are sensitive to small energy errors (Jiang et al., 2012).

4. Variational and dynamical-systems viewpoints

Outside PDE numerics, non-conservative error is also defined relative to the gap between a generalized invariant and the usual mechanical energy. In the generalized Lagrange formalism for stationary first-order dynamical systems,

$B(U)\partial_x U$6

with a Lagrangian linear in velocities,

$B(U)\partial_x U$7

time symmetry still yields a conserved quantity, $B(U)\partial_x U$8, even in non-conservative cases. The paper then defines non-conservative error in two complementary ways: first as the mismatch $B(U)\partial_x U$9, between the generalized time-symmetry integral and the usual mechanical energy $G=F+R$0; and second as the rate

$G=F+R$1

with generalized Poisson bracket $G=F+R$2. For the damped harmonic oscillator, $G=F+R$3, this reduces to the familiar dissipation law

$G=F+R$4

while $G=F+R$5 remains conserved. The resulting interpretation is that non-conservative error quantifies departure from classical energy behavior, not the absence of a conserved quantity altogether (Ushveridze, 2022).

A different dynamical meaning appears in accelerated flows under non-conservative vector fields. For Nesterov’s ODE driven by a field with Helmholtz decomposition $G=F+R$6, the destabilizing contribution to the Lyapunov derivative is identified explicitly as

$G=F+R$7

This term can dominate the damping and produce instability even for arbitrarily small non-conservative components in the linear case. The restart-based hybrid system removes the term at jumps by resetting $G=F+R$8, and explicit admissible reset periods

$G=F+R$9

yield uniform global exponential stability. Here non-conservative error is neither a flux defect nor an energy drift; it is the rotational cross-term that destroys Lyapunov descent under dynamic damping (Ochoa et al., 22 Jan 2025).

These two viewpoints share a structural message. A non-conservative system may still possess a conserved quantity generated by time symmetry, but that quantity need not coincide with mechanical energy, and rotational or skew-symmetric components can destabilize dynamics even when they are small. This suggests that “non-conservative” is best understood relative to the quantity one is trying to preserve: mechanical energy, a Noether integral, a Lyapunov function, or a weak-solution jump condition (Ushveridze, 2022, Ochoa et al., 22 Jan 2025).

5. Physics-informed learning, path consistency, and shock selection

Physics-informed neural networks reproduce the conservative-versus-non-conservative distinction almost exactly. For Burgers and compressible Euler, one can train on either the conservative residual

$p\,\partial_x A$0

or the non-conservative residual $p\,\partial_x A$1. In smooth regimes the two are equivalent; the Burgers experiment with smooth Gaussian initial data confirms that conservative and non-conservative formulations yield indistinguishable numerical and PINN solutions. At shocks, however, the non-conservative residual does not encode the correct distributional jump condition, and standard non-conservative numerical solvers can fail to propagate shocks or propagate them with severe smearing. The Adaptive Weight and Viscosity architecture introduces a trainable artificial viscosity $p\,\partial_x A$2, a viscosity penalty $p\,\partial_x A$3, and adaptive residual weighting such as $p\,\partial_x A$4, which stabilizes training and can recover correct shock propagation in Burgers and good shock resolution in 1D and 2D Euler even when training on primitive-variable residuals (Neelan et al., 27 Jun 2025).

That empirical success is not universal. A later analysis isolates the mechanism of failure for unsteady systems such as the Sod shock tube. If the primitive-variable non-conservative system is regularized as

$p\,\partial_x A$5

then transforming back to conservative variables introduces the extra source

$p\,\partial_x A$6

Inside a shock layer this source does not vanish and shifts the jump condition away from the classical Rankine–Hugoniot relation, producing wrong shock speeds. The paper documents that standard non-conservative PINNs can have small early-time errors in the Sod problem at $p\,\partial_x A$7 s yet show incorrect shock position by $p\,\partial_x A$8 s, whereas a path-integral PINN that adds the loss

$p\,\partial_x A$9

recovers the correct shock location (Neelan et al., 2 Apr 2026).

The same theme appears in the modified conservative PINN for the generalized Buckley–Leverett equation with discontinuous porosity. The original scalar law in non-conservative form, $A$0, is recast into the two-by-two resonant conservative system

$A$1

and the cPINN uses two subdomain networks: one for the standing $A$2-wave associated with $A$3, and one for the scalar $A$4-wave with an Oleinik-aware residual. Flux matching across the interface and rescaling near critical states suppress the oscillations that otherwise appear in non-conservative treatments near vacuum. The reported solutions are identical in conservative and non-conservative variables once the interface structure and entropy constraints are enforced (Quita et al., 2023).

A recurring misconception is that non-conservative PINNs fail only because of optimization. The cited papers distinguish two different issues. One is training stiffness, which AWV can mitigate. The other is structural inconsistency of the residual near discontinuities, which requires path-consistent or entropy-consistent constraints. The latter cannot be repaired by optimization alone when the learned residual enforces the wrong weak form (Neelan et al., 27 Jun 2025, Neelan et al., 2 Apr 2026, Quita et al., 2023).

6. Data-driven decomposition, drifting models, and recurring remedies

In data-driven mechanics, non-conservative error can be treated as an identifiable component of the force field. The Neural New-Physics Detector decomposes the observed force into

$A$5

with the conservative part represented by a Lagrangian Neural Network and the non-conservative part by a universal approximator network. Its objective

$A$6

combines force recovery error with a penalty on the magnitude of the predicted non-conservative component. The paper proves a universal phase transition at $A$7: for $A$8, the model recovers the exact decomposition with $A$9; for $\Xi$0, the non-conservative component is suppressed and the recovery error equals the magnitude of the true non-conservative part. Empirically this allows friction to be rediscovered from a damped double pendulum, Neptune from Uranus’ orbit, and gravitational waves from an inspiraling orbit (Liu et al., 2021).

An analogous distinction appears in one-step generative drifting. The conservative drifting velocity is a gradient field,

$\Xi$1

whereas the original displacement-based Laplace drift is non-conservative. The finite-particle analysis shows that the non-conservative Laplace method decomposes into a sharp-score mismatch plus an unavoidable scale-mismatch residual,

$\Xi$2

This residual does not disappear merely by increasing particle count $\Xi$3; it vanishes only when the local Laplace scales of data and model align. By contrast, the conservative KDE-gradient method yields finite-particle bounds driven by a reciprocal-KDE self-interaction term and quadrature errors, with root residual-velocity rate $\Xi$4 under an $\Xi$5-uniform quadrature regularity condition (Balasubramanian, 21 May 2026).

Across these literatures, the remedies are strikingly consistent. One class rewrites the problem in a genuinely conservative or balance-preserving form: conservative energy fluxes for self-gravity, global-flux reconstructions, resonant conservative augmentations, or score-gradient drifts (Jiang et al., 2012, Ciallella et al., 1 Jul 2025, Quita et al., 2023, Balasubramanian, 21 May 2026). A second class leaves the non-conservative formulation in place but adds a consistency mechanism: DLM path integrals, local energy-residual corrections, subcell finite-volume corrections, a posteriori level-set shifts, restart maps, or explicit penalties on the non-conservative component (Abgrall et al., 2017, Gaburro et al., 2024, Long et al., 2022, Ochoa et al., 22 Jan 2025, Liu et al., 2021). A plausible implication is that non-conservative error is best controlled not by generic stabilization alone but by restoring the exact structural object that the formulation has lost: a flux, a path, a balance law, a symmetry-generated invariant, or a gradient representation.