Convergence to equilibrium for the backward Euler scheme and applications

A result of Lojasiewicz [22, 23] asserts that if F is real analytic, then any bounded solution U(t) of (1) converges to a critical point of F as t tends to infinity. The proof relies on the so-called Lojasiewicz inequality for real analytic functions. If F is only assumed C and if d ≥ 2, the result may fail (see [1] and references therein). The Lojasiewicz inequality also allows to obtain optimal convergence rates to equilibrium. This result has many generalizations, in particular to infinite dimension [19, 27] and to gradient-like flows [15, 24]), and there is a growing literature extending it to other situations (see for instance [6, 14]). The aim of this paper is to show that these convergence results also apply to some stable discretizations of gradient flows, and first of all, to the backward Euler scheme. We recall that for the model example (1), the backward Euler scheme reads: U − U ∆t = −∇F (U), n ≥ 0, (2)

1. Introduction. In this paper, we are concerned with the asymptotic behaviour, as time goes to infinity, of the solution of the backward Euler scheme applied to gradient flows. As a model example, consider the following gradient flow: where U = (u 1 , . . . , u d ) t , F ∈ C 1,1 loc (R d , R) and A result of Lojasiewicz [22,23] asserts that if F is real analytic, then any bounded solution U (t) of (1) converges to a critical point of F as t tends to infinity. The proof relies on the so-called Lojasiewicz inequality for real analytic functions. If F is only assumed C ∞ and if d ≥ 2, the result may fail (see [1] and references therein). The Lojasiewicz inequality also allows to obtain optimal convergence rates to equilibrium. This result has many generalizations, in particular to infinite dimension [19,27] and to gradient-like flows [15,24]), and there is a growing literature extending it to other situations (see for instance [6,14]). The aim of this paper is to show that these convergence results also apply to some stable discretizations of gradient flows, and first of all, to the backward Euler scheme. We recall that for the model example (1), the backward Euler scheme reads: where ∆t > 0 is the time step. The implicit Euler scheme is widely used in the numerical or theoretical treatment of flows governed by parabolic semi-linear PDE's, because of its stability properties (see for instance [3,17]). The question of convergence to equilibrium for the backward Euler scheme is therefore a natural question, and to our knowledge, it seems to be still open (see however Remark 1.1 below). Notice that a related question has been answered for some explicit descent methods by Absil, Mahony and Andrews in a context of optimization [1]. These authors, however, assume a so-called "primary descent condition" which, in the case of the backward Euler scheme, is impossible to check a priori. In [11], Gajewski and Griepentrog also prove convergence to equilibrium for a descent method (distinct from the backward Euler scheme) in a Hilbert setting, using the Lojasiewicz inequality.
The starting point of our paper, Theorem 2.4 in Section 2, shows that if F is real analytic and coercive, then any solution (U n ) n≥0 of (2) converges to a critical point of F as n goes to infinity. We also obtain, in Proposition 2.5, optimal convergence rates, similarly to the continuous case. As a corollary of Theorem 2.4 and under reasonable assumptions, the backward Euler approximation of a gradient flow is shown to converge uniformly in time towards the exact solution, as the time step tends to 0. The rest of the paper is organized as follows: in Section 3, we show how the preceding results apply to a time and space discretization of the Cahn-Hilliard equation. In Section 4, we extend Theorem 2.4 to the θ-scheme, and in the last section, we extend Theorem 2.4 in infinite dimension to the case of the semilinear heat equation.
Remark 1.1. After this paper was complete, we learned that Bolte et al. recently proved convergence to equilibrium for the backward Euler scheme with variable step size in a Hilbert setting: see [4,Theorem 24]. Their proof is also based on a Lojasiewicz inequality.
2. Convergence to equilibrium for the backward Euler scheme. In this section, we will consider the backward Euler scheme as a descent method. In this regard, we assume the following coercivity condition on F : This assumption, together with the continuity of F , implies that F is bounded from below. For convenience, we define the implicit Euler scheme as follows: let U 0 ∈ R d and for all n ≥ 0, Here and in the following, V denotes the euclidean norm of a vector V ∈ R d . Assumption (3) guarantees that problem (4) has at least one solution for all n. If F ∈ C 1 (R d ), then clearly, any solution U n+1 of (4) satisfies (2) for all n. A solution of (4) may not be unique, but in many applications, uniqueness of solutions for (2) is guaranteed for ∆t > 0 small enough, assuming some additional assumptions on F (see Sections 3 and 4, for instance); in such cases, (4) and (2) are equivalent. Convergence to equilibrium for (4) will be a consequence of the Lojasiewicz inequality, which is defined as follows.
Definition 2.1. We say that F ∈ C 1 (R d , R) satisfies the Lojasiewicz inequality at a point U ∈ R d if there exist ν ∈ (0, 1/2], γ > 0 and σ > 0 such that In some (few) cases, the exponent ν is known explicitly, and it is called the Lojasiewicz exponent at U . It allows the computation of convergence rates. For instance, F : R d → R defined by F (U ) = U p , p ≥ 2, satisfies (5) atŪ = 0 with ν = 1/p. Another example is the following result, which is a consequence of the inverse mapping Theorem and of a Taylor's expansion (see for instance [16]): ,j≤d is non singular, then F satisfies the Lojasiewicz inequality at U with exponent ν = 1/2. However, most of the time, the Lojasiewicz inequality is a consequence of analyticity, and the exponent ν is not known. This is the fundamental result: [23]). If F : R d → R is real analytic in the neighborhood of a point U, then F satisfies the Lojasiewicz inequality at U .
We denote S = {V ∈ R d : ∇F (V ) = 0} the set of critical points of F . Our main result is: (3) and that the Lojasiewicz inequality (5) holds for every U ∈ S. Let (U n ) n≥0 be a sequence defined by (4). Then there exists U ∞ ∈ S such that U n → U ∞ as n → +∞.
Proof. By (4), Thus, the sequence (F (U n )) is nonincreasing, and since F is bounded from below, F (U n ) tends to some F ⋆ ∈ R. Without loss of generality, we assume F ⋆ = 0. By (3), (U n ) is bounded, so there exist U ∞ ∈ R d and a subsequence (U n k ) such that U n k → U ∞ as k → +∞. By (6), U n+1 − U n → 0, and letting n = n k tend to ∞ in (2), we obtain ∇F (U ∞ ) = 0, so U ∞ ∈ S. Since F satisfies the Lojasiewicz inequality at U ∞ , there exist ν ∈ (0, 1/2] and σ, γ > 0 such that Let n be such that U n+1 − U ∞ < σ. We consider two cases: Thus, In both cases, for all n such that U n+1 − U ∞ < σ, we have Now, letF > 0 be small enough so that We choosen large enough such that Un − U ∞ < σ/3, and F (Un) ≤F . Let N ≥n be the largest integer (including +∞) such that U n − U ∞ < 2σ/3 for all n ≤ n ≤ N . Assume by contradiction that N is finite. We deduce from (6) that So we may apply (8) to everyn ≤ n ≤ N and since (F (U n )) is nonincreasing, we obtain , and this contradicts the definition of N . So N = +∞ and the whole sequence converges.
We are also able to estimate the convergence rates to equilibrium: they are analogous to the continuous case (see for instance [16]).
Similarly, the Euler scheme, in spherical coordinates, reduces to For p = 2, we have (1 + 2∆t)r n+1 = r n , so that r n = r 0 (1 + 2∆t) −n , and we see that the geometric convergence proved in Proposition 2.5 is optimal. If p > 2, the solution of (10) with initial condition r(0) = r 0 > 0 is On the other hand, it is easily seen that if r n ≥ r(n∆t), then the solution of (11) satisfies r n+1 ≥ r((n + 1)∆t). By induction, the sequence uniquely defined by (11) and the initial condition r 0 = r(0) satisfies r n ≥ r(n∆t), so that for some constant C = C(r 0 , p, ∆t) > 0. Thus the exponent found in Theorem 2.5 is optimal (notice that 1/(p − 2) = ν/(1 − 2ν)).
Proof of Proposition 2.5. We argue as in the proof of Theorem 2.4. By (6), the sequence (F (U n )) is nonincreasing and converges to F (U ∞ ) = 0. There exists n ≥ 0 such that, for n ≥ñ, U n − U ∞ < σ, where σ is defined by the Lojasiewicz inequality (7). If F (Uñ) = 0 for someñ, then by (6), U n = Uñ = U ∞ for all n ≥ñ, and the result is obvious. From now on, we assume that F (U n ) > 0 for all n. Then, applying estimate (8) with n, n + 1, · · · and summing, we get Next, we estimate the rate of decay of the sequence (F (U n )). Let us define the function G : (0, +∞) → (0, +∞) by The sequence (G(F (U n ))) is nondecreasing and tends to +∞. Moreover, for n ≥ñ, In both cases, we have for n ≥ñ, Writing this relation forñ,ñ + 1, · · · , n − 1 and summing, we obtain Thus, and recalling (12), the proof is complete in the case ν = 1/2.
In [24, Theorem 2.2], Miranville and Rougirel noticed that the Lojasiewicz inequality can sometimes be used to show the continuity of the limiting state with respect to some parameters. Here, we adapt their remark to show a result concerning convergence of the discrete solution to the continuous solution as ∆t → 0. For this purpose, for ∆t > 0, we denote (U n ∆t ) a sequence generated by (4), and we define the continuous and piecewise affine function U ∆t : [0, +∞) → R d associated to (U n ∆t ) by We have: This assumption cannot be removed in general: if U ∞ ∈ S has a non trivial unstable set, the conclusion of Proposition 2.6 is hardly true.
In the next section, we give an application of Theorem 2.4 to a PDE. The following remark will be useful.
Remark 2.7. In Theorem 2.4 and Proposition 2.6, we can replace the usual euclidean norm (which arises through (4)) by any euclidean norm. More precisely, let A be a symmetric positive definite matrix and let V 1 , V 2 A = (AV 1 , V 2 ) (where (·, ·) is the usual scalar product on R d ) be the scalar product associated to A, with norm · A . In Theorem 2.4 and Proposition 2.6, we can replace assumption (4) by the following definition: The proofs are similar, replacing · by · A . Notice that in this case, U n+1 satisfies This is the implicit Euler scheme for the flow which is of course the gradient flow of F for the scalar product ·, · A .
3. Application to a time and space discretization of the Cahn-Hilliard equation. Let us consider the Cahn-Hilliard equation where Ω is a bounded domain of R d (1 ≤ d ≤ 3) with Lipschitz continuous boundary, γ > 0, and f is a polynomial of degree 2p + 1 with positive leading coefficient, p ≥ 1 if d = 1 or d = 2, and 1 ≤ p ≤ 2 if d = 3. Equation (17) is suplemented with an initial condition u 0 and Neumann or periodic boundary conditions (if periodic boundary conditions are considered, then Ω is a parallelepiped). The Cahn-Hilliard equation is mass conservative, and it can be seen as a gradient flow for the H −1 scalar product, for the energy where F (s) is an antiderivative of f (s). Convergence to equilibrium for solutions of the Cahn-Hilliard equation, as t → +∞, has been shown in [26]. This is all the more interesting since stationary solutions of (17) are neither unique nor isolated in general (see for instance [18,28]). We show below that this convergence still holds for a time and space discretization of (17), where the time discretization is the backward Euler scheme, and the space discretization is a finite element method. The variational formulation of (17), as introduced by Elliott et al [10], reads as follows. Let V = H 1 (Ω) (for Neumann boundary conditions) or V = H 1 per (Ω) (for periodic boundary conditions), and for a given u 0 ∈ V , find u, w : where (·, ·) denotes the L 2 (Ω)-scalar product. The L 2 (Ω)-norm is denoted · 0 . For the space semi-discretization, we replace V by V h , a P k finite element approximation approximation of V (see for instance [7,10] The total discretization by the backward Euler scheme was considered by Elliott [9]. It reads: let u h,0 ∈ V h , and for n ≥ 1, find (u h,n , w h, Since any constant function on Ω h belongs to V h , the mass is preserved, i.e. (u h,n , 1) = (u h,0 , 1), ∀n ≥ 0.
A steady state for this problem is a solution (ū h ,w h ) ∈ V h × R such that The assumption on the nonlinearity f implies that there exists c f ≥ 0 such that Elliott proved existence and uniqueness of a solution (u h,n , w h,n ) n≥0 for this scheme, for the nonlinearity f (s) = s 3 − β 2 s [9]. He also proved dynamical stability, which implies convergence to equilibrium under the assumption that the steady states are isolated. However, such an assumption is impossible to check in general. Here, we prove convergence to equilibrium without any knowledge on the structure of steady states.
Proof. The proof of uniqueness is the same as in [9] (it requires the smallness assumption made on ∆t). For the proof of existence and convergence to equilibrium, we apply Theorem 2.4 and Remark 2.7: we show below that the scheme (18)- (19) can be seen as the backward Euler scheme applied to a gradient flow which is real analytic.
Let (e j ) 1≤j≤N be an orthonormal basis of V h for the L 2 (Ω h )-scalar product, with e 1 ≡ cste. In matrix representation, (18)- (19) reads: where U n , (resp. W n ) is the vector of the coordinates of u h,n (resp. w h,n ), A is the matrix of the Laplacian, I is the identity matrix, and The matrix A is not invertible because of the constants. However, since the mass is preserved, we can eliminate the constants by writing, for any vector X ∈ R N , X = (x 1 ,Ẋ), withẊ ∈ R N −1 . Now,Ȧ is invertible, because it is the matrix of the scalar product (∇·, ∇·) on V h /R. The scheme in matrix form becomes This is equivalent toẆ = −Ȧ −1 (U n −U n−1 )/∆t, We recognize the backward Euler scheme for the floẇ This is a gradient flow for the function where, forU ∈ R N −1 , u h is defined by u h = u 0 1 e 1 + N i=2 u i e i (u 0 1 is the first component of u h,0 ). The assumptions on f imply that F is bounded from below, so E h (U ) → +∞ as ∇u h 0 → +∞. SinceU → ∇u h 0 is a norm on R N −1 , this means that E h satisfies assumption (3).
The assumption on f implies that E h is a polynomial of the variables (u 2 , . . . , u N ). The proof is complete.
A similar result holds for the discrete version of the viscous Cahn-Hilliard equation considered in [2]. Notice that in both cases, the convergence rates from Section 2 are valid (for some ν ∈ (0, 1/2] which is not known explicitly here).
4. Convergence to equilibrium for the θ-scheme. In this section, we extend the preceding results on the backward Euler scheme to a more general family of schemes, namely the θ-scheme. Recall that for a fixed θ ∈ [0, 1], the θ-scheme associated to the gradient flow (1) reads: For θ = 1 we recover the backward Euler scheme (and for θ = 0, the forward Euler scheme). The θ-scheme has order one, except when θ = 1/2, in which case it has order two. The 1/2-scheme is also known as the Crank-Nicholson scheme.
In [17], Humphries and Stuart proved that when θ ∈ [1/2, 1] (i.e. when the implicit part is dominant) and if F satisfies an additional natural assumption, then the scheme (21) is dynamically stable. For this purpose, they assume the following one-sided Lipschitz condition on ∇F : there exists c ≥ 0 such that Assumption (22) is equivalent to saying that ∇F + cI is a monotone operator [5], or equivalently, that the function (6) and (22), then for every U n , equation (21) has at least one solution U n+1 , obtained as a minimizer of the function which is bounded from below and coercive. Moreover, if ∆t < 1/(θc), then this solution is unique, by strict convexity.
Following [17], we now define, for ∆t > 0, The stability result of Humphries and Stuart reads: (3) and (22). Then, for every 0 < ∆t ≤ 1/c, and for every sequence (U n ) defined by (21), the sequence (F ∆t (U n )) satisfies In particular, the sequence (U n ) is bounded, and this can be used to prove that the ω-limit set of any initial condition U 0 is a subset of S. By showing a modified Lojasiewicz inequality for F ∆t , we are able to prove convergence to equilibrium. In the remainder of this section, the parameter θ is fixed in [1/2, 1]. (22) and that the Lojasiewicz inequality (5) holds for everyŪ ∈ S. Let U 0 ∈ R d . Then there exists ∆t ⋆ = ∆t ⋆ (U 0 ) ∈ (0, 1/c) such that for every ∆t ∈ (0, ∆t ⋆ ], if (U n ) n≥0 is the unique sequence defined by (21), then U n → U ∞ as n → +∞, for some U ∞ ∈ S.
Theorem 4.2 can be applied, for instance, to the space and time discretization of the viscous diffusion equation u t − ν∆u t = ∆f (u) (ν > 0, f (u) = u 3 − u) used in [25]. It could also be applied, as in Section 3, to the Cahn-Hilliard equation or to the viscous Cahn-Hilliard equation, by replacing the backward Euler scheme by the θ-scheme.
Remark 4.4. 1. It is possible to obtain convergence rates for the θ-scheme, by the same proof as in Section 2. 2. It is also possible to adapt Proposition 2.6 to the θ-scheme.

5.
A model problem in infinite dimension. Our aim, in this section, is to see that Theorem 2.4 can also be generalized, to some extent, in infinite dimension. For this purpose, we consider the following model problem: where Ω is a bounded domain of R d with smooth boundary. Equation (26) is completed with Neumann, Dirichlet, or periodic boundary conditions, and with an initial condition. In [19,27], it was proved that a solution of (26) converges to equilibrium, provided that f is real analytic and that the orbit of u is precompact in an adequate function space. This problem has also been considered by many authors (see the references in [19]). Here, we prove that under some assumptions which remind the ones made in Theorem 4.2, any solution of the backward Euler scheme applied to (26) converges to equilibrium. In order to be slightly more general, we first consider an abstract version of (26).

A general result.
Let V denote a closed subspace of H 1 (Ω) which is dense in L 2 (Ω); in particular, V is a Hilbert space for the H 1 (Ω)-scalar product. Identifying L 2 (Ω) with its dual, we have the continuous inclusions V ֒→ L 2 (Ω) ֒→ V ′ . Let a : V × V → R be a continuous bilinear form on V which is symmetric and coercive, i.e.
We consider a Caratheodory function f : Ω × R → R assumed to be of class C 1 in the variable s ∈ R and such that and we denote F (x, s) = s 0 f (x, σ)dσ. We assume the following growth condition: We assume in addition that (compare with (22)) Under the previous assumptions, by monotonicity, for every initial data u 0 ∈ L 2 (Ω), the abstract evolution equation has a unique weak solution u ∈ C 0 ([0, +∞); L 2 (Ω)) such that u(0) = u 0 (see [21] or [5, Example 2.3.7 and Theorem 3.17]). Equation (32) can be seen as a gradient flow, for the L 2 (Ω) scalar product, of the energy functional By the assumptions made above, E is of class C 2 on V , with Now, consider the backward Euler scheme applied to (32): let u 0 ∈ L 2 (Ω), and for n ≥ 0, let u n+1 be such that where δt > 0 is the time step. It is easy to see that for δt > 0 small enough, problem (34) has a unique solution u n+1 for all n. Indeed, by (27) and (31), where · 0 denotes the L 2 (Ω)-norm. Thus, if 1/δt ≥ c f , then for every u n ∈ L 2 (Ω), is α-convex on V : it has a global minimizer u n+1 in V , which is the unique solution of (34).
We will use the following definition, which generalizes Definition 2.1 in this context.
Definition 5.1. We say that E : V → R defined by (33) satisfies the Lojasiewicz inequality in V ′ at a pointū ∈ V if there exist ν ∈ (0, 1/2] and σ > 0, We also define the set of critical points of E We have: Theorem 5.2. Let δt > 0 such that 1/δt ≥ c f . Assume that (27)-(31) hold and that E satisfies the Lojasiewicz inequality at everyū ∈ S E . Let (u n ) be a sequence defined by (34).
Remark 5.3. Proposition 2.6 could also be generalized in a similar way. It is also possible to obtain convergence rates in L 2 (Ω)-norm, by the method used in section 2 and methods used in the continuous case (see for instance [16]).
We have Theorem 5.4. Assume that f : R → R is real analytic and satisfies (40)-(42). Let δt ≤ 1/c f , and let (u n ) be a sequence defined by the backward Euler scheme (44).
An important example of a function f which satisfies the required assumptions is f (s) = s 3 − s, when 1 ≤ d ≤ 3. In this case, equation (26) is known as the Allen-Cahn (or Ginzburg-Landau) equation.
Proof. We apply Theorem 5.2. We only need to prove that for everyū ∈ S E , the functional E satisfies the Lojasiewicz inequality, and that ∪ n≥1 {u n } is precompact in H 1 0 (Ω). Letū ∈ S E . By elliptic regularity [12],ū ∈ C ∞ (Ω). Corollary 4.5 in [6] shows that E satisfies the Lojasiewicz inequality atū. This proves the first point.
For the second point, we first notice that E(u 1 ) < +∞, because u 1 is a minimizer in V of the functional v → v − u 0 2 0 2 + E(v).
In the case d = 2, we directly obtain from the Sobolev imbedding theorem that (f (u n+1 )) is bounded in any L q (Ω), q < ∞ and we conclude similarly.
In the case d = 1, H 1 (Ω) is compactly embedded in C(Ω) so (f (u n+1 )) is uniformly bounded in L ∞ (Ω). We conclude as above, and this completes the proof.