Definition 11. An approximation-discretization to the Cauchy problem (7.1, 7.2) is numerically stable if there is some discrete norm in space and constants such that the corresponding approximation satisfies

for high enough resolutions, smooth initial data , and .Note:

- The previous definition applies both to the semi-discrete case (where space but not time is discretized) as well as the fully-discrete one. In the latter case, Eq. (7.3) is to be interpreted at fixed time. For example, if the timestep discretization is constant, then Eq. (7.3) needs to hold for fixed and arbitrarily large . In other words, the solution is allowed to grow with time, but not with the number of timesteps at fixed time when resolution is increased.
- The norm in general depends on the spatial approximation, and in Sections 8 and 9 we discuss some definitions for the finite difference and spectral cases.
- From Definition 11, one can see that an ill-posed problem cannot have a stable discretization, since otherwise one could take the continuum limit in (7.3) and reach the contradiction that the original system was well posed.
- As in the continuum, Eq. (7.3) implies uniqueness of the numerical solution .
- In Section 3 we discussed that if, in a well-posed homogeneous Cauchy problem, a forcing term is
added to Eq. (7.1),
then the new problem admits another estimate, related to the original one via Duhamel’s formula,
Eq. (3.23). A similar concept holds at the semi-discrete level, and the discrete estimates change
accordingly (in the fully-discrete case the integral in time is replaced by a discrete sum),
In other words, the addition of a lower-order term does not affect numerical stability, and without loss of generality one can restrict stability analyses to the homogeneous case.

- The difference between the exact solution and its numerical approximation satisfies an equation analogous to (7.5), where is related to the truncation error of the approximation. If the scheme is numerically stable, then in the linear and semi-discrete cases Eq. (7.6) implies If the approximation is consistent, the truncation error converges to zero as resolution is increased, and Equation (7.7) implies that so does the norm of the error . That is, stability implies convergence. The inverse is also true and this equivalence between convergence and stability is the celebrated Lax–Richtmyer theorem. The equivalence also holds in the fully-discrete case.
- In the quasi-linear case, one follows the principle of linearization, as described in Section 3.3. One linearizes the problem, and constructs a stable numerical scheme for the linearization. The expectation, then, is that the scheme also converges for the nonlinear scheme. For particular problems and discretizations this expectation can be rigorously proven (see, for example, [259]).

From here on denotes some discretization of space and time. This includes both finite difference and spectral collocation methods, which are the ones discussed in Sections 8 and 9, respectively. In addition, we use the shorthand notation

In order to gain some intuition into the general problem of numerical stability we start with some examples of simple, low-order approximations for a test problem. Consider uniform grids both in space and time

and the advection equation, on a periodic domain with , and smooth periodic initial data. Then the solution can be represented by a Fourier series: where and the stability of the following schemes can be analyzed in Fourier space.Example 33. The one-sided Euler scheme.

Eq. (7.10) is discretized with a one-sided FD approximation for the spatial derivative and evolved in time
with the forward Euler scheme,

Using Parseval’s identity, we find

and therefore, we see that the inequality (7.3) can only hold for all if For , this is the case if and only if the CFL factor satisfies and in this case the well-posedness estimate (7.3) holds with and . The upper bound in condition (7.19) for this example is known as the CFL limit, and (7.18) as the von Neumann condition. If , , while for the scheme is unconditionally unstable even though the underlying continuum problem is well posed.Next we consider a scheme very similar to the previous one, but which turns out to be unconditionally unstable for , regardless of the direction of propagation.

Example 34. A centered Euler scheme.

Consider first the semi-discrete approximation to Eq. (7.10),

The semi-discrete centered approximation (7.20) and the fully-discrete centered Euler scheme (7.21) constitute the simplest example of an approximation, which is not fully-discrete stable, even though its semi-discrete version is. This is related to the fact that the Euler time integration is not locally stable, as discussed in Section 7.3.2.

The previous two examples were one-step methods, where can be computed in terms of . The following is an example of a two-step method.

Example 35. Leap-frog.

A way to stabilize the centered Euler scheme is by approximating the time derivative by a centered
difference instead of a forward, one-sided operator:

In the above example the amplification matrix can be diagonalized through a transformation that is uniformly bounded:

with , , and The eigenvalues are of unit modulus, . In addition, the norms of and its inverse are Therefore, the condition number of can be bounded for all : provided that and it follows that the Leap-frog scheme is stable under the condition (7.30).The previous examples were explicit methods, where the solution (or ) can be explicitly computed from the one at the previous timestep, without inverting any matrices.

Example 36. Crank–Nicholson.

Approximating Eq. (7.10) by

Example 37. Iterated Crank–Nicholson.

Approximating the Crank–Nicholson scheme through an iterative scheme with a fixed number of iterations
is usually referred to as the Iterated Crank–Nicholson (ICN) method. For Eq. (7.10) it proceeds as
follows [414]:

- First iteration: an intermediate variable is calculated using a second-order-in-space centered difference (7.32) and an Euler, first-order forward-time approximation, Next, a second intermediate variable is computed through averaging, The full time step for this first iteration is
- Second iteration: it follows the same steps. Namely, the intermediate variables
- Further iterations proceed in the same way.

The resulting discretization is numerically stable for and iterations, and unconditionally unstable otherwise. In the limit the ICN scheme becomes the implicit, unconditionally-stable Crank–Nicholson scheme of the previous example. For any fixed number of iterations, though, the method is explicit and stability is contingent on the CFL condition . The method is unconditionally unstable for because the limit of the amplification factor approaching one in absolute value [cf. Eq. (7.34)] as increases is not monotonic. See [414] for details and [380] for a similar analysis for “theta” schemes.

Similar definitions to the one of Definition 11 are introduced for the IBVP. For simplicity we explicitly discuss the semi-discrete case. In analogy with the definition of a strongly–well-posed IBVP (Definition 9) one has

Definition 12. A semi-discrete approximation to the linearized version of the IBVP (5.1, 5.2, 5.3) is numerically stable if there are discrete norms at and at and constants and such that for high-enough resolution the corresponding approximation satisfies

for all . If the constant can be chosen strictly positive, the problem is called strongly stable.In addition, the semi-discrete version of Definitions 6 and 7 lead to the concepts of strong stability in the generalized sense and boundary stability, respectively, which we do not write down explicitly here. The definitions for the fully-discrete case are similar, with time integrals such as those in Eq. (7.39) replaced by discrete sums.

Living Rev. Relativity 15, (2012), 9
http://www.livingreviews.org/lrr-2012-9 |
This work is licensed under a Creative Commons License. E-mail us: |