I try to write tests to the implementations of numerical methods. The best way to do this is to study the observed order of accuracy and check that it matches with the theoretical order of accuracy.

However, I noticed that the observed order of accuracy lies in a band $$ [p_{\mathrm{theor}} - \Delta, p_{\mathrm{theor}} + \Delta], $$ where $p_{\mathrm{theor}}$ is the theoretical order of accuracy and $\Delta$ is some absolute error that must be chosen for any particular method. For some of the methods I implement I can take $\Delta=0.05$, for other $\Delta=0.1$ and so forth.

So, my question is: how should one choose $\Delta$ to be sure that there are no mistakes in the implementation?


From practical experience, I find plotting the results of discretization refinement studies best to figure out whether or not the observed order of accuracy is "close enough". However, writing an automated test of a method's order of accuracy requires some sort of regression to calculate the observed order of accuracy from data, and then defining an appropriate $\Delta$ -- a comparison tolerance -- such that if the observed order of accuracy is within $\Delta$ of the theoretical order of accuracy, the test "passes".

I don't know of any theory that suggests an appropriate $\Delta$ or an appropriate way to bracket the asymptotic order of convergence region a priori. I have always written these kinds of tests by reverse-engineering them from plots, so I typically choose $\Delta = .2$ or $\Delta = .3$, because I see the order of accuracy deviate by at most that amount.

In certain cases (e.g., discontinuous initial conditions, or stiff problems), theory can predict how the order of accuracy of a method will degrade; if this additional theory isn't taken into account, then you might need to choose a larger $\Delta$ to compensate. For instance, in Section 8.7 of "Finite-Volume Methods for Hyperbolic Problems", Randy LeVeque discusses how convergence theory for a first-order upwind method assumes smooth solutions. If instead, a weak solution is discontinuous, the order of accuracy in space (in the $L^{1}$ norm) is one-half instead of one. So in a case like that, if you assume that the true order should be one, you might need to make $\Delta$ as large as 0.6, but in some cases, the degradation of the error is known, and you can adjust the theoretical order of accuracy accordingly.

  • 2
    $\begingroup$ "I don't know of any theory that suggests..." -- in reality, it's also not quite clear of course how you compute the convergence rate. You could use only the last two data points, or fit to all of them. The point being that the convergence rate of an implementation is an asymptotic statement that can not be evaluated in practice, only approximated by a series of data points. $\endgroup$ – Wolfgang Bangerth Sep 4 '15 at 21:55
  • $\begingroup$ Agreed; that's why the sentence you cite ends with "...or an appropriate way to bracket the asymptotic order of convergence region a priori." Indeed, one could use as many (or as few) points as one sees fit, provided at least two points are used, and the points are plausibly in the asymptotic region (which I would check using a plot). $\endgroup$ – Geoff Oxberry Sep 5 '15 at 1:30

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.