Dissipation

Roughly speaking, numerically-induced energy loss. If it is worse that the dissipation inherent in the PDE, it will lead to incompatible approximations to the real solution. The requirement that a scheme have $\vert g(h\xi)\vert\le 1$ ensures stability. If $\vert g(h\xi)\vert=1$ for all $h\xi$, it means that the amplitude of the mode $\xi$ is not affected by time-stepping: it does not get diminished. However, if $\vert h(\xi)\vert<1$ for some or all $h\xi$, a slight or large loss in amplitude is incurred. This is dissipatioin (numerical) and is an artifice of the scheme rather than of the equation being approximated.

For evolutionary PDE's, as can be imagined, if a scheme has more numerical dissipation than is inherent in the PDE being approximated, the solution will eventually be different in amplitude and phase to the approximation solution provided by a scheme with less dissipation. The phase phenomenom must be considered too, since the scheme will in general have different rates of dissipation for different modes.

In order to measure to dissipation concretely, we have to agree on a definition of dissipation. This is one possibility:

definition: Let $\theta= \xi h$. A scheme is dissipative of order $2r$ if there exists a positive contant $c$, independent of $h$ and $k$, such that each amplification factor $g_{\nu}(\theta)$ satisfies

$$\vert g(\theta)\vert^2\le 1-c\sin^{2r}\left(\frac{1}{2}\theta\right)$$

Other definitions will replace $\sin^{2r}\left(\displaystyle \frac12\theta\right)$ by $\vert\theta\vert^{2r}$.

Example

This happens at boundary conditions. Suppose we are solving a problem over a space grid indexed by $m=0, 1\cdots M$. Let $u^n_{m}$ be an approximation to $U(x_m, t_n)$. Hence the edge variables are $u^n_0$ and $u^n_M$.

Lax-Wendroff:

$$\vert g(\theta)\vert^2=1-4a^{e2}\lambda^2(1-a^2\lambda^2)\sin^4\frac12\theta$$

For $\vert a\lambda\vert=1$ we have $\vert g(\theta)\vert=1$, nondissipative. But for $0<\vert a\lambda\vert<1$ the scheme is of order 4 in dissipation.

Example: Show that

Leapfrog and Crank Nicholson $\rightarrow$ both non-dissipative since their amplificaiton factors are identically 1 in magnitude.

Example

Lax-Friedrichs: show that $\vert g(\theta)\vert=1$ for $\theta=0$ and $\pi$, but less than 1 for their values $\Rightarrow$ dissipative.

Remark: Sometimes dissipation is good. It may also be added to schemes in order to stabilize them. For example, adding

Dissipation may also be added to schemes in order to stabilize them. For example, adding

$$\frac{\varepsilon }{2k}\left(\frac{1}{2}h\delta\right)^4v_m^{n-1}$$

where $\delta^2v_{m}^{n-1}\equiv v_{m+1}-2v_m+v_{m=1}$

and $\varepsilon \ll 1$ leads to

$$g_{\pm}=-a\lambda\sin\theta\pm\sqrt{1-a^2\lambda^2\sin^2\theta-\varepsilon \sin^4\frac12\theta}$$

$$\mbox {if }\varepsilon <1-a^2\lambda^2 \mbox {scheme is stable and of } {\cal O}(4) \mbox { in dissipation.}$$

Dispersion There are PDE's that have dispersive terms (KdV, Nonlinear Schrodinger Equation, etc). In hyperbolic problems, these dispersive terms force each Fourier mode to travel at different speeds. Hence, if a wave that was compact at some time is subjected to dispersion (and is not balanced by other other effects, such as could be possible with nonlinearity or dissipation), will eventually spread out in space and time. An example of a real physical system in which dispersive effects are readily observed is: throw a rock into a pond and watch the concentric waves propagate out of the center of impact. Far from the center you see that waves of different wavelengths will separate. This would not happen if the surface of the lake, which is capable of supporting waves, were not dispersive. In the absence of dispersion the initial disturbance set up by the rock would propagate out as a single and compact ring of waves.

Dispersion can also be caused unwittingly by certain numerical approximations to equations. If it is unwanted, it is a form of distortion and it turns out a fairly important one. Suppose we were solving the one-way wave equation with constant speed. In Figure 25, which would correspond to the approximation with a numerical method with zero dispersion, we have the initial data, which can be thought of as a superposition of waves (via Fourier methods) of wavenumber $\kappa$ and corresponding frequency $\omega$ all traveling at speed $c$, constant. Hence, the dispersion relation $\omega = \kappa c$, where $c$ is constant. In Figure 26 we would have the same initial data with each wave of component traveling at speed $c(\kappa)$ and the initial data would then spread and distort.

Figure 25: Non-dispersive approximation to the one-way wave equation with wave speed $c$, constant. The dispersion relation is a straight line with slope $c$

Figure 26: Dispersive approximation to the one-way wave equation with wave speed $c$ no longer a constant, hence the dispersion relation is $\omega = c(\kappa )$. The local slope of this function below is the speed, which shows that the higher frequencies are traveling faster than the slower ones.

$c$ is fixed constant speed. $k$ is the wave number and $\omega$ the frequency.

Take

$$U(t,x)=\frac{1}{\sqrt{2\pi}}\int^{\infty}_{-\infty}e^{i\omega x} e^{-i\omega at}\hat U_0(\omega)d\omega$$ (131)

This is a solution to

$$U_t+aU_x =$$ 0

$$U(x,0) = U_0(x)$$

Here $a$ is constant.

From (132) we conclude that

$$\hat U(t+k, \omega)=e^{-i\omega k}\widehat U(t,\omega).$$ (132)

A one-step finite difference scheme gives

$$\hat U^{n+1}=g(h\xi)\hat U^n$$ (133)

Comparing (133) and (134) we expect $g(h\xi)$ to be a good approximation to $e^{-1\xi ak}$. Let's write $g(h\xi)$ in term s of magnitude and phase:

$$g(h\xi)=\vert g(h\xi)\vert e^{-i \xi\alpha(h\xi)}k$$

``phase speed'' is the speed at which waves of wave number $\xi$ are propagated by the FD scheme.

If $\alpha(h\xi)=a$, constant, for all $\xi\Rightarrow$ waves would propagate at correct speed $\ldots$ thus is what we expect from analytical solution of $u_t+au_x=0$.

The ``phase error'' could be a measured by $a-\alpha(h\xi)$.

When the waves travel at different speeds, we say the FD scheme is dispersive and can use the phase error to assess how badly this affects the solution. When the waves travel at constant speed $\alpha_0=k'a$, where $k'$ is a constant and $a$ is the wave speed in $u_t+au_x=0$, we say the FD scheme is non-dispersive.

Example

Lax-Wendroff: $g=1-2(a\lambda)^2\sin^2\frac{1}{2}h\xi ia\lambda\sin h\xi$

$$\mbox {and so take }[\alpha(h\xi)\xi k]=\frac{a\lambda\sin h\xi} {1-2(a\lambda)^2\sin^2\frac12 h\xi}$$

But we want $$\frac{\tan^{-1}}{\xi k}\left[\displaystyle \frac{a\lambda\sin h\xi} {1-2(a\lambda)^2\sin^2\frac12h\xi}\right]=\alpha(h\xi)$$

Take the low wave number limit $\xi\rightarrow 0$

$$\alpha(h\xi)=a\left\{1-\frac{1}{6}(h\xi)^2 [1-(a\lambda)^2]+0(h\xi)^4\right\}$$

So for $h\xi$ small and $(a\lambda)<1, \alpha(h\xi)<a$. Also if $\vert a\lambda\vert\rightarrow 1$ then the dispersion is smaller.

Example Find dispersion behavior for larger $\xi\rightarrow h^{-1}\pi$.

Remarks

(1)

In general, for hyperbolic problems, we usually want to take $\vert a\lambda\vert$ close to stability limit: usually gives largest time steps and in general, for dissipative schemes, the least dissipation. But more importantly, most likely would yield the smallest dissipation and disperse errors.

(2)

The comments here apply broadly to all evolutionary PDE approximations.

(3)

The leap-frog is an example of a scheme with no dispersion error when $a\lambda=1$. It doesn't have dissipation either when $a\lambda=1$. However, the leap-frog method have a couple of notorious problems: a) it can have bad stability problems when coupled to certain boundary conditions. b) being that it is a 2nd order-in-line method, its approximation is made of 2 traveling wave solutions, traveling in opposite directions, generally. However, if the hyperbolic problem being approximated only admits 1-way wave solutions, a single wave, care must be exercised in either not exciting the spurious solution or in actively suppressing it $\ldots$ one popular suppressing technique is to revert to a single Euler every 100's or 1000's time steps and then reverting back to Leapfrog (doing too many Euler will defeat the purpose of using leapfrog and will generate a lot of dissipation). There are other techniques for this, such as using time- and/or space-averaging filters.

Group Velocity and Propagation of Wave Packets

The group velocity is the speed at which the energy in a wave packet travels at. It is a useful concept in nonlinear and dispersive equations. It can be used to explain some rather striking behavior of certain schemes, including the explanation of certain instabilities caused by boundary conditions (see Trefethen).

We've seen that dispersive FD schemes will cause a pure wave with wave number $\xi_0$ to travel with phase speed $\alpha(h\xi_0)$. We want to know what is the velocity of propagation of the center of mass of a wave packet.

The scheme group velocity is

$$\gamma(\theta)\equiv\frac{d}{d\theta}(\theta\alpha(\theta)),$$   where $$\quad \theta=h\xi$$

A wave packet example suppose initial data of the form

$$u(0,x)=e^{i\xi_0x}p(x)$$

where $p(x)$ is a relatively smooth function decaying rapidly about its centor of mass. For

$$u_t+au_x=0$$

the solution is $u(t,x)=e^{i\xi_0(x-at)}p(x-at)$

$p(x)$ is the envelope of the wave packet and $e^{i\xi_0x}$ is the carrier wave, see Figure 27

Figure 27: Wave Packet

For a finite difference approximation below, $\gamma$ is the group velocity:

$$u(t,x)=e^{i\xi_0(x-\alpha(h\xi_0)t)}p(x-\gamma(h\xi_0)t)$$

(in the case of zero dissipation). Note that since $\alpha(h\xi)\rightarrow a$ as $h\rightarrow 0$ we have that $\gamma(h\xi_0)\rightarrow a$ as $h\rightarrow 0$ and $u\rightarrow U$, the exact solution.