# Lax-Wendroff method for linear advection - Stability analysis

Question 1: Consider the wave equation $$u_t + c(x) u_x = 0 ,$$ where $x\in \Omega \subset \Bbb R$ and $c(x)$ is a function of $x$.

(a) Show that the Lax-Wendroff scheme for this PDE is given by $$u_j^{n+1} = u_j^n - c_j \Delta t \frac{D_x u_j^n}{2 \Delta x} + \frac{c_j^2 \Delta t^2}{2} \frac{\delta_x^2 u_j^n}{\Delta x^2} + \frac{c_j \Delta t^2}{8 \Delta x ^2} (D_x c_j)(D_x u_j^n) ,$$ where $D_x$ is the first central difference operator, $\delta_x^2$ is the second central difference operator, and $\Delta t$ and $\Delta x$ are the mesh-spacing in $t$ and $x$, respectively. The $j$ and $n$ are space and time indices, respectively, and $u_j^n$ is the grid function such that $u_j^n\approx u(x_j,t_n)$ and $c_j \approx c(x_j)$.

I would also like to answer the following question

Assuming $c(x)$ is a constant, and given that the initial value problem is well posed, carry out a von Neumann Stability analysis and hence show that the scheme is convergent provided $$\bigg \lvert \frac{c\Delta t}{ \Delta x} \bigg \rvert \leq1$$

I have tried to do the stability analysis but I am not getting anywhere near what the question requires. Could someone please check my working and provide some advice.

The scheme can be written as

$$u_{j}^{n+1} = u_j^n - \frac{p}{2}(u_{j+1}-u_{j-1}) + \frac{p^2}{2}(u_{j+1}^n-2u_j^n+u_{j-1}^n) + \frac{p^2}{8}(u_{j+1}-u_{j-1}), \quad \text{where p = \frac{c\Delta t}{5 \Delta x}}$$

I think this is how the scheme should be written. However, I am not sure about the $D_{x}c_{j}$ term, since it is a constant, should I just ignore the $D_x$ operator?

For the stability analysis, using a trial solution of the form $u_{j}^{n} = A \xi^{n}e^{i \omega j}$ gives

$$A \xi^{n+1} e^{i \omega j} = A \xi^{n} e^{i \omega j} - \frac{p}{2}\left[ A \xi^{n} e^{i \omega (j+1)} - A \xi^{n} e^{i \omega (j-1)} \right] + \frac{p^{2}}{2} \left[ A \xi^{n} e^{i \omega (j+1)} - 2 A \xi^{n} e^{i \omega j} + A \xi^{n} e^{i \omega (j-1)} \right] + \frac{p^2}{8}\left[ A \xi^{n} e^{i \omega (j+1)} + A \xi^{n} e^{i \omega (j-1)} \right]$$

Dividing through and by $u_j^n = A \xi^{n} e^{i \omega j}$ I get

$$\xi = 1 + \frac{p}{2}\left[ e^{i\omega} - e^{-i\omega} \right] - \frac{p^{2}}{2} \left[ e^{i\omega} - 2 + e^{-i\omega} \right] + \frac{p^{2}}{8} \left[ e^{i \omega} + e^{-i \omega} \right] = 0$$

Than, using Eulers identity

$$\implies \xi = 1 + \frac{p^{2}}{2}\cdot 2 \cos \omega - \frac{p^{2}}{2} (2 \cos \omega - 2) - \frac{p^{2}}{8}\cdot 2 \cos \omega$$

Can anyone help me from here?

• I tried to fix your post, but there were so many typos and errors that I'm not even sure I have fixed it correctly. Can you check to make sure the post is as you wanted it. – mattos Apr 19 '18 at 0:42
• Also, $e^{i \omega} - e^{-i \omega} = 2i \sin \omega$ not $2 \cos \omega$. – mattos Apr 19 '18 at 0:43
• Yes u edited it correctly, thank you so much! – italy Apr 19 '18 at 1:02

Recall how the Lax-Wendroff method is obtained in the constant-speed case :

• a Taylor series in time is written: $$u(x,t_{n+1}) = u(x,t_{n}) + \Delta t\, u_t(x,t_{n}) + \frac{1}{2}\Delta t^2\, u_{tt}(x,t_{n}) + \dots$$
• the time derivatives are eliminated using the PDE: $u_t = -c u_x$ and $u_{tt} = c^2 u_{xx}$.
• the spatial derivatives are replaced by central finite-difference approximations.

Thus, the following scheme is obtained: $$u_j^{n+1} = u_j^n - c\Delta t \frac{D_x u_j^n}{2 \Delta x} + \frac{c^2\Delta t^2}{2} \frac{\delta_x^2 u_j^n}{\Delta x^2} \, .$$ This method can be adapted to the variable velocity case.

Let us analyze the stability of the Lax-Wendroff scheme above (constant-speed case). Assuming a perturbation of the form $u_j^{n} = \xi^n \text{e}^{\text i k x_j}$, one has \begin{aligned} \xi &= 1 - c \Delta t \frac{\text{e}^{\text i k \Delta x} - \text{e}^{-\text i k \Delta x}}{2 \Delta x} + \frac{c^2 \Delta t^2}{2} \frac{\text{e}^{\text i k \Delta x} - 2 + \text{e}^{-\text i k \Delta x}}{\Delta x^2} \\ &= 1 - \text{i} \kappa \sin (k\Delta x) + \kappa^2 \left(\cos (k\Delta x) - 1\right) . \end{aligned} where $\kappa = c\frac{\Delta t}{\Delta x}$ is the Courant number. Thus, the squared modulus of the amplification factor is \begin{aligned} |\xi|^2 &= (1 - \kappa^2 (1-\cos(k\Delta x)))^2 + \kappa^2 (1-\cos^2(k\Delta x))\\ &= 1 - 2\kappa^2 (1-\cos(k\Delta x)) + \kappa^2 (1-\cos^2(k\Delta x)) + \kappa^4 (1-\cos(k\Delta x))^2 \\ &= 1 - \kappa^2 (1-\kappa^2) (1-\cos(k\Delta x))^2 \\ &= 1 - 4\kappa^2 (1-\kappa^2) \sin^4\left(\tfrac{1}{2}k\Delta x\right) . \end{aligned} Finally, the Lax-Wendroff scheme is Neumann-stable provided that $|\xi|^2 \leq 1$, which implies that the CFL condition $|\kappa| \leq 1$ is satisfied.

 R.J. LeVeque, Finite Volume Methods for Hyperbolic Problems, Cambridge University Press, 2002.

• Do you ignore the last term while doing the stability analysis ? – italy Apr 24 '18 at 15:48
• Oh okay, thank you – italy Apr 24 '18 at 16:17
• I am trying the get this amplification factor but I don’t seem to get the same answer... $1 -pi sin \omega - p^{2} (1-cos \omega )$ from here I know I would have to take the complex conjugate was of this, $\xi ^{2} = 1 +p^{2}(1-cos^{2} \omega)-2p^{2}(1-cos \omega) p^{4}(1-cos \omega)^{2}$ I have tried to simplify this but I get anything meaningful can you suggest anything from here? Thank you – italy Apr 24 '18 at 19:07
• The amplification factor I am looking for should be $\xi = 1-2p^{2} sin^{2} \frac{\omega}{2} -2ip sin \frac{\omega}{2} cos\frac{\omega}{2}$ – italy Apr 24 '18 at 19:16
• @italy Thanks for the catch: the square was not at the right place in my answer. Now, a few intermediate steps have been added. Note that $\cos u = \cos 2\frac{u}{2} = 1 - 2\sin^2 \frac{u}{2}$ and $\sin u = \sin 2\frac{u}{2} = 2\sin \frac{u}{2}\cos \frac{u}{2}$. – Harry49 Apr 24 '18 at 22:23