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?
Best 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.