# How should a numerical solver treat conserved quantities?

Since the good old Noether Theorem ultimately states that any physical system will exhibit conserved quantities, how should they be treated best in numerical solvers?

On the one hand, one can observe the numerical propagation of such quantities and use their deviation from conservation as indicator on the solution's quality, e.g in order to adapt a time step or mesh size.

On the other hand, one can try and craft a solver such that the conserved quantities are implicitly conserved (under more or less all circumstances). That takes away one indicator on the solution's quality, but there seems to be a justified hope that in exchange the solver is simply "better".

Are there any more detailed discussions on which of these possibilities is "better" and why?

• If you "enforce" conservation, there is at least hope that the numerical solution is a "shadow" solution, i.e. there exist slightly modified initial conditions such that their exact trajectory stays close to the numerical solution forever. – Hagen von Eitzen Mar 11 '13 at 16:37
• @HagenvonEitzen Not necessarily - take the simple harmonic oscillator $\ddot x + x =0$ plus a timestep too large and you'll end up with a wrong frequency, getting absolutely out of sync with the real solution. And if you increase the timestep even further, you'll even end up with imaginary coordinates/velocities... – Tobias Kienzler Mar 12 '13 at 5:59
• OK, but that would be close to an exact solution of $\ddot x+ax=0$ with $a\approx 1$. So not only initial conditions, but also other constants might be "slightly off". – Hagen von Eitzen Mar 12 '13 at 7:21
• @HagenvonEitzen Yes, exactly. And thereby indirectly other equations relating conserved quantities to these "constants", e.g. the relation between frequency and energy is then also off. So that might actually yield another indicator for the solution's proximity, i.e. a conservative scheme will still be able to provide a time step / mesh size adjuster, it's just less obvious than using the conserved quantities themselves – Tobias Kienzler Mar 12 '13 at 7:51
• @TobiasKienzler You can also stabilize constraints (or conservative quantities) with Baumgarte method. – AnilB Mar 14 '13 at 23:25

In numerical PDE, especially dealing with conservation laws $\partial_t u - \nabla \cdot\mathbf{F}(u) = g$, there are a class of numerical methods having local conservative properties. For example: Finite Volume Methods (FVM), and FVM's higher order counterpart Discontinuous Galerkin Methods, or Mixed Finite Element Methods.

• "local conservation of flux" of Finite Volume method: for simple stationary diffusion equation $-\Delta u = f$, the numerical solution $u_h$ satisfies: $$\int_{\partial K} \frac{\partial u_h}{\partial n} dS = \int_K f \tag{1}$$ where $K$ is a control volume of the original Delaunay triangulation. The control volume is the Voronoi tessellation dual to Delaunay triangulation (dotted line polygon in the following figure):

• If we add time in the equation: $\partial_t u - \Delta u = f$, the local conservation of flux (1) will still be satisfied at every time step if you finite volume method. Due to its simple implementation of this nice property, FVM is very popular in computational fluid dynamics dealing with incompressible flow. A reference of FVM would be this note.

• Discontinuous Galerkin methods are very popular in recent years among engineering applications. In that it is the generalization of FVM to higher order numerical approximation.

• Mixed Finite Element method (MFEM) embraces the other type of idea toward conservation of flux, for $-\Delta u = f$. MFEM decompose the second order PDE into a first order system: $$\begin{cases}\boldsymbol{\sigma} +\nabla u = 0 \\ \nabla \cdot \boldsymbol{\sigma} = f \end{cases}\tag{2}$$ and solve for the approximated flux $\boldsymbol{\sigma}$ in a discrete space directly, so that the numerical flux $\boldsymbol{\sigma}_h$ has weak conservation of flux, in terms of: $$\int_T \nabla \cdot \boldsymbol{\sigma} \phi = \int_T f \phi$$ in a simplex of the Delaunay triangulation where $\phi$ is a test function in the approximation space of $u$.

• Other numerical PDE methods like Continuous Galerkin Finite Element methods (CGFEM), Finite Difference Methods (FDM) do not have any of these kinds of conservation.

To sum up: FVM and some schemes in Discontinuous Galerkin methods address the conservation directly by solving the conservation equation in the approximation space. MFEM imposes conservation weakly. CGFEM doesn't have any form of conservation.

• Which one is better? There is not a single definitive answer. All is about trading different amenable and amiable aspects of each approach. DGM has local conservation of flux, but loses continuity and needs extra stabilization terms in the weak formulation. CGFEM has continuity and is stable, but no conservation. MFEM has somewhat both, but needs quite a lot extra computational resources because you have to solve an extra quantity.

• Which one to choose? For PDEs in which the discontinuity exists, and propagates at a finite speed, for example the linear transport along a flow field $\mathbf{a}$: $\partial_t u - \nabla\cdot (\mathbf{a}u) = 0$. Choosing a scheme that addresses conservation would be wiser opposing to choose something like CGFEM. For PDEs where the discontinuities propagate at an infinite speed, and is smoothened out immediately, like diffusion equation $\partial_t u - \nabla\cdot(\alpha \nabla u) = f$. Choosing CGFEM would be better in the sense that it is stable.