Measuring the convergence order of a numerical scheme for PDE

To obtain the convergence order of a numerical scheme, the following formula is used

$$R = \frac{ \log_2 \| e_{\text{new}} \| - \log_2 \| e_{\text{old}} \| }{ \log_2 \| \Delta x_{\text{new}} \| - \log_2 \| \Delta x_{\text{old}} \| }$$

where $$e$$ is the error in the approximation of the numerical solution of the PDE and the exact solution and $$\Delta x$$ is the mesh size in the spatial direction (see this post).

For simplicity, let us take $$v_t + v_x = 0$$ and suppose we want to investigate the convergence order using $$N=200$$ mesh points. Consider $$x \in [0,1]$$ and $$\phi(x) = u(x,0)$$ initial cauchy data. Say using FTBS scheme. The exact solution is $$u(x,t) = \phi(x-t)$$. Suppose we want to study it at $$t=t_{\text{max}}$$. So we have

$$e_{\text{old}} = | \text{approx}_{N=200} - F(x-t_{\text{max}}) |$$

where $$\text{approx}_{N= n }$$ means approximation using $$n$$ mesh intervals and $$\Delta x_{\text{old}} = 1/N = 1/200$$.

To compare it with need $$e_{\text{new}}$$, my question is

Is there a rule of thumb into what size of mesh intervals to pick for the $$e_{\text{new}}, \Delta x_{\text{old}}$$?

Say if we take $$M=100$$ then

$$e_{\text{new}} = |\text{approx}_{N=100} - F(x'-t_{\text{max}}) |$$

and $$\Delta x_{\text{new}} = 1/100$$. This gives a value, but if we choose $$M=199$$, say we have a negligible number. So, what is the rule of thumb to pick new number of mesh intervals?

• Usually people want to use these numbers to demonstrate numerically observed order of convergence. Often people like to half the mesh size so that, e.g., the first order of convergence would mean that the errors are decreasing with a factor of two. It’s just easier to divide by two or four or a power of two rather than by a factor of 199/200 Feb 17, 2019 at 18:41

1 Answer

For the case of $$v_t +v_x = 0$$, the time step is related to the mesh size via a fixed value of the Courant number $$\Delta t/\Delta x$$. Therefore, error is analyzed in terms of the mesh size only. As specified in the comments, people like to use powers of two for the number of points $$N$$. This gives the following type of error plot:

The $$x$$-axis is the mesh size $$\Delta x$$ in log scale, the $$y$$-axis is the error $$e$$ in log scale. If $$e = C (\Delta x)^R$$, then we have $$\log e = R \log\Delta x + \log C$$, i.e. the error plot is a line with slope $$R = \frac{\text d \log e}{\text d \log \Delta x}$$ in log-log coordinates. For example, one can use the $$L^2$$-error \begin{aligned} e = \| v_N - v_\text{ref}\| &= \sqrt{\int_0^1 (v_N(x) - v_\text{ref}(x))^2\, \text d x} \\ &\simeq \sqrt{ \Delta x \sum_{j=1}^{N-1} (v_j - v_\text{ref}(x_j))^2 } \end{aligned} but other choices are possible. On the figure, we can see that the slope of the error plot converges towards a fixed value as $$\Delta x \to 0$$, which is the order of convergence $$R \approx 1.57$$ measured numerically. To reach the theoretical value of the convergence order, it is important to perform the error measurements by using smooth signals! Then, the estimation is all the more precise as the mesh size is small. A too small mesh size could lead to round-off errors and too expensive computations. Therefore, a compromise is made.