### 7.4 Strict or time-stability

Consider a well-posed problem. For the sake of simplicity, we assume it is a linear initial-value one; a similar discussion holds for linear IBVPs. According to Definition 3, there are constants and such that smooth solutions satisfy

Definition 15. For a numerically-stable semi-discrete approximation, there are resolution independent constants such that for all initial data

The approximation is called strict or time stable if, for finite differences,

Similar definitions hold in the fully-discrete case and/or when the spatial approximation is not a finite difference one. Essentially, (7.115) attempts to capture the notion that the numerical solution should not have, at a fixed resolution, growth in time, which is not present at the continuum. However, the problem with the definition is that it is not useful if the estimate (7.113) is not sharp, since neither will be the estimate (7.114), and the numerical solution can still exhibit artificial growth.

Example 41. Consider the problem (drawn from [278])

where , for which the solution is the stationary one

Defining

it follows that
This energy estimate implies that the spatial norm of cannot grow faster than . However, in principle this need not be a sharp bound. In fact, Eq. (7.116) is, in disguise, the advection equation for which the general solution does not grow in time. Therefore, a numerical scheme for the problem (7.116, 7.117, 7.118) whose solutions are allowed to grow exponentially at the rate is strictly stable according to Definition 15 but the classification of the scheme as such is not of much use if the growth does take place. In order to illustrate this, we show the results for two schemes. In the first one, spurious growth takes place although the scheme is strictly stable, and in the second case one obtains a strictly-stable scheme with respect to a sharp energy estimate, which does not exhibit such growth.

If the system (7.116, 7.117, 7.118) is approximated by

where is a finite difference operator satisfying summation by parts (SBP) (Section 8.3) and the boundary condition is imposed through a projection method (Section 10), it can be shown – as discussed in Section 8 – that the following semi-discrete estimate holds, at least for analytic boundary conditions, with the discrete SBP scalar product defined by Eq. (8.21) in the next section,
Technically, the semi-discrete approximation (7.122) is strictly stable, since the continuum (7.121) and semi-discrete (7.123) estimates agree. However, this does not preclude spurious growth, as the bounds are not sharp. The left panel of Figure 3 shows results for the SBP operator (see the definition in Example 51 below), boundary conditions through an orthogonal projection, and third-order Runge–Kutta time integration. At any time, the errors do converge to zero with increasing resolution as expected since the scheme is numerically stable. However, at any fixed resolution there is spurious growth in time.

On the other hand, discretizing the system as

is, in general, not equivalent to (7.122) because difference operators do not satisfy the Leibnitz rule exactly. In fact, defining
it follows that
which, being an equality, is as sharp as an estimate can be.

The right panel of Figure 3 shows a comparison between discretizations (7.122) and (7.124), as well as (7.122) with the addition of numerical dissipation (see Section 8.5), in all cases at the same fixed resolution. Even though numerical dissipation does stabilize the spurious growth in time, the strictly-stable discretization (7.124) is considerably more accurate. Technically, according to Definition 15, the approximation (7.122) is also strictly stable, but it is more useful to reserve the term to the cases in which the estimate is sharp. The approximation (7.124), on the other hand, is (modulo the flux at boundaries, discussed in Section 10) energy preserving or conservative.

In order to construct conservative or time-stable semi-discrete schemes, one essentially needs to write the approximation by grouping terms in such a way that when deriving at the semi-discrete level what would be the conservation law at the continuum, the need of using the Leibnitz rule is avoided. In addition, the numerical imposition of boundary conditions also plays a role (see Section 10).

In many application areas, conservation or time-stability play an important role in the design of numerical schemes. That is not so much (at least so far) the case for numerical solutions of Einstein’s equations, because in general relativity there is no gauge-invariant local notion of conserved energy unlike many other nonlinear hyperbolic systems (most notably, in Newtonian or special relativistic Computational Fluid Dynamics); see, however, [400]. In addition, there are no generic sharp estimates for the growth of the solution that can guide the design of numerical schemes. However, in simpler settings such as fields propagating on some stationary fixed-background geometry, there is a notion of conserved local energy and accurate conservative schemes are possible. Interestingly, in several cases such as Klein–Gordon or Maxwell fields in stationary background spacetimes the resulting conservation of the semi-discrete approximations follows regardless of the constraints being satisfied (see, for example, [278]). A local conservation law in stationary spacetimes can also guide the construction of schemes to guarantee stability in the presence of coordinate singularities [105, 375, 225, 310], as discussed in Section 7.6.

In addition, there has been work done on variational, symplectic or mimetic integration techniques for Einstein’s equations, which aim at exactly or approximately preserving the discrete constraints, while solving the discrete evolution equations. See, for example, [304, 139, 201, 200, 76, 110, 359, 173, 358, 174, 360, 357].