SIAM J. NUMER. ANAL. Vol. 31, No. 4, August 1994

\begin{matrix}c \\ ◯ \end{matrix}1994 Society for Industrial and Applied Mathematics

1000–1018

THE BEHAVIOR OF FINITE ELEMENT SOLUTIONS OF SEMILINEAR PARABOLIC PROBLEMS NEAR STATIONARY POINTS

S. LARSSON^{†} and J.-M. SANZ-SERNA^{‡}

Abstract. We study the qualitative behavior of spatially semidiscrete finite element solutions of a semilinear parabolic problem near an unstable hyperbolic equilibrium ū. We show that any continuous trajectory is approximated by an appropriate discrete trajectory, and vice versa, as long as they remain in a sufficiently small neighborhood of ū. Error bounds of optimal order in the L_{2} and H^{1} norms hold uniformly over arbitrarily long time intervals. In particular, the local stable and unstable manifolds of the discrete problem converge to their continuous counterparts. Therefore, the discretized dynamical system has the same qualitative behavior near ū as the continuous system.

Key words. Semilinear parabolic problem, unstable hyperbolic stationary point, finite element method, error estimate, stable, unstable manifold, shadowing

AMS subject classifications. 65M15, 65M60

1. Introduction. Classical error estimates for approximate solutions of nonlinear evolution problems on a time interval τ < t < T involve an error constant that grows exponentially with the length T − τ of the interval. In general this is the best that can be expected, because trajectories of initial value problems may diverge from each other at an exponential rate. On the other hand, exponentially growing bounds become meaningless even for moderate values of T − τ in situations where trajectories actually contract as time increases.

A typical instance of “contracting” trajectories occurs in the neighborhood of a stable equilibrium, a case that has been considered for ordinary differential equations by Stetter [12, Chapters 3.5 and 4.6] and for nonlinear parabolic partial differential equations by Heywood and Rannacher [8], Larsson [9] and Sanz-Serna and Stuart [11]. It turns out that, under suitable technical assumptions, the exponential growth of the error predicted by the classical bounds does not materialize: numerical methods provide approximations that are accurate uniformly in t near a stable equilibrium.

Numerical methods cannot be expected to do very well in a neighborhood V of an unstable equilibrium ū due to the divergence of trajectories. Figure 1 depicts a situation, which is typical for an unstable hyperbolic equilibrium (cf. Section 2 below): there is a local stable manifold M_{S} composed of solutions that remain in V for t > 0 and approach the equilibrium as t → ∞. The local unstable manifold M_{U} is composed of solutions that are in V for t < 0 and approach the equilibrium as t → −∞. Trajectories not contained in the stable or unstable manifolds remain in V during a finite time interval τ < t < T; the length of this interval depends on the individual trajectory and may be arbitrarily large. By considering initial values on opposite sides of M_{S}, it is easily seen that the distance between trajectories may grow exponentially with t, even if they are close initially. Therefore we are in a situation,

† Department of Mathematics, Chalmers University of Technology and University of Göteborg, S–412 96 Göteborg, Sweden (E-mail: stig@math.chalmers.se).

‡ Departamento de Matemática Aplicada y Computación, Facultad de Ciencias, Universidad de Valladolid, Valladolid, Spain (E-mail: sanzserna@cpd.uva.es). Partially supported by “Dirección General de Investigación Científica y Técnica” under project PB89–0351.

Fig. 1. Trajectories in the neighborhood of an unstable hyperbolic equilibrium.

where the exponential error growth predicted by classical bounds is not pessimistic and numerical methods cannot be expected to approximate individual trajectories well over long time intervals.

In a recent paper [3], restricted to ordinary differential equations, Beyn showed that numerical methods can nevertheless be useful for long-time integration near an unstable hyperbolic equilibrium: any continuous trajectory is approximated by an appropriate discrete trajectory, and vice versa, to the correct order of approximation and uniformly in t as long as they remain in a sufficiently small neighborhood V of ū. The conclusion is that, when a numerical method is used to simulate the phase portrait in V through the generation of a number of numerical trajectories, we can be sure that the result is accurate: there are trajectories of the original problem that are close, uniformly in t, to the computed trajectories. The difference with the case of a stable equilibrium is that there the discrete trajectory approximates the continuous trajectory starting from the same initial value, while here a given discrete trajectory approximates a continuous trajectory with a different (and a priori unknown) initial value.

The purpose of the present work is to extend Beyn’s result to the case of partial differential equations. We consider a model situation, where a semilinear parabolic problem is discretized in space by piecewise linear finite elements. Several generalizations are possible but we have decided not to deal with them in order to gain in clarity. Our main result, Theorem 1, is analogous to that of Beyn: there is an H^{1}-neighborhood V of the unstable hyperbolic stationary point ū, such that for any exact solution u(t) contained in V for τ < t < T there is a numerical solution u^{h}(t), which approximates u(t) accurately for τ < t < T. We emphasize that the initial value u^{h}(τ) is a priori unknown: a numerical solution starting from some a priori prescribed approximation of u(τ) (as in the standard error analysis) may deviate from u(t) at an exponential rate as noted before. Conversely, for any numerical solution u^{h}(t) contained in V for τ < t < T there exists an exact solution u(t), which is close to it. The distance between u(t) and u^{h}(t) is measured in the L_{2} and H^{1} norms, and we prove error bounds of optimal order with error constants that are independent of T − τ. In a similar manner, we show in Theorem 2 that the discrete problem has local stable and unstable manifolds which converge at an optimal rate to their continuous counterparts. The conclusion is that the discretized dynamical system has the same qualitative behavior near ū as the continuous system.

Theorem 1 is related to the concept of shadowing in dynamical systems theory, but the classical shadowing lemma is not directly applicable in the present situation, see, for example, [6] and [10].

Our work is in a sense complementary to that of Alouges and Debussche [1], which deals with discretization in time but not in space. However it is not trivial to combine our result with theirs to obtain a result for a completely discrete scheme. We plan to study space and time discretization simultaneously in a future work.

Section 2 is devoted to the presentation and analysis of the partial differential equation to be solved and Section 3 deals with the discretization in space. Our main results, Theorems 1 and 2, are stated at the end of Section 3.

2. The continuous problem. Throughout this work Ω is a bounded convex polygonal domain in R^{d}, d = 1, 2 or 3, and we let (·, ·) and ∥ · ∥ denote the usual inner product and norm in L_{2} = L_{2}(Ω). The norms in the standard Sobolev spaces H^{s} = H^{s}(Ω), s ≥ 0, are denoted by ∥ · ∥_{s}. H_{0}^{1} = H_{0}^{1}(Ω) is the space of functions v ∈ H^{1} satisfying the Dirichlet boundary condition v|_{∂Ω} = 0, and H^{−1} = H^{−1}(Ω) is the dual space of H_{0}^{1} with norm ∥v∥_{−1} = sup_{χ∈H_{0}^{1}} |(v, χ)|/∥χ∥_{1}. For u ∈ H_{0}^{1} we let B(ρ, u) denote the closed ball of radius ρ > 0 centered at u, i.e., B(ρ, u) = {v ∈ H_{0}^{1} : ∥v − u∥_{1} ≤ ρ}. The symbols C and C(ρ) are used to denote generic constants whose values may change from one occurrence to the next.

We consider the model semilinear parabolic problem

(2.1)  u_{t} − Δu = f(u), x ∈ Ω;  u = 0, x ∈ ∂Ω,

where u = u(x, t), Δu = ∑_{i=1}^{d} ∂^{2}u/∂x_{i}^{2}, and u_{t} = ∂u/∂t. We assume that f ∈ C^{2}(R) and, if d = 2 or 3, we assume in addition that

(2.2)  |f^{(j)}(u)| ≤ C(1 + |u|^{δ−j}),  j = 0, 1, 2, u ∈ R,

with δ = 3 if d = 3, and 2 ≤ δ < ∞ if d = 2. It then follows that f is locally Lipschitz from H_{0}^{1} into L_{2} (cf. Lemma 2.2 below), which, by standard techniques (see, e.g., [7]), implies local existence and uniqueness of solutions to the initial value problem for (2.1).

We further assume that (2.1) has a stationary solution ū ∈ H_{0}^{1} ∩ H^{2}, i.e., ū is a solution of

−Δū = f(ū), x ∈ Ω;  ū = 0, x ∈ ∂Ω.

Let A = −Δ − aI, where a(x) = f′(ū(x)), denote the linearized operator with domain D(A) = H_{0}^{1}∩H^{2} and let {λ_{i}}_{i=1}^{∞} be its eigenvalues numbered in nondecreasing order, with {ϕ_{i}}_{i=1}^{∞} the corresponding L_{2}-orthonormal eigenfunctions. Note that, by standard embeddings, ū and hence a are continuous functions in the closure Ω̄. The stationary point ū is assumed to be hyperbolic (i.e., we assume that 0 is not one of the eigenvalues λ_{i}) and unstable (i.e., there is at least one negative eigenvalue). Hence there is a positive integer q such that λ_{i} < 0 for 1 ≤ i ≤ q, and λ_{i} > 0 for q+1 ≤ i < ∞.

We then let P_{1} denote the orthogonal projection of L_{2} onto X_{1} = span{ϕ_{i}}_{i=1}^{q} and set P_{2} = I −P_{1}, X_{2} = X_{1}^{⊥} (orthogonal complement in L_{2}). Finally we define E(t) = e^{−tA}, A_{j} = A|_{X_{j}}, E_{j}(t) = e^{−tA_{j}} for j = 1, 2, i.e.,

E_{1}(t)P_{1}v = ∑_{i=1}^{q} e^{−tλ_{i}}(v, ϕ_{i})ϕ_{i}, E_{2}(t)P_{2}v = ∑_{i=q+1}^{∞} e^{−tλ_{i}}(v, ϕ_{i})ϕ_{i}, v ∈ L_{2}.

Thus E_{1} and E_{2}, respectively, correspond to the unstable and stable parts of the semigroup E generated by −A. The following lemma contains bounds for these operators. Although the result is valid in greater generality, we state it with precisely those assumptions that will be needed in the sequel.

Lemma 2.1. There are positive numbers M and β such that

\begin{matrix}∥E_{1}(t)P_{1}v∥_{γ} ≤ Me^{βt}∥v∥_{α}, & t ≤ 0, \\ ∥E_{2}(t)P_{2}v∥_{γ} ≤ Mt^{− \frac{γ−α}{2}} e^{−βt}∥v∥_{α}, & t > 0, \end{matrix}

for v ∈ H_{0}^{1}, α = 0, 1, γ = 0, 1, 2 with α ≤ γ.

Proof. Let B = −Δ with domain D(B) = D(A) = H_{0}^{1} ∩ H^{2}. We first note the equivalence of norms

(2.3)  ∥v∥_{γ} ≈ ∥A_{2}^{γ/2}v∥,  v ∈ X_{2} ∩ D(B^{γ/2}), γ = 0, 1, 2.

This is trivial when γ = 0. Since B − A = aI is bounded in L_{2} and ∥Av∥ = ∥A_{2}v∥ ≥ λ_{q+1}∥v∥ for v ∈ X_{2}, we obtain

∥Bv∥ ≤ ∥(B − A)v∥ + ∥Av∥ ≤ C∥Av∥, v ∈ X_{2} ∩ D(B),

and also an analogous inequality with A and B interchanged. In view of the equivalence of norms ∥Bv∥ ≈ ∥v∥_{2} for v ∈ D(B), this proves the case γ = 2 of (2.3). The remaining case γ = 1 follows by interpolation.

Since the spectrum of A_{2} is positive and bounded away from 0, we thus have

∥E_{2}(t)P_{2}v∥_{γ} ≤ Ct^{− \frac{γ−α}{2}} e^{−βt}∥P_{2}v∥_{α}, t > 0,

and the second bound in the lemma follows, if P_{2} is bounded with respect to the H^{α} norm. For α = 0 this is obvious. For α = 1 it follows from the orthogonality of P_{2} with respect to the indefinite bilinear form

(2.4)  A(u, v) = (∇u, ∇v) − (au, v),

(following a standard practice we use the same letter A to denote both the linear operator A and the corresponding bilinear form A(·, ·)) by the following estimation:

\begin{aligned}c∥P_{2}v∥_{1}^{2} & ≤ ∥A_{2}^{1/2}P_{2}v∥^{2} = A(P_{2}v, P_{2}v) \\ & = A(P_{2}v, v) ≤ C∥P_{2}v∥_{1}∥v∥_{1}, v ∈ H_{0}^{1}. \end{aligned}

Finally, the bound for E_{1}(t) is clear, since X_{1} is finite dimensional and A_{1} is bounded and negative definite. ⊔ ⊓

Using the linearized operator A we may now write equation (2.1) as

(2.5)  u_{t} + Au = F (u);  F (u) = f(u) − au.

We are interested in the behavior of solutions of (2.5) in the neighborhood of ū. It is therefore convenient to introduce a dependent variable z(t) = u(t) − ū, which satisfies

(2.6)  z_{t} + Az = G(z);  G(z) = F (u) − F (ū).

We shall need the fact that the mappings F, G : H_{0}^{1} → L_{2} are locally Lipschitz with a constant that can be made arbitrarily small in a neighborhood of ū and 0, respectively. More precisely, we have the following result.

Lemma 2.2. If u_{i} ∈ B(ρ, ū), and z, z_{i} ∈ B(ρ, 0) for i = 1, 2, then for j = 0, 1 we have the bounds

(2.7)  ∥F (u_{1}) − F (u_{2})∥_{−j} ≤ k(ρ)∥u_{1} − u_{2}∥_{1−j},
(2.8)  ∥G(z_{1}) − G(z_{2})∥_{−j} ≤ k(ρ)∥z_{1} − z_{2}∥_{1−j},
(2.9)  ∥G′(z)v∥_{−j} ≤ k(ρ)∥v∥_{1−j},
(2.10)  ∥G(z)∥_{1} ≤ k(ρ)∥z∥_{2},
(2.11)  ∥G′(z)v∥_{1} ≤ k(ρ)∥v∥_{2} + C(ρ)∥z∥_{2}∥v∥_{1},

where k(ρ) = O(ρ), C(ρ) = O(1) as ρ → 0.

Proof. Since F (u_{1}) − F (u_{2}) = G(z_{1}) − G(z_{2}) (with z = u − ū), it is sufficient to show the bounds for G. The inequality (2.8) follows readily from (2.9). We prove (2.9) for d = 3 only; the remaining cases can be proved in a similar way. We have, by Hölder’s and Sobolev’s inequalities,

∥G′(z)v∥ ≤ ∥G′(z)∥_{L_{3}}∥v∥_{L_{6}} ≤ C∥G′(z)∥_{L_{3}}∥v∥_{1}

(with a slight abuse of notation we let G′(z) denote both the linear operator G′(z) : H_{0}^{1} → L_{2} and the related function x ↦ G′(z(x)) in L_{3}(Ω)). Moreover, by (2.2) (recall that δ = 3 if d = 3),

∥G′(z)∥_{L_{3}} = ∥∫_{0}^{1} f′′(ū + tz)z dt∥_{L_{3}} ≤ ∫_{0}^{1} ∥f′′(ū + tz)∥_{L_{6}} dt ∥z∥_{L_{6}}
≤ C (1 + ∥ū∥_{L_{6}} + ∥z∥_{L_{6}}) ∥z∥_{L_{6}}
≤ C (1 + ∥ū∥_{1} + ∥z∥_{1}) ∥z∥_{1} ≤ C (1 + ρ)ρ,

which proves the case j = 0 of (2.9). The case j = 1 is obtained in a similar way, using the inequalities

∥v∥_{−1} = sup_{χ∈H_{0}^{1}} \frac{|(v, χ)|}{∥χ∥_{1}} ≤ sup_{χ∈H_{0}^{1}} \frac{∥v∥_{L_{6/5}}∥χ∥_{L_{6}}}{∥χ∥_{1}} ≤ C∥v∥_{L_{6/5}},

∥G′(z)v∥_{L_{6/5}} ≤ ∫_{0}^{1} ∥f′′(ū + tz)∥_{L_{6}} dt ∥z∥_{L_{6}} ∥v∥_{L_{2}}.

The remaining bounds (2.10) and (2.11) are proved by the same techniques. ⊔ ⊓ Let −∞ < τ < t < T < ∞. A solution of (2.5) defined for t ∈ [τ, T] satisfies

u(t) = E(t − τ)u(τ) + ∫_{τ}^{t} E(t − s)F (u(s)) ds,

and hence

P_{2}u(t) = E_{2}(t − τ)P_{2}u(τ) + ∫_{τ}^{t} E_{2}(t − s)P_{2}F (u(s)) ds.

Similarly,

P_{1}u(T) = E_{1}(T − t)P_{1}u(t) + ∫_{t}^{T} E_{1}(T − s)P_{1}F (u(s)) ds,

so that, by application of E_{1}(t − T), we get

P_{1}u(t) = E_{1}(t − T)P_{1}u(T) − ∫_{t}^{T} E_{1}(t − s)P_{1}F (u(s)) ds.

Substitution of these expressions into u(t) = P_{1}u(t) + P_{2}u(t) yields

\begin{aligned}u(t) & = E_{1}(t − T)v − ∫_{t}^{T} E_{1}(t − s)P_{1}F (u(s)) ds \\ & + E_{2}(t − τ)w + ∫_{τ}^{t} E_{2}(t − s)P_{2}F (u(s)) ds, \end{aligned} (2.12)

where v = P_{1}u(T), w = P_{2}u(τ). Conversely, a solution of (2.12) (see, e.g., [7] for the appropriate definition of solution) satisfies (2.5). The corresponding equation for a solution z of (2.6) is

\begin{aligned}z(t) & = E_{1}(t − T)P_{1}z(T) − ∫_{t}^{T} E_{1}(t − s)P_{1}G(z(s)) ds \\ & + E_{2}(t − τ)P_{2}z(τ) + ∫_{τ}^{t} E_{2}(t − s)P_{2}G(z(s)) ds. \end{aligned} (2.13)

Note also that ū satisfies (2.12), i.e.,

\begin{aligned}ū & = E_{1}(t − T)P_{1}ū − ∫_{t}^{T} E_{1}(t − s)P_{1}F (ū) ds \\ & + E_{2}(t − τ)P_{2}ū + ∫_{τ}^{t} E_{2}(t − s)P_{2}F (ū) ds. \end{aligned} (2.14)

The representation in (2.12) holds for bounded time intervals [τ, T]. If we let τ = 0 and T → ∞ in (2.12) under the assumption that ∥u(T)∥_{1} remains bounded, then we obtain the equation

(2.15) u(t) = E_{2}(t)w − ∫_{t}^{∞} E_{1}(t − s)P_{1}F (u(s)) ds + ∫_{0}^{t} E_{2}(t − s)P_{2}F (u(s)) ds,

which is satisfied by solutions of (2.5) defined in [0, ∞). Solutions defined in (−∞, 0] can also be accommodated by means of a similar device.

Consider now a solution u = u(t) of (2.5), which enters a small ball B(ρ, ū) at time τ and exits B(ρ, ū) at time T. As pointed out in the introduction, the initial value problem, where u(τ) is prescribed, is ill posed due to the unstable character of the stationary point ū. However, equation (2.12) shows that the boundary value problem, where P_{2}u(τ) and P_{1}u(T) are prescribed, is well posed. This is the key idea of the present work. We make this precise in the following lemma, where M is the constant in Lemma 2.1.

Lemma 2.3. There is a positive number ρ such that, for any real numbers τ, T with τ < T and any v ∈ X_{1}, w ∈ X_{2} with

(2.16)  ∥v − P_{1}ū∥_{1} + ∥w − P_{2}ū∥_{1} ≤ \frac{ρ}{2M} ,

equation (2.12) has a unique solution u such that u(t) ∈ B(ρ, ū) for t ∈ [τ, T].

Proof. We shall apply Banach’s fixed point theorem in the space C = C([τ, T] : H_{0}^{1}) normed by ∥u∥_{C} = sup_{τ≤t≤T} ∥u(t)∥_{1}. Equation (2.12) can be written as u = S(v, w) + T (u), where

(2.17)  \begin{matrix}S(v, w)(t) = E_{1}(t − T)v + E_{2}(t − τ)w, \\ T (u)(t) = − ∫_{t}^{T} E_{1}(t − s)P_{1}F (u(s)) ds + ∫_{τ}^{t} E_{2}(t − s)P_{2}F (u(s)) ds. \end{matrix}

Let B = {u ∈ C : ∥u − ū∥_{C} ≤ ρ}. We want to choose ρ such that the operator T is a contraction in B. For future reference we introduce

(2.18)  \begin{matrix}J(t) = ∫_{t}^{∞} e^{β(t−s)}(1 + s^{− \frac{1}{2}} ) ds + ∫_{0}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}(1 + s^{− \frac{1}{2}} ) ds, \\ K = sup_{t≥0} J(t). \end{matrix}

It is easy to show that K is finite, so that we may choose ρ > 0 such that

(2.19)  Mk(ρ)K ≤ \frac{1}{2}.

For u, z ∈ B we then have, by Lemma 2.1 and Lemma 2.2,

∥T (u)(t) − T (z)(t)∥_{1} ≤ ∫_{t}^{T} ∥E_{1}(t − s)P_{1}[F (u(s)) − F (z(s))]∥_{1} ds
+ ∫_{τ}^{t} ∥E_{2}(t − s)P_{2}[F (u(s)) − F (z(s))]∥_{1} ds
≤ M ∫_{t}^{T} e^{β(t−s)}∥F (u(s)) − F (z(s))∥ ds
+ M ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}∥F (u(s)) − F (z(s))∥ ds

≤ Mk(ρ)(∫_{t}^{T} e^{β(t−s)} ds + ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)} ds)∥u − z∥_{C}

≤ Mk(ρ)J(t − τ)∥u − z∥_{C} ≤ Mk(ρ)K∥u − z∥_{C}
≤ \frac{1}{2}∥u − z∥_{C}, t ∈ [τ, T],

so that T is a contraction. It remains to check that the operator u ↦ S(v, w) + T (u) maps B into itself. This will be achieved if ∥S(v, w) − ū + T (ū)∥_{C} ≤ \frac{1}{2}ρ, or

∥S(v, w) − S(P_{1}ū, P_{2}ū)∥_{C} ≤ \frac{1}{2}ρ,

since in view of (2.14) ū = S(P_{1}ū, P_{2}ū)+T (ū). But, in view of Lemma 2.1 and (2.16), we have

\begin{aligned}∥S(v, w)(t) − S(P_{1}ū, P_{2}ū)(t)∥_{1} & ≤ ∥E_{1}(t − T)(v − P_{1}ū)∥_{1} + ∥E_{2}(t − τ)(w − P_{2}ū)∥_{1} \\ & ≤ M(∥v − P_{1}ū∥_{1} + ∥w − P_{2}ū∥_{1}) ≤ \frac{1}{2}ρ, \end{aligned}

for t ∈ [τ, T]. Thus equation (2.12) has a unique solution u ∈ B. ⊔ ⊓

Note that the lemma also holds, with essentially the same proof, for solutions on the semi-infinite time interval [0, ∞). More precisely, there is ρ > 0 such that equation (2.15) has a unique solution u(t) ∈ B(ρ, ū) for t ∈ [0, ∞) for all w ∈ X_{2} with ∥w − P_{2}ū∥_{1} ≤ ρ/2M. This implies that the local stable manifold of ū, defined by

(2.20) M_{S}(ρ) = \{u_{0} ∈ H_{0}^{1} : ∥P_{2}(u_{0} − ū)∥_{1} ≤ \frac{ρ}{2M} and u(t; u_{0}) ∈ B(ρ, ū) for t ≥ 0\},

is homeomorphic to the ball B(ρ/2M, P_{2}ū) in X_{2}. Here u(t; u_{0}) denotes the solution u(t) of (2.5) satisfying u(0) = u_{0}. Moreover, one can show that M_{S}(ρ) is tangent to X_{2} at ū, and that u(t; u_{0}) → ū exponentially in H_{0}^{1} as t → ∞ for all u_{0} ∈ M_{S}(ρ). Similar considerations hold for solutions defined in (−∞, 0], leading to the construction of the local unstable manifold M_{U}(ρ). We refer to Henry [7, Theorem 5.2.1] for the details.

Our next lemma concerns the regularity of a solution u(t) ∈ B(ρ, ū) for t ∈ [τ, T]. In our error analysis below we will need bounds for certain derivatives of u. Since T − τ may be arbitrarily large, it is crucial that these bounds are independent of τ and T.

Lemma 2.4. There are positive numbers ρ and C such that, for any τ, T with τ < T, and for any solution u of (2.12) with u(t) ∈ B(ρ, ū) for t ∈ [τ, T], we have the bounds

(2.21)  ∥D_{t}^{l}u(t)∥_{m} ≤ C(1 + (t − τ)^{−l− \frac{m−1}{2}} ),  t ∈ (τ, T],

for l = 0, 1, m = 1, 2.

Proof. Let ρ be given by (2.19). It is convenient to estimate z(t) = u(t)−ū, which satisfies (2.6) and (2.13). The desired bounds for u then follow, since ū ∈ H_{0}^{1} ∩ H^{2}. The reason why this is convenient is that G(z) ∈ H_{0}^{1}, so that we may employ Lemma 2.1 with v = G(z) and α = 1. This is not possible with F (u), which may be nonzero on ∂Ω.

Using Lemma 2.1, (2.10), (2.19) and equation (2.13) we get

∥z(t)∥_{2} ≤ ∥E_{1}(t − T)P_{1}z(T)∥_{2} + ∥E_{2}(t − τ)P_{2}z(τ)∥_{2}

+ ∫_{t}^{T} ∥E_{1}(t − s)P_{1}G(z(s))∥_{2} ds + ∫_{τ}^{t} ∥E_{2}(t − s)P_{2}G(z(s))∥_{2} ds

≤ M∥z(T)∥_{1} + M(t − τ)^{− \frac{1}{2}} ∥z(τ)∥_{1}

+ M ∫_{t}^{T} e^{β(t−s)}∥G(z(s))∥_{1} ds + M ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}∥G(z(s))∥_{1} ds

≤ Mρ(1 + (t − τ)^{− \frac{1}{2}} ) + Mk(ρ)(∫_{t}^{∞} e^{β(t−s)}(1 + (s − τ)^{− \frac{1}{2}} ) ds

+ ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}(1 + (s − τ)^{− \frac{1}{2}} ) ds) sup_{τ<s<T}(φ(s − τ)∥z(s)∥_{2})

≤ Mρ(1 + (t − τ)^{− \frac{1}{2}} ) + Mk(ρ)J(t − τ) sup_{τ<s<T}(φ(s − τ)∥z(s)∥_{2})

≤ Mρ(1 + (t − τ)^{− \frac{1}{2}} ) + Mk(ρ)K sup_{τ<s<T}(φ(s − τ)∥z(s)∥_{2})

≤ Mρ(1 + (t − τ)^{− \frac{1}{2}} ) + \frac{1}{2} sup_{τ<s<T}(φ(s − τ)∥z(s)∥_{2}),

where φ(s) = (1 + s^{− \frac{1}{2}} )^{−1} ≤ 1. Hence

φ(t − τ)∥z(t)∥_{2} ≤ Mρ + \frac{1}{2} sup_{τ<s<T}(φ(s − τ)∥z(s)∥_{2}),

so that

(2.22)  ∥z(t)∥_{2} ≤ C(ρ)(1 + (t − τ)^{− \frac{1}{2}} ),

which implies the case l = 0, m = 2 of (2.21). The case l = 0, m = 1 is trivial.

We now estimate v = u_{t} = z_{t}. It can be shown by standard techniques that u is differentiable with respect to t and that v satisfies v_{t} + Av = G′(z)v, and hence

\begin{aligned}v(t) & = E_{1}(t − T)P_{1}v(T) − ∫_{t}^{T} E_{1}(t − s)P_{1}(G′(z(s))v(s)) ds \\ & + E_{2}(t − τ)P_{2}v(τ) + ∫_{τ}^{t} E_{2}(t − s)P_{2}(G′(z(s))v(s)) ds. \end{aligned}

We refer to Henry [7, Theorem 3.5.2] for the details. It remains to prove the required bounds for v with a constant C(ρ) independent of τ and T. (The bound given in [7] is obtained by Gronwall’s lemma and therefore not applicable in the present situation, where we require bounds that do not deteriorate as T − τ → ∞.) We already have, by (2.10) and (2.22),

(2.23)  ∥v(t)∥ ≤ ∥Az(t)∥ + ∥G(z(t))∥ ≤ C(ρ)∥z(t)∥_{2} ≤ C(ρ)(1 + (t − τ)^{− \frac{1}{2}} ).

Using Lemma 2.1, (2.9) and (2.18), we obtain

∥v(t)∥_{1} ≤ ∥E_{1}(t − T)P_{1}v(T)∥_{1} + ∥E_{2}(t − τ)P_{2}v(τ)∥_{1}

+ ∫_{t}^{T} ∥E_{1}(t − s)P_{1}(G′(z(s))v(s))∥_{1} ds

+ ∫_{τ}^{t}∥E_{2}(t − s)P_{2}(G′(z(s))v(s))∥_{1} ds

≤ M∥v(T)∥ + M(t − τ)^{− \frac{1}{2}} ∥v(τ)∥

+ Mk(ρ)(∫_{t}^{T} e^{β(t−s)}∥v(s)∥_{1} ds + ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}∥v(s)∥_{1} ds)

≤ M(1 + (t − τ)^{− \frac{1}{2}} )(∥v(T)∥ + ∥v(τ)∥)

+ Mk(ρ)K sup_{τ<s<T}(φ(s − τ)∥v(s)∥_{1}).

In view of (2.19) this shows that

φ(t − τ)∥v(t)∥_{1} ≤ 2M(∥v(T)∥ + ∥v(τ)∥),

which together with (2.23) with τ replaced by (t − τ)/2 implies

(2.24)  ∥v(t)∥_{1} ≤ C(ρ)(1 + (t − τ)^{−1}),

which is the case l = 1, m = 1 of (2.21). Similarly we may show

(2.25)  ∥v(t)∥_{1} ≤ 2M(∥v(T)∥_{1} + ∥v(τ)∥_{1}).

For the proof of the remaining bound (l = 1, m = 2) we use (2.11), (2.22) and (2.25) in a similar way as above:

∥v(t)∥_{2} ≤ M(1 + (t − τ)^{− \frac{1}{2}} )(∥v(T)∥_{1} + ∥v(τ)∥_{1})

+ Mk(ρ)K sup_{τ<s<T}(φ(s − τ)∥v(s)∥_{2})

+ MC(ρ)(∫_{t}^{T} e^{β(t−s)}∥z(s)∥_{2}∥v(s)∥_{1} ds

+ ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}∥z(s)∥_{2}∥v(s)∥_{1} ds)

≤ M(1 + (t − τ)^{− \frac{1}{2}} + 2MC(ρ)K)(∥v(T)∥_{1} + ∥v(τ)∥_{1})

+ \frac{1}{2} sup_{τ<s<T}(φ(s − τ)∥v(s)∥_{2}).

Together with (2.24) this leads to ∥v(t)∥_{2} ≤ C(ρ)(1 + (t − τ)^{− \frac{3}{2}} ) and the proof is complete. ⊔ ⊓

3. The semidiscrete problem. We now proceed to discuss the spatially semidiscrete approximation of equation (2.1) by the standard piecewise linear finite element method. Thus we denote by S^{h} the subspace of H_{0}^{1} that consists of piecewise polynomials of degree ≤ 1 with respect to a “triangulation” of the convex polygonal domain Ω with maximum mesh size h. The semidiscrete solutions u^{h}(t) ∈ S^{h} satisfy

(3.1)  (u_{t}^{h}, χ) + (∇u^{h}, ∇χ) = (f(u^{h}), χ)  ∀χ ∈ S^{h}.

Before pursuing the discussion of u^{h} further, we collect some basic results and assumptions concerning the finite element method. Let α = max_{x∈Ω} a(x), so that the bilinear form

Ã(ψ, χ) = (∇ψ, ∇χ) + ((α − a)ψ, χ)

is H_{0}^{1}-elliptic. Then there is an operator G : L_{2} → H_{0}^{1} such that

(3.2)  Ã(Gf, χ) = (f, χ)  ∀χ ∈ H_{0}^{1}, f ∈ L_{2}.

In other words G = (A + αI)^{−1}. It is easy to see that G is a selfadjoint, positive definite, compact, linear operator on L_{2}, and by the standard regularity theory for elliptic problems, we have the inequality

(3.3)  ∥Gf∥_{2} ≤ C∥f∥,  f ∈ L_{2}.

We also define an “elliptic projection” operator R^{h} : H_{0}^{1} → S^{h} by

(3.4)  Ã(R^{h}v − v, χ) = 0  ∀χ ∈ S^{h}.

Under the usual regularity assumptions on the triangulation, and in view of (3.3), the standard error analysis for elliptic problems yields

(3.5)  ∥R^{h}v − v∥ + h∥R^{h}v − v∥_{1} ≤ Ch^{m}∥v∥_{m},  v ∈ H_{0}^{1} ∩ H^{m}, m = 1, 2,

see for example [4, Chapter 3]. This error bound will form the basis for our error analysis. For convenience of the presentation of our main result (cf. Remark 2 below) we assume, in addition, that the triangulation is such that the orthogonal projection P^{h} of L_{2} onto S^{h} is bounded (uniformly in h) with respect to the H^{1} norm. It is easy to see that this is true under an inverse assumption. For a more general discussion of the H^{1} boundedness of P^{h} the reader may consult [5].

We now resume our discussion of the semidiscrete problem (3.1). It is convenient to linearize the equation around ū. We thus introduce the operator A^{h} : S^{h} → S^{h}, defined by

(A^{h}v, χ) = A(v, χ) ∀χ ∈ S^{h},

where A(·, ·) is the bilinear form corresponding to the operator A, see (2.4). Equation (3.1) now becomes

(3.6)  u_{t}^{h} + A^{h}u^{h} = P^{h}F (u^{h});  F (u) = f(u) − au.

Let {λ_{i}^{h}}_{i=1}^{N^{h}} be the eigenvalues of A^{h} numbered in nondecreasing order, and let {ϕ_{i}^{h}}_{i=1}^{N^{h}} be the corresponding L_{2}-orthonormal eigenfunctions. It is well known [2] that the eigenvalues of A^{h} converge to those of A as h → 0. Hence there is h_{0} > 0 such that, for h < h_{0}, we have λ_{i}^{h} < 0 for 1 ≤ i ≤ q, and λ_{i}^{h} > 0 for q + 1 ≤ i ≤ N^{h}. We then set X_{1}^{h} = span{ϕ_{i}^{h}}_{i=1}^{q}, X_{2}^{h} = span{ϕ_{i}^{h}}_{i=q+1}^{N^{h}}, and denote by P_{j}^{h} the orthogonal projections of L_{2} onto X_{j}^{h} for j = 1, 2. Note that P_{j}^{h} = P_{j}^{h}P^{h} and P_{1}^{h} = P^{h} − P_{2}^{h}.

Next we define A_{j}^{h} = A^{h}|_{X_{j}^{h}} and evolution operators E_{j}^{h}(t) = e^{−tA_{j}^{h}}, i.e.,

E_{1}^{h}(t)P_{1}^{h}v = ∑_{i=1}^{q} e^{−tλ_{i}^{h}}(v, ϕ_{i}^{h})ϕ_{i}^{h}, E_{2}^{h}(t)P_{2}^{h}v = ∑_{i=q+1}^{N^{h}} e^{−tλ_{i}^{h}}(v, ϕ_{i}^{h})ϕ_{i}^{h}, v ∈ L_{2}.

The next result provides bounds for these operators. Note that there is no loss of generality in assuming that the constants M and β in the following lemma are the same as in Lemma 2.1.

Lemma 3.1. There are positive numbers h_{0}, M and β such that for h < h_{0} we have

\begin{matrix}∥E_{1}^{h}(t)P_{1}^{h}v∥_{1} ≤ Me^{βt}∥P_{1}^{h}v∥_{1}, & t ≤ 0, \\ ∥E_{2}^{h}(t)P_{2}^{h}v∥_{1} ≤ Me^{−βt}∥P_{2}^{h}v∥_{1}, & t ≥ 0, \end{matrix}

and

\begin{matrix}∥D_{t}^{j}E_{1}^{h}(t)P_{1}^{h}v∥_{γ} ≤ Me^{βt}∥v∥_{α}, & t ≤ 0, \\ ∥D_{t}^{j}E_{2}^{h}(t)P_{2}^{h}v∥_{γ} ≤ Mt^{−j− \frac{γ−α}{2}} e^{−βt}∥v∥_{α}, & t > 0, \end{matrix}

for v ∈ H_{0}^{1}, j = 0, 1, α = −1, 0, γ = 0, 1.

Proof. We first establish the equivalence of norms (uniform in h)

(3.7)  ∥v∥_{1} ≈ ∥(A_{2}^{h})^{1/2}v∥,  v ∈ X_{2}^{h}.

Since λ_{q+1}^{h} ≥ λ_{q+1}, we have

(3.8)  A(v, v) ≥ λ_{q+1}∥v∥^{2},  v ∈ X_{2}^{h}.

Hence

(3.9) ∥v∥_{1}^{2} ≤ C∥∇v∥^{2} = C(A(v, v) + (av, v)) ≤ CA(v, v) ≤ C∥v∥_{1}^{2},  v ∈ X_{2}^{h},

which is the desired result, in view of the identity ∥(A_{2}^{h})^{1/2}v∥^{2} = A(v, v). Using (3.7) and (3.8) we may now prove

∥D_{t}^{j}E_{2}^{h}(t)P_{2}^{h}v∥_{γ} ≤ Ct^{−j− \frac{γ−α}{2}} e^{−βt}∥(A_{2}^{h})^{α/2}P_{2}^{h}v∥, t > 0.

This implies the second bound in the lemma, because of the above equivalence of norms, and also the fourth bound for α = 0, because of the boundedness of P_{2}^{h} in L_{2}. For the case when α = −1 we note that

(3.10)  A((A_{2}^{h})^{−1}P_{2}^{h}v, χ) = (v, χ)  ∀χ ∈ X_{2}^{h},

which together with (3.9) implies

∥(A_{2}^{h})^{−1/2}P_{2}^{h}v∥ = sup_{χ∈X_{2}^{h}} \frac{(v, χ)}{A(χ, χ)^{\frac{1}{2}}} ≤ C sup_{χ∈X_{2}^{h}} \frac{(v, χ)}{∥χ∥_{1}} ≤ C sup_{χ∈H_{0}^{1}} \frac{(v, χ)}{∥χ∥_{1}} = C∥v∥_{−1}.

This completes the proof of the bounds for E_{2}^{h}(t). The bounds for E_{1}^{h}(t) are easily obtained by noting that X_{1}^{h} is finite dimensional with dimension q, and that the spectrum of A_{1}^{h} is negative and bounded away from 0 (uniformly in h). ⊔ ⊓

Equation (3.6) is now equivalent to

\begin{aligned}u^{h}(t) & = E_{1}^{h}(t − T)v^{h} − ∫_{t}^{T} E_{1}^{h}(t − s)P_{1}^{h}F (u^{h}(s)) ds \\ & + E_{2}^{h}(t − τ)w^{h} + ∫_{τ}^{t} E_{2}^{h}(t − s)P_{2}^{h}F (u^{h}(s)) ds, \end{aligned} (3.11)

where v^{h} = P_{1}^{h}u^{h}(T), w^{h} = P_{2}^{h}u^{h}(τ), cf. (2.12). Again the cases T = ∞ and τ = −∞ can be catered for after minor changes, cf. (2.15).

In our next lemma we will show that equation (3.6) has a unique stationary solution in a neighborhood of ū. In other words, there is ū^{h} ∈ S^{h} satisfying

(3.12)  A^{h}ū^{h} = P^{h}F (ū^{h}),

which is implies

\begin{aligned}ū^{h} & = E_{1}^{h}(t − T)P_{1}^{h}ū^{h} − ∫_{t}^{T} E_{1}^{h}(t − s)P_{1}^{h}F (ū^{h}) ds \\ & + E_{2}^{h}(t − τ)P_{2}^{h}ū^{h} + ∫_{τ}^{t} E_{2}^{h}(t − s)P_{2}^{h}F (ū^{h}) ds. \end{aligned} (3.13)

Lemma 3.2. There are positive numbers ρ, h_{0} and C such that, for h < h_{0}, equation (3.12) has a unique solution ū^{h} ∈ B(ρ, ū) ∩ S^{h}. Moreover, ∥ū^{h} − ū∥_{1} ≤ Ch.

Proof. We first note that the operator A^{h} is invertible and that there is a constant C_{1} such that

∥(A^{h})^{−1}P^{h}v∥_{1} ≤ C_{1}∥v∥, v ∈ L_{2},

if h is sufficiently small. To see this we use (3.9) and (3.10) with χ = (A_{2}^{h})^{−1}P_{2}^{h}v to find

\begin{aligned}c∥(A_{2}^{h})^{−1}P_{2}^{h}v∥_{1}^{2} & ≤ A((A_{2}^{h})^{−1}P_{2}^{h}v, (A_{2}^{h})^{−1}P_{2}^{h}v) = (v, (A_{2}^{h})^{−1}P_{2}^{h}v) \\ & ≤ ∥v∥ ∥(A_{2}^{h})^{−1}P_{2}^{h}v∥, \end{aligned}

which shows ∥(A_{2}^{h})^{−1}P_{2}^{h}v∥_{1} ≤ C∥v∥. The inequality ∥(A_{1}^{h})^{−1}P_{1}^{h}v∥_{1} ≤ C∥v∥ holds by finite dimensionality (uniform in h). This proves the desired bound.

In view of (3.12) we are seeking a fixed point of the operator Q^{h} defined by Q^{h}(u) = (A^{h})^{−1}P^{h}F (u). For u, v ∈ B(ρ, ū) we have, in view of (2.7),

\begin{aligned}∥Q^{h}(u) − Q^{h}(v)∥_{1} & = ∥(A^{h})^{−1}P^{h}(F (u) − F (v))∥_{1} ≤ C_{1}∥F (u) − F (v)∥ \\ & ≤ C_{1}k(ρ)∥u − v∥_{1}. \end{aligned}

Hence Q^{h} is a contraction on B(ρ, ū), if we choose ρ so that C_{1}k(ρ) ≤ \frac{1}{2}. It remains to show that Q^{h} maps B(ρ, ū) into itself, if h is sufficiently small. Using the identity

(3.14)  A^{h}R^{h} = P^{h}(A − α(R^{h} − I)),

which follows easily from the definitions of A^{h} and R^{h}, see (3.4), we get

\begin{aligned}A^{h}(Q^{h}(u) − R^{h}ū) & = P^{h}F (u) − A^{h}R^{h}ū = P^{h}(F (u) − Aū + α(R^{h}ū − ū)) \\ & = P^{h}(F (u) − F (ū) + α(R^{h}ū − ū)), \end{aligned}

so that

Q^{h}(u) − R^{h}ū = Q^{h}(u) − Q^{h}(ū) + α(A^{h})^{−1}P^{h}(R^{h}ū − ū).

In view of (3.5) we thus have

(3.15) ∥Q^{h}(u) − ū∥_{1} ≤ \frac{1}{2}∥u − ū∥_{1} + (1 + |α|C_{1})∥R^{h}ū − ū∥_{1} ≤ \frac{1}{2}∥u − ū∥_{1} + Ch,

so that ∥Q^{h}(u) − ū∥_{1} ≤ \frac{1}{2}ρ + Ch ≤ ρ for u ∈ B(ρ, ū), if h is sufficiently small. Hence Q^{h} maps B(ρ, ū) into itself, and we conclude that Q^{h} has a unique fixed point ū^{h} ∈ B(ρ, ū). The fixed point belongs to S^{h} because the range of Q^{h} lies in S^{h}. The error bound follows immediately from (3.15). ⊔ ⊓

We can now prove an existence result for equation (3.11), analogous to Lemma 2.3 for the continuous problem.

Lemma 3.3. There are positive numbers ρ and h_{0} such that, for any h < h_{0}, for any real numbers τ, T with τ < T and for any v^{h} ∈ X_{1}^{h}, w^{h} ∈ X_{2}^{h} with

(3.16)  ∥v^{h} − P_{1}^{h}ū^{h}∥_{1} + ∥w^{h} − P_{2}^{h}ū^{h}∥_{1} ≤ \frac{ρ}{4M} ,

equation (3.11) has a unique solution u^{h} such that u^{h}(t) ∈ B(ρ, ū) for t ∈ [τ, T].

Proof. The argument is a slight modification of the proof of Lemma 2.3. Equation (3.11) is written as a fixed point equation u^{h} = S^{h}(v^{h}, w^{h}) + T^{h}(u^{h}), where the operators on the right hand side are defined in the obvious way, cf. (2.17). We want to solve this equation in the same ball B as in the proof of Lemma 2.3. Clearly T^{h} is a contraction on B with the same choice of ρ. It remains to check that the operator u ↦ S^{h}(v^{h}, w^{h}) + T^{h}(u) maps B into itself for small h. In view of (3.13), (3.16) and Lemma 3.2, this follows from

∥S^{h}(v^{h}, w^{h})(t) + T^{h}(ū) − ū∥_{1}
= ∥S^{h}(v^{h}, w^{h})(t) − S^{h}(P_{1}^{h}ū^{h}, P_{2}^{h}ū^{h})(t) + T^{h}(ū) − T^{h}(ū^{h}) − (ū − ū^{h})∥_{1}
≤ ∥E_{1}^{h}(t − T)(v^{h} − P_{1}^{h}ū^{h})∥_{1} + ∥E_{2}^{h}(t − τ)(w^{h} − P_{2}^{h}ū^{h})∥_{1} + \frac{3}{2}∥ū^{h} − ū∥_{1}
≤ M(∥v^{h} − P_{1}^{h}ū^{h}∥_{1} + ∥w^{h} − P_{2}^{h}ū^{h}∥_{1}) + Ch ≤ \frac{1}{4}ρ + Ch ≤ \frac{1}{2}ρ,

for t ∈ [τ, T]. This completes the proof. ⊔ ⊓

A variant of Lemma 3.3 can be used to construct the local unstable and stable manifolds M_{U}^{h}(ρ) and M_{S}^{h}(ρ) of ū^{h} in the same way as for the continuous problem. In this context we define (cf. (2.20))

(3.17) M_{S}^{h}(ρ) = \{u_{0}^{h} ∈ S^{h} : ∥P_{2}^{h}(u_{0}^{h} − ū^{h})∥_{1} ≤ \frac{ρ}{4M} , u^{h}(t; u_{0}^{h}) ∈ B(ρ, ū) for t ≥ 0\}.

In our next result we estimate the difference between two solutions u^{h}(t) and u(t) that remain in a small ball B(ρ, ū) for t ∈ [τ, T]. It is important that the error constant is independent of τ and T, because T − τ may be arbitrarily large. In the first error estimate there is a weak singularity at t = τ due to the possible lack of regularity of u(τ); note that we assume only that u(τ) ∈ H_{0}^{1}. The second estimate holds uniformly as t → τ, but the rate of convergence is correspondingly lower.

Lemma 3.4. There are positive numbers ρ, h_{0} and C such that, for any h < h_{0}, for any τ, T with τ < T, for any solutions u^{h} of (3.6) and u of (2.5) with u^{h}(t), u(t) ∈ B(ρ, ū) for t ∈ [τ, T], and for j = 0, 1, t ∈ (τ, T] we have

\begin{aligned}∥u^{h}(t) − u(t)∥_{j} & ≤ C(1 + (t − τ)^{− \frac{1}{2}} ) \\ & × (∥P_{1}^{h}(u^{h}(T) − u(T))∥ + ∥P_{2}^{h}(u^{h}(τ) − u(τ))∥ + h^{2−j}), \end{aligned}

and

∥u^{h}(t) − u(t)∥ ≤ C(∥P_{1}^{h}(u^{h}(T) − u(T))∥ + ∥P_{2}^{h}(u^{h}(τ) − u(τ))∥ + h).

Proof. Choose ρ as in (2.19), so that the conclusion of Lemma 2.4 holds. Following a standard practice we write

e(t) ≡ u^{h}(t) − u(t) = (u^{h}(t) − R^{h}u(t)) + (R^{h}u(t) − u(t)) ≡ θ^{h}(t) + ρ^{h}(t),

where R^{h} is the projection defined in (3.4), so that by (3.5) and Lemma 2.4, we have

(3.18)  ∥D_{t}^{l}ρ^{h}(t)∥_{j} ≤ Ch^{m−j}∥D_{t}^{l}u(t)∥_{m} ≤ Ch^{m−j}(1 + (t − τ)^{−l− \frac{m−1}{2}} ),

for l, j = 0, 1, m = 1, 2. In particular, we have the required bounds for ρ^{h}(t). In order to estimate θ^{h}(t) = u^{h}(t) − R^{h}u(t) we note that, in view of (3.14) and the equations (3.6) and (2.5) satisfied by u^{h} and u,

θ_{t}^{h} + A^{h}θ^{h} = P^{h}(F (u^{h}) − F (u) + αρ^{h} − ρ_{t}^{h}),

and hence that

θ^{h}(t) = E_{1}^{h}(t − T)P_{1}^{h}θ^{h}(T) + E_{2}^{h}(t − τ)P_{2}^{h}θ^{h}(τ)
− ∫_{t}^{T} E_{1}^{h}(t − s)P_{1}^{h}(F (u^{h}(s)) − F (u(s)) + αρ^{h}(s) − ρ_{t}^{h}(s)) ds
+ ∫_{τ}^{t} E_{2}^{h}(t − s)P_{2}^{h}(F (u^{h}(s)) − F (u(s)) + αρ^{h}(s) − ρ_{t}^{h}(s)) ds.

It is convenient to divide this expression into two parts θ^{h} = θ_{1}^{h} + θ_{2}^{h}, where

\begin{aligned}θ_{1}^{h}(t) & = E_{1}^{h}(t − T)P_{1}^{h}θ^{h}(T) + E_{2}^{h}(t − τ)P_{2}^{h}θ^{h}(τ) \\ & + ∫_{t}^{T} E_{1}^{h}(t − s)P_{1}^{h}ρ_{t}^{h}(s) ds − ∫_{τ}^{t} E_{2}^{h}(t − s)P_{2}^{h}ρ_{t}^{h}(s) ds, \end{aligned}

and

\begin{aligned}θ_{2}^{h}(t) & = − ∫_{t}^{T} E_{1}^{h}(t − s)P_{1}^{h}(F (u^{h}(s)) − F (u(s)) + αρ^{h}(s)) ds \\ & + ∫_{τ}^{t} E_{2}^{h}(t − s)P_{2}^{h}(F (u^{h}(s)) − F (u(s)) + αρ^{h}(s)) ds. \end{aligned}

We rewrite θ_{1}^{h} by integration by parts,

θ_{1}^{h}(t) = E_{1}^{h}(t − T)P_{1}^{h}e(T) + E_{2}^{h}(t − τ)P_{2}^{h}e(τ) − P_{1}^{h}ρ^{h}(t)
− E_{2}^{h}((t − τ)/2)P_{2}^{h}ρ^{h}((t − τ)/2) − ∫_{t}^{T} D_{s}E_{1}^{h}(t − s)P_{1}^{h}ρ^{h}(s) ds
+ ∫_{τ}^{(t−τ)/2} D_{s}E_{2}^{h}(t − s)P_{2}^{h}ρ^{h}(s) ds − ∫_{(t−τ)/2}^{t} E_{2}^{h}(t − s)P_{2}^{h}ρ_{t}^{h}(s) ds.

Hence, by Lemma 3.1 and (3.18) with j = 0, m = 1,

∥θ_{1}^{h}(t)∥_{1} ≤ M(1 + (t − τ)^{− \frac{1}{2}} )(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + ∥ρ^{h}(t)∥ + ∥ρ^{h}((t − τ)/2)∥)
+ M ∫_{t}^{T} e^{β(t−s)}∥ρ^{h}(s)∥ ds + M ∫_{τ}^{(t−τ)/2}(t − s)^{− \frac{3}{2}} e^{−β(t−s)}∥ρ^{h}(s)∥ ds
+ M ∫_{(t−τ)/2}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}∥ρ_{t}^{h}(s)∥ ds
≤ C(1 + (t − τ)^{− \frac{1}{2}} )(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + h)
+ Ch(∫_{t}^{T} e^{β(t−s)} ds + ∫_{τ}^{(t−τ)/2}(t − s)^{− \frac{3}{2}} e^{−β(t−s)} ds
+ ∫_{(t−τ)/2}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}(1 + (s − τ)^{−1}) ds)
≤ C(1 + (t − τ)^{− \frac{1}{2}} )(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + h).

For the remaining term θ_{2}^{h} we apply (2.7) with j = 0 to get

∥θ_{2}^{h}(t)∥_{1} ≤ M ∫_{t}^{T} e^{β(t−s)}(k(ρ)∥e(s)∥_{1} + C∥ρ^{h}(s)∥) ds

+ M ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}(k(ρ)∥e(s)∥_{1} + C∥ρ^{h}(s)∥) ds

≤ Mk(ρ)K sup_{τ<s<T}(φ(s − τ)∥θ^{h}(s)∥_{1}) + CKh,

where φ(s) = (1+s^{− \frac{1}{2}} )^{−1} and K is defined in (2.18). Since Mk(ρ)K ≤ \frac{1}{2} and φ(s) ≤ 1, we may conclude that

φ(t − τ)∥θ^{h}(t)∥_{1} ≤ C(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + h) + \frac{1}{2} sup_{τ<s<T}(φ(s − τ)∥θ^{h}(s)∥_{1}),

and the H^{1} norm error estimate follows. For the L_{2} norm estimates we first use (3.18) with j = 0, m = 1, 2 to get

∥θ_{1}^{h}(t)∥ ≤ M(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + ∥ρ^{h}(t)∥ + ∥ρ^{h}((t − τ)/2)∥)

+ M ∫_{t}^{T} e^{β(t−s)}∥ρ^{h}(s)∥ ds + M ∫_{τ}^{(t−τ)/2}(t − s)^{−1}e^{−β(t−s)}∥ρ^{h}(s)∥ ds

+ M ∫_{(t−τ)/2}^{t} e^{−β(t−s)}∥ρ_{t}^{h}(s)∥ ds
≤ C(1 + (t − τ)^{− \frac{m−1}{2}} )(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + h^{m})
+ Ch^{m}(∫_{t}^{T} e^{β(t−s)}(1 + (s − τ)^{− \frac{m−1}{2}} ) ds
+ ∫_{τ}^{(t−τ)/2}(t − s)^{−1}e^{−β(t−s)}(1 + (s − τ)^{− \frac{m−1}{2}} ) ds
+ ∫_{(t−τ)/2}^{t} e^{−β(t−s)}(1 + (s − τ)^{−1− \frac{m−1}{2}} ) ds)
≤ C(1 + (t − τ)^{− \frac{m−1}{2}} )(∥P_{1}^{h}e(T)∥ + ∥P_{2}^{h}e(τ)∥ + h^{m}).

For the estimate of ∥θ_{2}^{h}(t)∥ we use (2.7) with j = 1 and obtain

∥θ_{2}^{h}(t)∥ ≤ M ∫_{t}^{T} e^{β(t−s)}(k(ρ)∥e(s)∥ + C∥ρ^{h}(s)∥) ds

+ M ∫_{τ}^{t}(t − s)^{− \frac{1}{2}} e^{−β(t−s)}(k(ρ)∥e(s)∥ + C∥ρ^{h}(s)∥) ds

≤ Mk(ρ)K sup_{τ<s<T}(φ(s − τ)∥θ^{h}(s)∥) + CKh^{m},

where now φ(s) = (1 + s^{− \frac{m−1}{2}} )^{−1} and K is defined in (2.18) as before. The proof can now be completed in the same way as above. ⊔ ⊓

We may also, with essentially the same argument, prove analogous error bounds over the semi-infinite time interval [0, ∞), that is, bounds of the difference between solutions of (2.15) and its discrete version. For example, we have

(3.19) ∥u^{h}(t) − u(t)∥ ≤ C(1 + t^{− \frac{m−1}{2}} )(∥P_{2}^{h}(u^{h}(0) − u(0))∥ + h^{m}),  t ∈ [0, ∞),

for m = 1, 2, whenever u^{h}(t), u(t) ∈ B(ρ, ū) for t ∈ [0, ∞). Error bounds over (−∞, 0] are obtained in a similar way. In fact, we have

(3.20)  ∥u^{h}(t) − u(t)∥ ≤ C(∥P_{1}^{h}(u^{h}(0) − u(0))∥ + h^{2}),  t ∈ (−∞, 0],

since there is no initial singularity in this case.

As a final preparation for the proof of our main result we give a bound for the difference between the projections P_{1}^{h} and P_{1}.

Lemma 3.5. There are positive numbers h_{0} and C such that, for h < h_{0}, we have

∥(P_{1}^{h} − P_{1})f∥ ≤ Ch^{2}∥f∥, f ∈ L_{2}.

Proof. The proof is a modification of a standard argument on abstract spectral approximation that can be found for example in [2]. We let G = (A + αI)^{−1} and R^{h} be the operators defined in (3.2) and (3.4), and define G^{h} = R^{h}G, i.e.,

Ã(G^{h}f, χ) = (f, χ) ∀χ ∈ S^{h}, f ∈ L_{2},

or, equivalently, G^{h} = (A^{h} + αI)^{−1}P^{h}. It follows that G has positive eigenvalues μ_{i} = 1/(λ_{i} + α) and the same eigenfunctions ϕ_{i} as A. Similarly, G^{h} has eigenvalues μ_{i}^{h} = 1/(λ_{i}^{h} + α) and the same eigenfunctions ϕ_{i}^{h} as A^{h}. Applying the error bound (3.5) for v = Gf together with (3.3), we obtain

(3.21)  ∥(G^{h} − G)f∥ ≤ Ch^{2}∥f∥,  f ∈ L_{2}.

Let Γ be a positively oriented circle in the complex plane such that {μ_{i}}_{i=1}^{q}, but no other eigenvalues of G, lie inside Γ. Then for small h the eigenvalues {μ_{i}^{h}}_{i=1}^{q} of G^{h}, but no others, lie inside Γ and we have the representations

P_{1} = \frac{1}{2πi} ∫_{Γ}(zI − G)^{−1} dz, P_{1}^{h} = \frac{1}{2πi} ∫_{Γ}(zI − G^{h})^{−1} dz.

Hence

\begin{aligned}P_{1}^{h} − P_{1} & = \frac{1}{2πi} ∫_{Γ}((zI − G^{h})^{−1} − (zI − G)^{−1}) dz \\ & = \frac{1}{2πi} ∫_{Γ}(zI − G^{h})^{−1}(G^{h} − G)(zI − G)^{−1} dz, \end{aligned}

and the required bound follows by (3.21), since clearly ∥(zI −G^{h})^{−1}∥ and ∥(zI −G)^{−1}∥ are bounded uniformly for z ∈ Γ and small h. ⊔ ⊓

We can now state and prove our main result.

Theorem 1. There are positive numbers ρ_{0}, h_{1} and C such that, for any h < h_{1}, and for any τ, T with τ < T, the following holds: if u is a solution of (2.5) with u(t) ∈ B(ρ_{0}, ū) for t ∈ [τ, T], then there is a solution u^{h} of (3.6) such that

(3.22)  ∥u^{h}(t) − u(t)∥_{j} ≤ C(1 + (t − τ)^{− \frac{1}{2}} )h^{2−j},  t ∈ (τ, T], j = 0, 1.

Conversely, if u^{h} is a solution of (3.6) with u^{h}(t) ∈ B(ρ_{0}, ū) for t ∈ [τ, T], then there is a solution u of (2.5) such that (3.22) holds.

Proof. Let ρ and h_{0} be such that the conclusions of all the previous lemmas hold. If u(t) ∈ B(ρ_{0}, ū) for t ∈ [τ, T] and h < h_{1}, then we choose

(3.23)  v^{h} = P_{1}^{h}u(T), w^{h} = P_{2}^{h}u(τ),

in equation (3.11). In order to apply Lemma 3.3 we must check condition (3.16). By Lemma 3.2 we obtain

∥v^{h} − P_{1}^{h}ū^{h}∥_{1} + ∥w^{h} − P_{2}^{h}ū^{h}∥_{1} = ∥P_{1}^{h}(u(T) − ū^{h})∥_{1} + ∥P_{2}^{h}(u(τ) − ū^{h})∥_{1}

≤ C (∥u(T) − ū^{h}∥_{1} + ∥u(τ) − ū^{h}∥_{1})

≤ C (∥u(T) − ū∥_{1} + ∥u(τ) − ū∥_{1} + 2∥ū − ū^{h}∥_{1})

≤ 2C (ρ_{0} + h_{1}) ≤ \frac{ρ}{4M} ,

for ρ_{0} and h_{1} sufficiently small. Here we used the boundedness of P_{i}^{h} in H^{1}, which follows from the H^{1} boundedness of the projections P^{h} (by assumption) and P_{i}^{h}|_{S^{h}} (cf. the corresponding statement about P_{i} in the proof of Lemma 2.1). Now Lemma 3.3 shows the existence of u^{h}(t) ∈ B(ρ, ū) for t ∈ [τ, T], and Lemma 3.4 yields the error bound (3.22), since P_{1}^{h}(u^{h}(T) − u(T)) = 0, P_{2}^{h}(u^{h}(τ) − u(τ)) = 0.

Conversely, if u^{h}(t) ∈ B(ρ_{0}, ū) for t ∈ [τ, T], then we choose v = P_{1}u^{h}(T) and w = P_{2}u^{h}(τ) in equation (2.12). In order to apply Lemma 2.3 we must check condition (2.16). Using the H^{1} boundedness of P_{i}, we obtain

∥v − P_{1}ū∥_{1} + ∥w − P_{2}ū∥_{1} = ∥P_{1}(u^{h}(T) − ū)∥_{1} + ∥P_{2}(u^{h}(τ) − ū)∥_{1}

≤ C(∥u^{h}(T) − ū∥_{1} + ∥u^{h}(τ) − ū∥_{1})

≤ 2Cρ_{0} ≤ \frac{ρ}{2M} ,

for ρ_{0} sufficiently small. Lemma 2.3 now yields the existence of u(t) ∈ B(ρ, ū) for t ∈ [τ, T]. In order to apply Lemma 3.4 we use the identities

P_{1}^{h}P_{2} = P_{1}^{h}(I − P_{1}) = P_{1}^{h}(P_{1}^{h} − P_{1}),
P_{2}^{h}P_{1} = (P^{h} − P_{1}^{h})P_{1} = P^{h}(P_{1} − P_{1}^{h})P_{1},

and the boundedness of P_{i}^{h}, P^{h} and P_{i} in L_{2} together with Lemma 3.5, to get

∥P_{1}^{h}(u^{h}(T) − u(T))∥ = ∥P_{1}^{h}P_{2}(u^{h}(T) − u(T))∥ ≤ Ch^{2}∥u^{h}(T) − u(T)∥ ≤ Cρ_{0}h^{2},

and, similarly, ∥P_{2}^{h}(u^{h}(τ) − u(τ))∥ = ∥P_{2}^{h}P_{1}(u^{h}(τ) − u(τ))∥ ≤ Cρ_{0}h^{2}. ⊔ ⊓

Remark 1. If u(τ) ∈ H_{0}^{1}∩H^{2} (smooth initial data) in the first part of the theorem, then the factor (t − τ)^{− \frac{1}{2}} is not needed in (3.22), as can be shown by the appropriate modifications of Lemma 3.4 and Lemma 2.4.

Remark 2. The assumption that P^{h} is bounded with respect to the H^{1} norm was used only in the first part of the previous proof. We made this assumption only in order to streamline the formulation of our main result. In fact, if P^{h} is not bounded in H^{1}, then we may choose v^{h} = P_{1}^{h}R^{h}u(T), w^{h} = P_{2}^{h}R^{h}u(τ), instead of (3.23), and condition (3.16) is checked in the same way as before, using now the H^{1} boundedness of P_{i}^{h}R^{h}. (Both P_{i}^{h} and R^{h} are bounded in H^{1}: R^{h} is bounded because it is the orthogonal projection with respect to the inner product Ã(·, ·), see (3.4), which is equivalent to the usual inner product in H_{0}^{1}, and P_{i}^{h} are bounded because the eigenfunctions ϕ_{i}^{h} are orthogonal with respect to Ã(·, ·)). Lemma 3.4 is applicable, since

∥P_{1}^{h}(u^{h}(T) − u(T))∥ = ∥P_{1}^{h}(R^{h} − I)u(T)∥ ≤ ∥(R^{h} − I)u(T)∥ ≤ Ch^{2}∥u(T)∥_{2},

and similarly ∥P_{2}^{h}(u^{h}(τ) − u(τ))∥ ≤ Ch^{2}∥u(τ)∥_{2}, if u(T), u(τ) ∈ H_{0}^{1} ∩ H^{2}.

Remark 3. Our analysis clearly applies also in the case of a stable hyperbolic equilibrium, i.e., when q = 0 and P_{1}^{h} = 0, P_{2}^{h} = P^{h}. Lemma 3.4 then gives long-time error bounds similar to those of [8], [9] and [11].

Finally we show that the local stable and unstable manifolds of ū^{h} converge to their continuous counterparts. More precisely, the following result shows that M_{S}^{h}(ρ_{0}) (as defined in (3.17)) lies in an O(h)-neighborhood of M_{S}(ρ) (defined in (2.20)) with respect to the L_{2} norm for some radii ρ_{0} < ρ, and vice versa. This rate of convergence is the best that can be expected, because the stable manifold contains nonsmooth elements of H_{0}^{1}, cf. (3.5) with m = 1. The unstable manifold, on the other hand, is smooth and the rate of convergence is O(h^{2}).

It is convenient to express the result in terms of the semidistance δ(A, B) = sup_{a∈A} inf_{b∈B} ∥a−b∥ between two subsets A, B of L_{2}. Thus δ(A, B) < ϵ if and only if A lies in an ϵ-neighborhood of B in L_{2}. It would be desirable to have error bounds in the Hausdorff metric d(A, B) = max (sup_{a∈A} inf_{b∈B} ∥a − b∥, sup_{b∈B} inf_{a∈A} ∥a − b∥), but we are not able to achieve this due to a lack of symmetry in our argument.

Theorem 2. There are positive numbers ρ_{0}, ρ, h_{1} and C with ρ_{0} < ρ such that, for any h < h_{1}, the quantities δ(M_{S}(ρ_{0}), M_{S}^{h}(ρ)) and δ(M_{S}^{h}(ρ_{0}), M_{S}(ρ)) are bounded by Ch, and δ(M_{U}(ρ_{0}), M_{U}^{h}(ρ)) and δ(M_{U}^{h}(ρ_{0}), M_{U}(ρ)) are bounded by Ch^{2}.

Proof. Let ρ_{0}, ρ and h_{1} be as in the proof of Theorem 1. If u_{0} ∈ M_{S}(ρ_{0}), then there is a solution u(t) = u(t; u_{0}) ∈ B(ρ_{0}, ū) of equation (2.5) for t ∈ [0, ∞) such that u(0) = u_{0}, see (2.20). Set w^{h} = P_{2}^{h}u_{0}. As in the previous proof we check that ∥w^{h}−P_{2}^{h}ū^{h}∥_{1} ≤ ρ/4M, and a variant of Lemma 3.3 with τ = 0, T = ∞ shows that the discrete version of equation (2.15) has a solution u^{h}(t) ∈ B(ρ, ū) for t ∈ [0, ∞) with P_{2}^{h}u^{h}(0) = w^{h}. It follows that u_{0}^{h} ≡ u^{h}(0) ∈ M_{S}^{h}(ρ), see (3.17). A variant of Lemma 3.4 with τ = 0, T = ∞, cf. (3.19), now yields ∥u^{h}(t) − u(t)∥ ≤ Ch for t ∈ [0, ∞). In particular, using this result with t = 0, we may conclude that

δ(M_{S}(ρ_{0}), M_{S}^{h}(ρ)) = sup_{u_{0}∈M_{S}(ρ_{0})} inf_{u_{0}^{h}∈M_{S}^{h}(ρ)} ∥u_{0}^{h} − u_{0}∥ ≤ Ch.

Arguing as in the second part of Theorem 1 we obtain an analogous inequality with M_{S} and M_{S}^{h} interchanged, which proves the statements about the stable manifolds.

the behavior of finite element solutions 1019

The corresponding statements about the unstable manifolds are proved in a similar way, see (3.20). ⊔ ⊓

references

[1] F. Alouges and A. Debussche, On the qualitative behavior of the orbits of a parabolic partial differential equation and its discretization in the neighborhood of a hyperbolic fixed point, Numer. Funct. Anal. Optim. 12 (1991), pp. 253–269.

[2] I. Babuška and J. E. Osborn, Finite element-Galerkin approximation of the eigenvalues and eigenvectors of selfadjoint problems, Math. Comp. 52 (1989), pp. 275–297.

[3] W.-J. Beyn, On the numerical approximation of phase portraits near stationary points, SIAM J. Numer. Anal. 24 (1987), pp. 1095–1113.

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

[5] M. Crouzeix and V. Thomée, The stability in L_{p} and W_{p}^{1} of the L_{2}-projection onto finite element function spaces, Math. Comp. 48 (1987), pp. 521–532.

[6] J. Guckenheimer and P. Holmes, Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields, Springer-Verlag, 1983.

[7] D. Henry, Geometric Theory of Semilinear Parabolic Equations, Lecture Notes in Mathematics, Vol. 840, Springer-Verlag, Berlin, Heidelberg, New York, 1981.

[8] J. G. Heywood and R. Rannacher, Finite element approximation of the nonstationary Navier-Stokes problem II. Stability of solutions and error estimates uniform in time, SIAM J. Numer. Anal. 23 (1986), pp. 750–777.

[9] S. Larsson, The long-time behavior of finite element approximations of solutions to semilinear parabolic problems, SIAM J. Numer. Anal. 26 (1989), pp. 348–365.

[10] J.-M. Sanz-Serna and S. Larsson, Shadows, chaos, and saddles, Appl. Numer. Math. 13 (1993), pp. 181–190.

[11] J.-M. Sanz-Serna and A. M. Stuart, A note on uniform in time error estimates for approximations to reaction-diffusion equations, IMA J. Numer. Anal. 12 (1992), pp. 457–462.

[12] H. J. Stetter, Analysis of Discretization Methods for Ordinary Differential Equations, Springer Tracts in Natural Philosophy, Volume 23, Springer-Verlag, New York, 1973.