# Numerical scheme with energy conservation?

I have a set of equations to integrate something in time $t$. At each time step I compute a scalar field $\phi(t)$ and a potential $V(\phi)$. I should also control the conservation of energy with an equation:

$E=\frac{1}{2}\left(\dot\phi^2+V(\phi)\right)$

Or numerically:

$E=\frac{1}{2}\left(\left(\frac{\phi_{n}-\phi_{n-1}}{t_{n}-t_{n-1}}\right)^2+V(\phi)\right)$

My question is: is it $V(\phi_{n})$ or $V(\phi_{n-1})$ ?

• It should be: $E=\frac{1}{2}\left(\left(\frac{\phi_{n+1/2}-\phi_{n-1/2}}{t_{n}-t_{n-1}}\right)^2+V(\phi_n)\right)$ – fibonatic Nov 16 '13 at 16:21
• The answer depends. Are you using $E$ to passively measure nonconservation? Then it doesn't matter. Are you using $E$ to control your integration? Then any formula works, but you should be aware of 1st-order vs. 2nd-order accuracy of the formulas. – user3224 Nov 16 '13 at 17:21
• You may also be interested in symplectic integrators - if you're solving a Hamiltonian system, these can tend to conserve key features of that system. en.wikipedia.org/wiki/Symplectic_integrator. See also energy drift for a discussion of where time integrators tend to be useful (en.wikipedia.org/wiki/Energy_drift) – Jesse Chan Nov 17 '13 at 21:28

It depends on how you want to solve it: via implicit or explicit methods.

The explicit method is certainly the easiest to code but constrains $dt$ for stability; the implicit method is a bit more difficult to code but has no such constraints on $dt$.

For example, if we consider a simple advection equation $$\frac{\partial w}{\partial t}+u\frac{\partial w}{\partial x}=0$$ Then the explicit method solves this as (using your notation, I think most people do $n\to n+1$, rather than $n-1\to n$) $$\frac{w^{n}_j+w^{n-1}_j}{dt}=u\left[\frac{w_{j+1}^{n-1}-w_{j-1}^{n-1}}{2\cdot dx}\right]$$ which is more commonly written as $$w^n_j=w^{n-1}_j+\chi u\left[w^{n-1}_{j+1}-w_{j-1}^{n-1}\right]$$ where $\chi=dt/2dx$. The implicit method solves this as $$\frac{w^{n}_j+w^{n-1}_j}{dt}=u\left[\frac{w_{j+1}^{n}-w_{j-1}^{n}}{2\cdot dx}\right]$$ which is more commonly written as $$-\chi uw_{j+1}^n+w^n_j-\chi uw^n_{j-1}=-w_j^{n-1}$$ which requires solving a tri-diagonal matrix (fortunately there is a relatively simple algorithm to solve this).

I generally would recommend trying the explicit method first and seeing if the time-stepping is good enough to be useful; if it's not, then it's time to use the implicit method.

I will just sum up the useful remarks made in the comments. Both approximations

$E=\frac{1}{2}\left(\left(\frac{\phi_{n}-\phi_{n-1}}{t_{n}-t_{n-1}}\right)^2+V(\phi_n)\right)$

and

$E=\frac{1}{2}\left(\left(\frac{\phi_{n}-\phi_{n-1}}{t_{n}-t_{n-1}}\right)^2+V(\phi_{n-1})\right)$

will give you first-order accurate approximation of the energy. If you wish to get higher than 1st order accuracy, you should use a more accurate finite difference approximation of $\dot\phi$. For instance, you could get second order with

$E=\frac{1}{2}\left(\left(\frac{\phi_{n+1/2}-\phi_{n-1/2}}{t_{n}-t_{n-1}}\right)^2+V(\phi_{n})\right)$.

The definition of the derivative is $$\dot{\phi}_n=\lim_{\delta t\rightarrow 0}\frac{\phi_{n+1}-\phi_n}{\delta t}$$ where $\delta t$ is $t_n-t_{n-1}$. Therefore, the derivative as you are taking it is calculated for the $n-1$-th step. Hence, I'd say you should go for $V(\phi_{n-1})$.