Preservation of Strong Stability Associated with Analytic Semigroups

Stig Larsson^{∗}

1 Introduction

We begin by recalling the elementary example,

y′ + λy = 0, t > 0; y(0) = y_{0},  (1)

that was discussed in the Introduction. The solution y(t) = e^{−λt}y_{0} satisfies the stability estimate,

|y(t)| = e^{−λt}|y_{0}| ≤ |y_{0}|, t ≥ 0,  (2)

which is uniform with respect to λ > 0. Recall that we showed that this stability estimate is preserved by the implicit Euler method. However, the explicit Euler method is only stable provided that Δt < 2/λ, so its stability is not uniform in λ. One consequence is that it is less efficient to compute with the explicit Euler method over long time intervals.

Note that in this example, we also have a stability estimate for the derivative. Namely,

|y′(t)| = λe^{−λt}|y_{0}| ≤ Ct^{−1}|y_{0}|, t > 0,  (3)

where C = max_{s≥0}(se^{−s}) = e^{−1}. We call this property strong stability. Note that it is also uniform with respect to λ.

When we generalize (1) to a system of equations y′ + Ay = 0, where A is a square matrix and y is a vector, the situation becomes quite a bit more complicated. Now the eigenvalues of A play the role of λ in (1). A system of equations y′+Ay = 0 is dissipative, if the eigenvalues of A have positive real parts, which is guaranteed, for example, if A is symmetric and positive definite. In this case, we seek stability bounds on norms of y and its derivatives analogous to the bounds that hold for (1). The uniformity of such bounds with respect to the eigenvalues of A is crucial because we often expect the eigenvalues to vary considerably in size.

∗Department of Mathematics, Chalmers University of Technology, SE–412 96 Göteborg, Sweden. Email: stig@math.chalmers.se.

Published in “Collected Lectures on Preservation of Stability under Discretization” (D. Estep and S. Tavener, eds.), SIAM, 2002, pp. 11–24.

In this chapter we generalize (1) even further to parabolic partial differential equations written as

u′ + Au = 0, t > 0; u(0) = u_{0},  (4)

where A is a differential operator with respect to the spatial variable. The prototypical example is the heat equation with A being minus the Laplacian. We can write the solution formally as u(t) = exp(−tA)u_{0}, where exp(−tA) is the solution operator. We are interested in the strong stability mentioned above, which, in this case, means that strong norms (involving space and time derivatives) of the solution can be bounded in terms of a weak norm (involving no derivatives) of the initial value. It turns out that (4) has this property precisely when the solution operator exp(−tA) is a so called analytic semigroup.

The strong stability, alias analytic semigroup property, is preserved by certain numerical methods, which makes it possible to prove error bounds with optimal and low regularity requirements. In this chapter we review the theory of analytic semigroups and discuss how it can be used in the error analysis of numerical methods that preserve this structure.

2 The heat equation

We consider the homogeneous initial-boundary value problem for the heat equation,

\begin{aligned}u_{t}(x, t) − Δu(x, t) & = 0, & x ∈ Ω, t > 0, \\ u(x, t) & = 0, & x ∈ ∂Ω, t > 0, \\ u(x, 0) & = u_{0}(x), & x ∈ Ω, \end{aligned}  (5)

where Ω is a bounded domain in R^{d}, u = u(x, t), u_{t} = ∂u/∂t, Δu = ∑_{i=1}^{d} ∂^{2}u/∂x_{i}^{2}. For simplicity, we assume in this presentation that Ω is a convex polygon, so that we have access to the elliptic regularity theory and so that finite element meshes can be fitted exactly to the domain.

We use the Hilbert space L_{2}(Ω), with its standard norm and inner product

∥v∥ = ( ∫_{Ω} |v|^{2} dx)^{1/2}, (v, w) = ∫_{Ω} vw dx.  (6)

The norms in the Sobolev spaces H^{m}(Ω), m ≥ 0, are denoted by

∥v∥_{H^{m}} = ( ∑_{|α|≤m} ∥D^{α}v∥^{2})^{1/2}.  (7)

The space H_{0}^{1}(Ω) consists of the functions in H^{1}(Ω) that vanish on ∂Ω.

We define the unbounded operator A = −Δ on L_{2}(Ω) with domain of definition D(A) = H^{2}(Ω) ∩ H_{0}^{1}(Ω). Then A is a closed, densely defined, and self-adjoint positive definite operator on L_{2}(Ω) with compact inverse. Using this notation we may write (5) as an initial-value problem in L_{2}(Ω), namely,

u′ + Au = 0, t > 0; u(0) = u_{0},  (8)

where u : R^{+} → L_{2}(Ω), u′ = du/dt, and u_{0} ∈ L_{2}(Ω). The solution of (8) is given by u(t) = E(t)u_{0}, with the solution operator, E(t) = exp(−tA), defined by

E(t)v = ∑_{j=1}^{∞} e^{−tλ_{j}}(v, φ_{j})φ_{j}, v ∈ L_{2}(Ω),  (9)

where λ_{j} > 0 and φ_{j} denote the eigenvalues and an orthonormal basis of eigenvectors of A.

Using Parseval’s formula,

∥v∥^{2} = ∑_{j=1}^{∞}(v, φ_{j})^{2},  (10)

we easily obtain the stability estimate

∥E(t)v∥ ≤ ∥v∥, t ≥ 0,  (11)

In order to estimate the derivative of the solution, we use max_{s≥0}(s^{2}e^{−2s}) = e^{−2} < 1 to get

∥E′(t)v∥^{2} = ∑_{j=1}^{∞} λ_{j}^{2}e^{−2tλ_{j}}(v, φ_{j})^{2} = t^{−2} ∑_{j=1}^{∞}(tλ_{j})^{2}e^{−2tλ_{j}}(v, φ_{j})^{2} ≤ t^{−2}∥v∥^{2}. (12)

Hence, we obtain the strong stability estimate:

∥E′(t)v∥ = ∥AE(t)v∥ ≤ t^{−1}∥v∥, t > 0.  (13)

The term strong stability refers to the fact that the norms on the left sides are stronger than the one on the right side.

Since we assume that the domain Ω is a convex polygon, we may take advantage of the regularity of elliptic equations,

∥v∥_{H^{2}} ≤ C∥Av∥, ∀v ∈ D(A) = H^{2}(Ω) ∩ H_{0}^{1}(Ω),  (14)

and (13) to get a related estimate of the second order spatial derivatives:

∥E(t)v∥_{H^{2}} ≤ Ct^{−1}∥v∥, t > 0.  (15)

The strong stability estimate (13) reflects the fact the solution operator, E(t) = exp(−tA), is an analytic semigroup on L_{2}(Ω). We proceed to discuss this concept in the following section.

3 Analytic semigroups

Let X be a Banach space with norm ∥·∥. A family {E(t)}_{t≥0} of bounded linear operators on X is called a semigroup of bounded linear operators, if

• E(0) = I (identity operator),

• E(t + s) = E(t)E(s) ∀ s, t ≥ 0 (semigroup property).

The semigroup is called strongly continuous, if

lim_{t→0^{+}} E(t)x = x ∀x ∈ X.  (16)

The infinitesimal generator of the semigroup is the linear operator G defined by

Gx = lim_{t→0^{+}} \frac{E(t)x − x}{t},  (17)

its domain of definition D(G) being the space of all x ∈ X for which the limit exists. One can prove that G is a closed and densely defined operator on X and that the function u(t) = E(t)u_{0} solves the initial-value problem

u′ = Gu, t > 0; u(0) = u_{0},  (18)

provided that u_{0} ∈ D(G). We therefore often write E(t) = exp(tG). We refer to [12, Chapter 1] for the proofs of these statements. In Section 2 the semigroup E(t) = exp(−tA) on L_{2}(Ω) is generated by −A = Δ.

On the other hand we might want to start with a linear operator A and seek conditions under which G = −A generates a semigroup E(t) = exp(−tA) such that u(t) = E(t)u_{0} solves

u′ + Au = 0, t > 0; u(0) = u_{0}.  (19)

Theorem 1 below provides such a condition.

The semigroup E(t) is called analytic or holomorphic, if it can be extended to a complex analytic function E(z) for z in a sector in the complex plane containing the positive t-axis. The following theorem, which we quote from [12, Chapter 2.5], gives two characterizations of an analytic semigroup; one in terms of the resolvent of its generator and the other one in terms of its derivative.

Recall that the resolvent set ρ(A) of the linear operator A is the set of complex numbers λ such that λI − A is invertible and the resolvent (λI − A)^{−1} is a bounded linear operator on X.

Theorem 1. Let E(t) be a strongly continuous semigroup with generator G = −A. Assume that E(t) is uniformly bounded, i.e., for some C,

∥E(t)∥ ≤ C, t ≥ 0,  (20)

and 0 ∈ ρ(A). Then the following are equivalent:

(a) The semigroup E(t) can be extended to a complex analytic function E(z) in a sector {z ∈ C : | arg(z)| ≤ δ}.

(b) There is an angle φ ∈ (0, π/2) and a positive number M such that ρ(A) ⊃ Σ_{φ} = {λ ∈ C : φ ≤ | arg(λ)| ≤ π} ∪ {0} and

∥(λI − A)^{−1}∥ ≤ \frac{M}{|λ|}, ∀λ ∈ Σ_{φ}, λ ≠ 0.  (21)

(c) E(t) is differentiable for t > 0 and there is a constant C such that

∥AE(t)∥ = ∥E′(t)∥ ≤ Ct^{−1}, t > 0.  (22)

Here ∥·∥ is the usual norm of bounded linear operators on X.

The differentiability in property (c) implies that the function u(t) = E(t)u_{0} solves (18) for all u_{0} ∈ X, not just for u_{0} ∈ D(G) = D(A).

The proof of the implication (b) ⇒ (c) is based on the formula

E(t) = \frac{1}{2πi} ∫_{Γ_{φ}} e^{−tλ}(λI − A)^{−1} dλ,  (23)

where the contour Γ_{φ} = {λ ∈ C : | arg(λ)| = φ} is the boundary of the sector in (21) oriented so that the imaginary part decreases along Γ_{φ}. The resolvent estimate (21) may then be used to show that the integral in (23) is convergent and defines a uniformly bounded analytic semigroup.

Instead of starting with an analytic semigroup with generator G = −A defined by (17), we may therefore start with a densely defined operator A with property (b) and define E(t) by (23). This will then be a uniformly bounded analytic semigroup with generator G = −A.

In Section 2 we showed that the semigroup discussed there has the property (c), see (13). This shows that it is analytic on X = L_{2}(Ω) and that the resolvent estimate in (b) holds for some angle φ ∈ (0, π/2). On the other hand, it is not hard to prove property (b) with arbitrary φ ∈ (0, π/2) for a self-adjoint, positive definite, linear operator A.

More generally, one can prove that a strongly elliptic partial differential operator of second order in a smooth bounded domain Ω together with the homogeneous Dirichlet boundary condition generates an analytic semigroup in L_{p}(Ω) for 1 < p < ∞. The approach is to prove the resolvent estimate (b) with respect to the L_{p}-norm, see [12, Chapter 7.3].

For such differential operators it is also possible to prove the resolvent estimate (b), and hence the strong stability (c), with respect to the L_{∞}-norm, see [16]. For example, the semigroup (9) satisfies

∥E(t)v∥_{L_{∞}} ≤ ∥v∥_{L_{∞}}, t ≥ 0 (maximum principle),
∥E′(t)v∥_{L_{∞}} ≤ Ct^{−1}∥v∥_{L_{∞}}, t > 0. (24)

However, the resulting semigroup is not analytic (not even strongly continuous) in L_{∞}(Ω), because the domain of definition of the generator is not dense in L_{∞}(Ω). Instead one has to work in the Banach space X = {u ∈ C(Ω̄) : u = 0 on ∂Ω}.

Unlike the case X = L_{p}(Ω), this space includes a boundary condition. This is often not desirable, because when the semigroup is used to solve the non-homogeneous equation u′ + Au = f, the boundary condition is imposed on f. There is also a theory of semigroups generated by operators with non-dense domain, [13].

Since our main interest here is the preservation of the strong stability under discretization, we shall not dwell on this, but proceed to the next section, where we study spatial discretization by the finite element method.

4 Spatial discretization

We now return to the heat equation as described in Section 2 and consider the non-homogeneous version of (8):

u′ + Au = f, t > 0; u(0) = u_{0}.  (25)

If the function f : R^{+} → L_{2}(Ω) is continuous, then the solution is given by

u(t) = E(t)u_{0} + ∫_{0}^{t} E(t − s)f(s) ds,  (26)

where E(t) = exp(−tA), defined in (9), is the analytic semigroup generated by −A. The weak formulation of (25) is: find u(t) ∈ H_{0}^{1}(Ω) such that

\begin{aligned}(u′, v) + a(u, v) & = (f, v), ∀v ∈ H_{0}^{1}(Ω), t > 0, \\ u(0) & = u_{0}, \end{aligned}  (27)

where a(u, v) = (∇u, ∇v) = (−Δu, v) = (Au, v) is the bilinear form associated with the operator A.

Let {V_{h}}_{0<h<1} be a family of finite dimensional subspaces of H_{0}^{1}, where each V_{h} consists of continuous piecewise polynomials of degree ≤ 1 with respect to a triangulation T_{h} of Ω with maximal mesh size h. We refer to, e.g., [4] and [5] for results about the finite element method for elliptic problems. The general reference for parabolic finite element problems is [17].

The approximate solution u_{h}(t) ∈ V_{h} of (25) is defined by

\begin{aligned}(u′_{h}, χ) + a(u_{h}, χ) & = (f, χ), ∀χ ∈ V_{h}, t > 0, \\ (u_{h}(0), χ) & = (u_{0}, χ), ∀χ ∈ V_{h}. \end{aligned}  (28)

Introducing the linear operator A_{h} : V_{h} → V_{h} and the orthogonal projection P_{h} : L_{2}(Ω) → V_{h}, defined by

(A_{h}ψ, χ) = a(ψ, χ), (P_{h}g, χ) = (g, χ) ∀ψ, χ ∈ V_{h}, g ∈ L_{2}(Ω),  (29)

we may write (28) as

u′_{h} + A_{h}u_{h} = P_{h}f, t > 0; u_{h}(0) = P_{h}u_{0}.  (30)

The operator A_{h} is self-adjoint positive definite (uniformly with respect to h) and the solution of (30) is given by

u_{h}(t) = E_{h}(t)P_{h}u_{0} + ∫_{0}^{t} E_{h}(t − s)P_{h}f(s) ds,  (31)

with the solution operator, E_{h}(t) = exp(−tA_{h}), defined by

E_{h}(t)v_{h} = ∑_{j=1}^{N_{h}} e^{−tλ_{h,j}}(v_{h}, φ_{h,j})φ_{h,j}, v_{h} ∈ V_{h},  (32)

where λ_{h,j} > 0 and φ_{h,j} denote the eigenvalues and an orthonormal basis of eigenvectors of A_{h}. Due to the similarity of this formula with (9), we immediately obtain stability estimates analogous to (11), (13):

∥E_{h}(t)v_{h}∥ ≤ ∥v_{h}∥, t ≥ 0,  (33)
∥E′_{h}(t)v_{h}∥ = ∥A_{h}E_{h}(t)v_{h}∥ ≤ t^{−1}∥v_{h}∥, t > 0.  (34)

Note in particular that these hold uniformly with respect to h. We may therefore say that the strong stability associated with the analytic semigroup E(t) is preserved under the spatial discretization.

We next demonstrate how the strong stability estimate (34) may be used in the error analysis. We need the Ritz projection R_{h} : H_{0}^{1}(Ω) → V_{h} defined by

a(R_{h}v, χ) = a(v, χ), ∀χ ∈ V_{h}.  (35)

A simple calculation shows that

R_{h} = A_{h}^{−1}P_{h}A.  (36)

In fact, R_{h}v is the finite element solution of an elliptic equation whose exact solution is v. Under the usual regularity assumptions on the triangulation and using the elliptic regularity estimate (14), we obtain the following error bound from the theory of finite elements for elliptic problems (see [4], [5])

∥R_{h}v − v∥ ≤ Ch^{2}∥v∥_{H^{2}}, ∀v ∈ H_{0}^{1}(Ω) ∩ H^{2}(Ω).  (37)

Let u and u_{h} denote the solutions of (25) and (30), respectively. We divide the error into two parts:

u_{h}(t) − u(t) = (u_{h}(t) − P_{h}u(t)) + (P_{h}u(t) − u(t)).  (38)

For the last term it follows from the optimality of the orthogonal projection and (37) that

∥P_{h}u(t) − u(t)∥ ≤ ∥R_{h}u(t) − u(t)∥ ≤ Ch^{2}∥u(t)∥_{H^{2}},  (39)

and it remains to estimate the first part θ(t) = u_{h}(t) − P_{h}u(t). A simple calculation using (30), (25), and (36) yields

\begin{matrix}θ′ + A_{h}θ = A_{h}P_{h}(A_{h}^{−1}P_{h}A − I)u = A_{h}P_{h}(R_{h} − I)u = A_{h}P_{h}ρ, \\ θ(0) = 0, \end{matrix}  (40)

where we introduced the notation ρ = R_{h}u − u. Hence, by the solution formula (32), we get

θ(t) = ∫_{0}^{t} A_{h}E_{h}(t − s)P_{h}ρ(s) ds,  (41)

and it should now be clear that (34) is useful.

We will derive several error estimates. Applying (34) directly, we note that ∥A_{h}E_{h}(t − s)∥ ≤ (t − s)^{−1} has a non-integrable singularity. In order to handle this we assume that our family of finite element meshes is such that we have an inverse estimate

∥A_{h}v_{h}∥ ≤ Ch^{−β}∥v_{h}∥, ∀v_{h} ∈ V_{h},  (42)

for some positive number β. For example, if the mesh family is quasi-uniform then this holds with β = 2. Hence, using also (33) and (37), we obtain

∥θ(t)∥ ≤ ( ∫_{0}^{t−h^{β}}(t − s)^{−1} ds + Ch^{−β} ∫_{t−h^{β}}^{t} ds) sup_{s≤t} ∥ρ(s)∥

≤ Ch^{2}(1 + log(t/h^{β})) sup_{s≤t} ∥u(s)∥_{H^{2}}, t > h^{β}.

(43)

Another way to avoid the singularity is to note that A_{h}E_{h}(t−s) = D_{s}E_{h}(t−s) and integrate by parts:

\begin{aligned}θ(t) & = ∫_{0}^{t/2} A_{h}E_{h}(t − s)P_{h}ρ(s) ds \\ & + P_{h}ρ(t) − E_{h}(t/2)P_{h}ρ(t/2) − ∫_{t/2}^{t} E_{h}(t − s)P_{h}ρ′(s) ds. \end{aligned}  (44)

Hence, by (33), (34), and (37),

∥θ(t)∥ ≤ ( ∫_{0}^{t/2}(t − s)^{−1} ds + 2) sup_{s≤t} ∥ρ(s)∥ + ∫_{t/2}^{t} s^{−1} ds sup_{s≤t} (s∥ρ′(s)∥)

≤ Ch^{2} sup_{s≤t} (∥u(s)∥_{H^{2}} + s∥u′(s)∥_{H^{2}}).

(45)

A third variant is obtained by integrating the first term on the right side of (44) by parts. With ρ̃(t) = ∫_{0}^{t} ρ(s) ds we get

\begin{aligned}θ(t) & = A_{h}E_{h}(t/2)P_{h}ρ̃(t/2) − ∫_{0}^{t/2} A_{h}^{2}E_{h}(t − s)P_{h}ρ̃(s) ds \\ & + P_{h}ρ(t) − E_{h}(t/2)P_{h}ρ(t/2) − ∫_{t/2}^{t} E_{h}(t − s)P_{h}ρ′(s) ds, \end{aligned}  (46)

so that

∥θ(t)∥ ≤ (2t^{−1} + ∫_{0}^{t/2}(t − s)^{−2} ds) sup_{s≤t} ∥ρ̃(s)∥

+ 2t^{−1} sup_{s≤t} (s∥ρ(s)∥) + ∫_{t/2}^{t} s^{−2} ds sup_{s≤t} (s^{2}∥ρ′(s)∥)

≤ Ch^{2}t^{−1} sup_{s≤t} (∥ũ(s)∥_{H^{2}} + s∥u(s)∥_{H^{2}} + s^{2}∥u′(s)∥_{H^{2}}).

(47)

Adding the contribution from (39) we conclude the following.

Theorem 2. Let u and u_{h} denote the solutions of (25) and (30), respectively. Then, with ũ(t) = ∫_{0}^{t} u(s) ds,

∥u_{h}(t) − u(t)∥ ≤ Ch^{2}(1 + log(t/h^{β})) sup_{s≤t} ∥u(s)∥_{H^{2}}, t > h^{β},  (48)
∥u_{h}(t) − u(t)∥ ≤ Ch^{2} sup_{s≤t} (∥u(s)∥_{H^{2}} + s∥u′(s)∥_{H^{2}}), t ≥ 0,  (49)
∥u_{h}(t) − u(t)∥ ≤ Ch^{2}t^{−1} sup_{s≤t} (∥ũ(s)∥_{H^{2}} + s∥u(s)∥_{H^{2}} + s^{2}∥u′(s)∥_{H^{2}}), t > 0. (50)

The first error estimate, (48), requires the least regularity of the three. Except for the slowly growing logarithmic factor, it is also uniform in time, while the other bounds may accumulate slightly because of the factor s that appears in front of some terms.

In order to compare the last two error estimates we apply them to the solution of the homogeneous equation (8), i.e., u(t) = E(t)u_{0}. Then, according to (15),

∥u(t)∥_{H^{2}} ≤ Ct^{−1}∥u_{0}∥,  (51)

and similarly

∥u′(t)∥_{H^{2}} ≤ Ct^{−2}∥u_{0}∥,  (52)

and

∥ũ(t)∥_{H^{2}} ≤ C∥Aũ(t)∥ = C∥∫_{0}^{t} Au(s) ds∥ = C∥∫_{0}^{t} u′(s) ds∥ ≤ C∥u_{0}∥. (53)

Therefore (50) becomes

∥E_{h}(t)P_{h}u_{0} − E(t)u_{0}∥ ≤ Ch^{2}t^{−1}∥u_{0}∥, t > 0,  (54)

which means that the error converges with optimal order even if the initial value is only in L_{2}(Ω). On the other hand, if u_{0} ∈ D(A) = H^{2}(Ω) ∩ H_{0}^{1}(Ω), then ∥AE(t)u_{0}∥ ≤ ∥Au_{0}∥, so that

∥u(t)∥_{H^{2}} + t∥u′(t)∥_{H^{2}} ≤ C∥u_{0}∥_{H^{2}},  (55)

and (49) becomes

∥E_{h}(t)P_{h}u_{0} − E(t)u_{0}∥ ≤ Ch^{2}∥u_{0}∥_{H^{2}}, t ≥ 0.  (56)

Thus, the error bound is uniform as t → 0 if the initial value is smooth.

In Theorem 2 we used the strong stability (34) in order to prove three error estimates with low and optimal regularity requirement. Using only the stability estimate (33) we obtain instead

∥u_{h}(t) − u(t)∥ ≤ Ch^{2}( sup_{s≤t} ∥u(s)∥_{H^{2}} + ∫_{0}^{t} ∥u′(s)∥_{H^{2}} ds), t ≥ 0,  (57)

which is slightly weaker than (49), e.g., it does not imply (56). It is proved by integration by parts as in (44) but with t/2 replaced by t.

The previous discussion is based on the use of the spectral theorem for the self-adjoint operator A = −Δ in the Hilbert space X = L_{2}(Ω). When the operator A is not self-adjoint or the space X is a Banach space, we base the argument on the contour integral (23) and the resolvent estimate (21). We demonstrate this by analyzing the finite element approximation of (5) in L_{∞}(Ω). We assume from now on that the mesh family is quasi-uniform.

In [14] it was shown, in the case of two spatial variables, that the semigroup (32) satisfies

∥E_{h}(t)v_{h}∥_{L_{∞}} + t∥E′_{h}(t)v_{h}∥_{L_{∞}} ≤ C log(1/h)∥v_{h}∥_{L_{∞}}, t > 0.  (58)

The same bound, but with an additional factor log(1/h), was proved in [18] in the case of one spatial variable. This means that the strong stability in the maximum norm (24) is almost preserved by the finite element discretization. Combining this with the maximum norm analog of (37), see [4],

∥R_{h}v − v∥_{L_{∞}} ≤ Ch^{2} log(1/h)∥v∥_{W_{∞}^{2}}, ∀v ∈ H_{0}^{1}(Ω) ∩ W_{∞}^{2}(Ω),  (59)

we may prove error estimates in the maximum norm analogous to Theorem 2 and (54), (56), (57), with (non-important) logarithmic factors.

However, recent work in [15], [19] has removed the logarithm from (58) and we now know that

∥E_{h}(t)v_{h}∥_{L_{∞}} + t∥E′_{h}(t)v_{h}∥_{L_{∞}} ≤ C∥v_{h}∥_{L_{∞}}, t > 0,  (60)

holds in any number of spatial variables. The logarithm in (59) cannot be removed, [10], so the resulting error estimates still contain logarithmic factors.

In view of Theorem 1 we know that (60) implies that the resolvent estimate

∥(λI − A_{h})^{−1}∥_{L_{∞}→L_{∞}} ≤ \frac{M}{|λ|}, ∀λ ∈ Σ_{φ}, λ ≠ 0,  (61)

holds for some angle φ ∈ (0, π/2). However, the angle in (21) increases towards π/2, meaning that the estimate deteriorates, as the constant C in (22) increases.

We do not have precise information about the constant in (60), and therefore no information about the angle in (61). Note, by the way, that (58) would allow φ → π/2 as h → 0, so here is an example where the removal the logarithm makes a difference. Anyway, this is not satisfactory, since we know that the corresponding resolvent estimate for the differential operator A = −Δ holds for any angle. In the following section I will explain why it may be important to know the angle.

Therefore it seems advantageous to begin by proving the resolvent estimate and then deduce the strong stability. This approach was taken in, e.g., [20], [7], [1], [2]. We now know that the resolvent estimate (61) holds with arbitrary angle. This was first proved in [7] for the case of one spatial variable, and recently for the general case in [2].

5 Temporal discretization

In this subsection we show that the above program can be carried through also for certain completely discrete schemes. We begin by looking at the backward Euler method. We replace the time derivative in (30) by a backward difference quotient ∂_{t}U_{n} = (U_{n} − U_{n−1})/k, where k is a time step and U_{n} is the approximation of u_{n} = u(t_{n}) and t_{n} = nk. The discrete solution U_{n} ∈ V_{h} thus satisfies

∂_{t}U_{n} + A_{h}U_{n} = P_{h}f_{n}, t_{n} > 0; U_{0} = P_{h}u_{0}.  (62)

The solution of (62) is given by

U_{n} = E_{kh}^{n}P_{h}u_{0} + k ∑_{j=1}^{n} E_{kh}^{n−j−1}P_{h}f_{j},  (63)

where the solution operator, E_{kh}^{n} = (I + kA_{h})^{−n}, may be written

E_{kh}^{n}v_{h} = r(kA_{h})^{n} = ∑_{j=1}^{N_{h}} r(kλ_{h,j})^{n}(v_{h}, φ_{h,j})φ_{h,j}, v_{h} ∈ V_{h},  (64)

where r(s) = (1 + s)^{−1} and where λ_{h,j} > 0 and φ_{h,j} denote the eigenvalues and an orthonormal basis of eigenvectors of A_{h}. Since r(s) ≤ 1 and |sr(s)^{n}| ≤ Cn^{−1} for s ≥ 0, see [17, Lemma 7.3], we easily obtain stability estimates analogous to (11), (13):

∥E_{kh}^{n}v_{h}∥ ≤ ∥v_{h}∥, t_{n} ≥ 0,  (65)
∥∂_{t}E_{kh}^{n}v_{h}∥ = ∥A_{h}E_{kh}^{n}v_{h}∥ ≤ Ct_{n}^{−1}∥v_{h}∥, t_{n} > 0.  (66)

Note in particular that these hold uniformly with respect to h and k. We may therefore say that the strong stability associated with the analytic semigroup E(t) is preserved under this spatial and temporal discretization. More generally, this holds for rational functions that are strongly A-stable (see [17, Lemma 7.3]):

|r(s)| ≤ 1, ∀s ≥ 0; r(∞) = 0.  (67)

Imitating the proof of Theorem 2 we obtain, instead of (41),

θ_{n} = k ∑_{j=1}^{n} A_{h}E_{kh}^{n−j−1}(P_{h}ρ_{j} + A_{h}^{−1}P_{h}ω_{j}),  (68)

where ω_{n} = ∂_{t}U_{n} − u′(t_{n}). Using summation by parts as in (44) we may prove, e.g.,

\begin{aligned}∥U_{n} − u(t_{n})∥ & ≤ Ch^{2} sup_{s≤t} (∥u(s)∥_{H^{2}} + s∥u′(s)∥_{H^{2}}) \\ & + Ck sup_{s≤t} (∥A^{−1}u′′(s)∥ + s∥u′′(s)∥), t_{n} ≥ 0, \end{aligned}  (69)

which is analogous to (49).

Leaving the Hilbert space setting, we now consider a rational approximation of the semigroup E(t) of the form

E_{kh}^{n}v_{h} = r(kA_{h})^{n},  (70)

where r is an A(θ)-stable rational function, i.e., for some angle θ ∈ [0, π/2],

|r(λ)| ≤ 1, for | arg(λ)| ≤ θ.  (71)

The operator is defined by the Dunford-Taylor integral,

r(kA_{h}) = r(∞)I + \frac{1}{2πi} ∫_{Γ_{φ}} r(kλ)(λI − A_{h})^{−1} dλ,  (72)

where the contour Γ_{φ} = {λ ∈ C : | arg(λ)| = φ} is the boundary of the sector in the maximum-norm resolvent estimate (61) oriented so that the imaginary part decreases along Γ_{φ}. Recall from the previous section that the angle φ can be chosen arbitrarily and, in particular, does not depend on h. If θ ≥ φ, then we can use (71) together with (61) to prove the stability estimate, see, e.g., [17, Theorem 8.2], [11], [6],

∥E_{kh}^{n}v_{h}∥_{L_{∞}} ≤ C∥v_{h}∥_{L_{∞}}, t_{n} ≥ 0.  (73)

Thus, we need the analyticity of the discrete semigroup E_{h}(t), in the form of the resolvent estimate, already in order to prove the standard stability. In fact, (73) is not true in general if the semigroup is only strongly continuous, [3], [6].

Looking for strong stability, we need to assume that the rational function r is strongly A(θ)-stable, i.e., we assume in addition that r(∞) = 0. Then one can prove

∥A_{h}E_{kh}^{n}v_{h}∥_{L_{∞}} ≤ Ct_{n}^{−1}∥v_{h}∥_{L_{∞}}, t_{n} > 0,  (74)

see [9]. This opens the possibility of proving error estimates analogous to those of Theorem 2.

For a similar result for a particular class of time discretizations but with variable timesteps, see [8].

Bibliography

[1] N.Yu. Bakaev, Maximum norm resolvent estimates for elliptic finite element operators, BIT, 41, 2001, pp. 215–239.

[2] N.Yu. Bakaev, V. Thomée, and L.B. Wahlbin, Maximum-norm estimates for resolvents of elliptic finite element operators, Tech. Report 2001:58, Department of Mathematics, Chalmers University of Technology and Göteborg University, Göteborg, Sweden, 2001.

[3] P. Brenner and V. Thomée, On rational approximation of semigroups, SIAM J. Numer. Anal., 16, 1979, pp. 683–694.

[4] S.C. Brenner and L.R. Scott, The Mathematical Theory of Finite Element Methods, Texts in Applied Mathematics, vol. 15, Springer-Verlag, New York, 1994.

[5] P.G. Ciarlet, The Finite Element Method for Elliptic Problems, North-Holland, Amsterdam, 1978.

[6] M. Crouzeix, S. Larsson, S. Piskarev, and V. Thomée, The stability of rational approximations of analytic semigroups, BIT, 33, 1993, pp. 74–84.

[7] M. Crouzeix, S. Larsson, and V. Thomée, Resolvent estimates for elliptic finite element operators in one dimension, Math. Comp., 63, 1994, pp. 121–140.

[8] K. Eriksson, C. Johnson, and S. Larsson, Adaptive finite element methods for parabolic problems. VI. Analytic semigroups, SIAM J. Numer. Anal., 35, 1998, pp. 1315–1325 (electronic).

[9] A. Hansbo, Strong stability and non-smooth data error estimates L_{s} → L_{r}, s ≤ 2 ≤ r, Tech. Report 1999:53, Department of Mathematics, Chalmers University of Technology and Göteborg University, Göteborg, Sweden, 1999, to appear in BIT.

[10] R. Haverkamp, Eine Aussage zur L_{∞}-Stabilität und zur genauen Konvergenzordnung der H_{0}^{1}-Projektionen, Numer. Math., 44, 1984, pp. 393–405.

[11] C. Palencia, A stability result for sectorial operators in Banach spaces, SIAM J. Numer. Anal., 30, 1993, pp. 1373–1384.

[12] A. Pazy, Semigroups of Linear Operators and Applications to Partial Differential Equations, Applied Mathematical Sciences, vol. 44, Springer-Verlag, New York, 1983.

[13] G. Da Prato and E. Sinestrari, Differential operators with non dense domain, Ann. Scuola Norm. Sup. Pisa Cl. Sci., 14, 1987, pp. 285–344.

[14] A.H. Schatz, V. Thomée, and L.B. Wahlbin, Maximum norm stability and error estimates in parabolic finite element equations, Comm. Pure Appl. Math., 33, 1980, pp. 265–304.

[15] A.H. Schatz, V. Thomée, and L.B. Wahlbin, Stability, analyticity, and almost best approximation in maximum norm for parabolic finite element equations, Comm. Pure Appl. Math., 51, 1998, pp. 1349–1385.

[16] H.B. Stewart, Generation of analytic semigroups by strongly elliptic operators, Trans. Amer. Math. Soc., 199, 1974, pp. 141–162.

[17] V. Thomée, Galerkin Finite Element Methods for Parabolic Problems, Springer Series in Computational Mathematics, vol. 25, Springer-Verlag, Berlin, 1997.

[18] sc V. Thomée and L.B. Wahlbin, Maximum-norm stability and error estimates in Galerkin methods for parabolic equations in one space variable, Numer. Math., 41, 1983, pp. 345–371.

[19] V. Thomée and L.B. Wahlbin, Stability and analyticity in maximum-norm for simplicial Lagrange finite element semidiscretizations of parabolic equations with Dirichlet boundary conditions, Numer. Math., 87, 2000, pp. 373–389.

[20] L.B. Wahlbin, A quasioptimal estimate in piecewise polynomial Galerkin approximation of parabolic problems, in Numerical Methods, Proceedings, Dundee 1981, G. A. Watson, ed., Lecture Notes in Mathematics, vol. 912, Springer-Verlag, New York, 1981, pp. 230–245.