Say I'd like to implement a numerical model of a physical system (e.g. a simple wave equation). Of course the energy conservation law is valid for this system. Is there a way to validate that my implementation of the wave equation is implemented in accordance to the law of energy conservation in general?

I do not want to implement some test scenarios and check whether the energy is constant at all time. Instead I'd like to use some sort of symbolical solver to verify this property of my implementation. Something like Cryptol but instead for physics and of course GPL.

  • $\begingroup$ This question was migrated from Software Quality and Assurance Stack Exchange. $\endgroup$ Oct 7, 2013 at 23:21

2 Answers 2


The typical way ODE and PDE numerical methods are tested is using something like the method of manufactured solutions; a more extensive report on validating large physics codes can be found here.

For an ODE solver, the basic idea for testing a method for solving

\begin{equation} \dot{u}(t) = f(u(t)) \end{equation}

would be to pick a solution, say $u(t) = e^{-t}$, and then calculate what the source term $f$ should be using a symbolic algebra system. For this example, $f(u) = -e^{-t} = -u$. You can also calculate initial conditions and boundary conditions by direct substitution. Here, $u(0) = 1$.

There are GPL- and BSD-licensed symbolic manipulation programs out there that will take care of the algebra for you. In Python, for instance, there Sage and Sympy can each do symbolic manipulation.

Even though you say you'd prefer not to, testing explicitly for things like symmetry or energy conservation (to within a tolerance) is something also done in practice.

Alternative approaches exist, but they're typically impractical. For instance, you could use a computer assisted proof system like HOL Light or Coq to reason about your program, but the computer assisted proof would take an impractically long time to execute. The Flyspeck project aims to do a computer-assisted proof of Kepler's conjecture, and the computation is estimated to take roughly 2.5 years, which is longer in serial than some "hero" computations in high-performance computing. You could also use something like interval arithmetic to ensure that your results lie within a given range, but things like dependency effects would likely make the output intervals very large relative to the precision you want, and you could suffer from a lack of infrastructure (i.e., there aren't many libraries that do interval arithmetic, and these libraries don't necessarily implement enough of what an average nonspecialist might need to do, like, say, solve systems of nonlinear equations).

I gather from your remark about Cryptol that you might want to do some sort of computer-assisted proof. Based on the above, I wouldn't recommend it, but if it could be done efficiently, that would be a very interesting contribution to the numerical methods literature.


If your differential equations support energy conservation in some form then you should be able to prove it analytically first. So you write the energy of the system as some kind of volume integral, and take time derivative of it, and use your differential equation for time-evolution of involved quantities (a wave equation in your case) to show that the time-derivative of total energy is identically zero - that would be proving that the model obeys energy conservation. But even if it does, it is a separate issue if your discretized in space and time equations solved on a computer still satisfy this conservation law to machine accuracy (look up conservative finite-difference or finite-volume methods).

Of course the equations representing physical laws should be consistent with energy conservation but in formulating a numerical model there are all kinds of approximations made, so the actual equations in the model are by no means guaranteed to satisfy the "real" conservation laws. For example, you model the motion of a body in viscous medium, so there is a friction force $F = - k V$. You can write and solve the equation of motion for the body, $m \dot{V}=- k V$, but obviously there is no energy conservation here since the thermal energy of the medium is not in the model. This does not mean that the model of body motion is not "right", it is perfectly correct for calculating the body trajectory.

On the other hand, it may happen that a model may satisfy conservation of something that is not a full physical energy but some kind of proxy of it. A good example is the "energy" in plasma turbulence, $E = \int ( {n}^2 + |\nabla \phi|^2 ) dV$, where $n, \phi$ are plasma density and electric potential fluctuations, see papers on Hasegawa-Wakatani drift turbulence, e.g., http://iri.columbia.edu/~suzana/papers/camargo_biskamp_scott95.pdf.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.