Dominique Chapelle

Asven Gariah

Jacques Sainte-Marie

Galerkin approximation with Proper Orthogonal Decomposition: new error estimates and illustrative examples

We propose a numerical analysis of Proper Orthogonal Decomposition (POD) model reductions in which a priori error estimates are expressed in terms of the projection errors that are controlled in the construction of POD bases. These error estimates are derived for generic parabolic evolution PDEs, including with non-linear Lipschitz right-hand sides, and for wave-like equations. A specific projection continuity norm appears in the estimates andwhereas a general uniform continuity bound seems out of reach -we prove that such a bound holds in a variety of Galerkin bases choices. Furthermore, we directly numerically assess this bound -and the effectiveness of the POD approach altogether -for test problems of the type considered in the numerical analysis, and also for more complex equations. Namely, the numerical assessment includes a parabolic equation with super-linear reaction terms, inspired from the FitzHugh-Nagumo electrophysiology model, and a 3D biomechanical heart model. This shows that the effectiveness established for the simpler models is also achieved in the reduced-order simulation of these highly complex systems.

Introduction

In general, the simulation of partial differential equations resorts to discretization techniques such as finite differences, finite elements, or finite volumes. This typically results in discrete systems of large dimensions, hence the solution process can be rather costly, especially in situations when many computational iterations are required, as often occurs in design, control applications and inverse modeling.

In order to obtain reduced-order models, two main approaches are generally used. The first one consists in analyzing the dynamics operator of the system considered and retaining only the "most significant parts". Modal Analysis (linear or non-linear normal modes), but also the Moment Matching Method [START_REF] Daniel | A multiparameter momentmatching model-reduction approach for generating geometrically parameterized interconnect performance models[END_REF][START_REF] Astolfi | Model reduction by moment matching for linear and nonlinear systems[END_REF] and Balanced Truncation [START_REF] Gugercin | A survey of model reduction by balanced truncation and some new results[END_REF][START_REF] Stykel | Balanced truncation model reduction for semidiscretized Stokes equation[END_REF] belong to this first family. Unfortunately, for complex and large systems, these tools can be difficult to use in practice, since e.g. the eigenmodes are costly to obtain.

The second strategy is more data-oriented in the sense that it mainly uses snapshots of the system to perform its reduction. The Reduced Basis [START_REF] Maday | A priori convergence theory for reduced-basis approximations of single-parameter elliptic partial differential equations[END_REF][START_REF] Rovas | Reduced-basis output bound methods for parabolic problems[END_REF][START_REF] Prud'homme | A mathematical and computational framework for reliable real-time solution of parametrized partial differential equations[END_REF][START_REF] Veroy | Reduced-basis approximation of the viscous Burgers equation: rigorous a posteriori error bounds[END_REF][START_REF] Rozza | Reduced basis approximation and a posteriori error estimation for affinely parametrized elliptic coercive partial differential equations: application to transport and continuum mechanics[END_REF] and the Proper Orthogonal Decomposition (POD) [START_REF] Kunisch | Proper orthogonal decomposition for optimality systems[END_REF][START_REF] Kahlbacher | Galerkin proper orthogonal decomposition methods for parameter dependent elliptic systems[END_REF][START_REF] Kunisch | Galerkin proper orthogonal decomposition methods for parabolic problems[END_REF][START_REF] Hinze | Proper orthogonal decomposition surrogate models for nonlinear dynamical systems: Error estimates and suboptimal control[END_REF][START_REF] Willcox | Balanced model reduction via the proper orthogonal decomposition[END_REF] are two techniques belonging to this second family. This second approach consists in projecting the system onto subspaces of reduced sizes, albeit containing the major part of the expected dynamical solution. The aim is to obtain low-dimensional systems capturing the essence of the phenomena of interest.

Proper Orthogonal Decomposition, also known as Karhunen-Loève decomposition or principal component analysis, is a method initially introduced for analyzing multidimensional data. This method essentially provides an orthonormal basis for representing the given data in an optimal manner with respect to a quadratic criterion. The work in [START_REF] Kosambi | Statistics in function space[END_REF] has been pioneering in the development of the POD technique. In fluid mechanics, POD has been successfully used to access the coherent structures in turbulent flows [START_REF] Holmes | Turbulence, Coherent Structures, Dynamical Systems and Symmetry[END_REF], and it is now widely used in engineering in general.

Despite its relative simplicity of development and use, the POD technique has some limitations, since in particular it does not guarantee stability e.g. when parametric variations are considered [START_REF] Amsallem | Interpolation method for adapting reduced-order models and application to aeroelasticity[END_REF]. Moreover, existing error estimates are expressed with respect to quantities which are not controlled in the construction of the POD basis [START_REF] Hinze | Error estimates for abstract linear-quadratic optimal control problems using proper orthogonal decomposition[END_REF]. This latter important issue is our primary concern here.

In this article, we propose new error estimates for the POD-based Galerkin approximation of the solutions of some classical and widely used PDE systems. First, we briefly recall the foundations of the POD decomposition. Then we derive the estimates for linear and non-linear parabolic equations, and also for linear hyperbolic systems. Finally, the theoretical results are confronted with numerical tests in various situations including a complex 3D biomechanical heart model.

Classical principles of POD reduction

In this section we briefly summarize the general principles and construction rules for POD reductions.

Construction of the POD basis

Let V be a separable Hilbert space with scalar product •, • and norm • . Let z(t), t ∈ [0, T ], be a function with regularity z ∈ L 2 (0, T ; V ).

We introduce Cov : V → V the covariance operator defined by Cov ϕ = T 0 z(t), ϕ z(t) dt.

Performing the POD (time-continuous here) of rank l of z over [0, T ] consists in finding the orthogonal projector π l of rank l solution of min πl z(t) -πl z(t) L 2 (0,T ;V ) .

The integer l is called the POD rank. The time-discrete POD consists in solving the problem

min πl N i z(t i ) -πl z(t i ) 2 V , (2)
and for a detailed presentation of the solution of [START_REF] Astolfi | Model reduction by moment matching for linear and nonlinear systems[END_REF], see e.g. [START_REF] Kunisch | Proper orthogonal decomposition for optimality systems[END_REF][START_REF] Kunisch | Galerkin proper orthogonal decomposition methods for a general equation in fluid dynamics[END_REF] and references therein. In the continuous case, i.e. for the problem [START_REF] Amsallem | Interpolation method for adapting reduced-order models and application to aeroelasticity[END_REF], the results and their proofs are, to some extent, similar but some aspects need to be specified. The four following propositions are the cornerstones of the solution of (1). The complete proofs of these propositions rely on straightforward adaptations of results contained in [START_REF] Kunisch | Galerkin proper orthogonal decomposition methods for a general equation in fluid dynamics[END_REF].

Proposition 1. There exists a unique sequence (λ i ) i∈I with I finite or countable, such that

λ i > 0, ∀i ∈ I, λ 1 ≥ λ 2 ≥ • • • ≥ λ N , if I finite (I = {1, 2, . . . , N }), λ 1 ≥ λ 2 ≥ • • • ≥ λ i ≥ . . . , λ i -→ i→∞ 0, otherwise (I = N * ),
and an orthonormal sequence

(ϕ i ) i∈I of V satisfying Cov ϕ i = λ i ϕ i , ∀i ∈ I,
such that (ϕ i ) i∈I is total in the orthogonal complement of the kernel of Cov, i.e.

V = Ker Cov ⊥ ⊕ Span{ϕ i , i ∈ I}. (3)
To understand (3) it is helpful to characterize the kernel of Cov with respect to z. Then we have the classical result.

Proposition 3. For all 1 ≤ l ≤ Card I, a solution π l of Problem ( 1) is determined by

V l = Im π l = Span(ϕ 1 , . . . , ϕ l ).
Moreover,

z -π l z L 2 (0,T ;V ) = Card I i=l+1 λ i 1/2 . ( 4
)
The POD reduction is interesting when the sequence of eigenvalues λ i tends rapidly to zero. Indeed, for small values of l we then have the approximation

z(t) ≈ l i=1 z(t), ϕ i ϕ i .
In practice, the dimension of V can be very large and it is costly to use the operator Cov. So it is convenient to introduce Cov the L 2 (0, T ) → L 2 (0, T ) operator defined by Thus there exists an orthonormal sequence (v i ) i∈I of L 2 (0, T ) eigenvectors of Cov, in finite number for each non-zero eigenvalue,

Cov v (s) = T 0 z(t), z(s) v(t) dt = T 0 v(t)z(t) dt, z(s) ,
Cov v i = λ i v i , ∀i ∈ I, such that (v i ) i∈I is total in the orthogonal complement of the kernel of Cov, i.e. L 2 (0, T ) = Ker Cov ⊥ ⊕ Span{v i , i ∈ I}.
We define ∀i ∈ I the V element

ϕ i = 1 √ λ i T 0 v i (t)z(t) dt.
Then, (ϕ i ) i∈I is an orthonormal sequence of eigenvectors of Cov, with the same sequence of corresponding eigenvalues Cov ϕ i = λ i ϕ i , ∀i ∈ I.

Reduced-order modeling

Considering z(x, t) the solution of a PDE problem, the POD-based reduced order modeling, or more simply POD reduction, consists in building a spatial Galerkin approximation z l (x, t) of z(x, t) in the POD space V l = Span(ϕ 1 , . . . , ϕ l ). Then the key point is to be able to control the reduction error, namely z -z l L 2 (0,T ;V ) . We tackle this problem in the following section.

New estimates for the POD reduction error

In this section, our objective is to derive POD-reduction error estimates bounded by approximation terms which can be conveniently controlled in the construction of the POD basis.

For the sake of generality and homogeneity with the existing literature, we introduce the classical abstract mathematical framework. Nevertheless, to fix the ideas the reader can keep in mind that in the examples considered, the abstract spaces H and V will typically correspond to L 2 (Ω) and H 1 0 (Ω), respectively. Let (V, ((•, •)), • ) and (H, (•, •), | • |) be two separable Hilbert spaces with continuous and dense embedding V → H, i.e.

|v| ≤ C Ω v , ∀v ∈ V.

We choose H as the pivot space -namely, we perform the identification of H with its dual space H -and then

V → H → V .
Let a be a symmetric bilinear form on V , continuous, and coercive, namely,

a(v, w) ≤ C a v w , ∀v, w ∈ V, a(v, v) ≥ c a v 2 , ∀v ∈ V.
Then a also defines a scalar product on V and we denote by • a the associated norm.

We point out that in our estimations we use C to denote a generic positive constant, independent of all discretization parameters, and that may take different values at various occurrences, including in the same equation.

Galerkin estimates for linear parabolic problems with H-orthogonal projectors

We formally introduce the abstract parabolic equation

d dt (u(t), v) + a(u(t), v) = (f (t), v), ∀v ∈ V, (5)
u(0) = u 0 . (6)
Equation ( 5) is to be understood in the sense of distributions in time. We then have the following existence and uniqueness result [8, XVIII, §3.2, Th. 1, §3.3, Th. 2].

Proposition 5. Assume f ∈ L 2 (0, T ; H) and u 0 ∈ H. Then there exists a unique solution u of Eqs. ( 5)-( 6) such that

u ∈ L 2 (0, T ; V ) ∩ C([0, T ]; H), du dt ∈ L 2 (0, T ; V ).
Considering now a finite-dimensional subspace V l , we formally introduce the spatial Galerkin approximation u l of u

u l (t) ∈ V l , (7) d dt (u l (t), v l ) + a(u l (t), v l ) = (f (t), v l ), ∀v l ∈ V l , (8)
u l (0) = u l 0 . (9)
Since V l is a finite-dimensional space, it is easy to prove the following result [8, XVIII, §3.3.1, Lemma 1].

Proposition 6. Assume f ∈ L 2 (0, T ; H) and u l 0 ∈ V l . Then there exists a unique solution u l of Eqs. ( 8)-( 9) such that

u l ∈ C([0, T ]; V l ), du l dt ∈ L 2 (0, T ; V l ).
Note that we have more regularity here than for the continuous solution only because we are considering a finite dimensional problem, and of course the corresponding estimates are not uniform with respect to the discretization.

Finally, let π l H and π l V respectively denote the H-orthogonal and V -orthogonal projectors of V onto the reduction space

V l . For all v ∈ V |v -π l H v| = inf v l ∈V l |v -v l |, v -π l V v = inf v l ∈V l v -v l .
With a view to estimating the reduction error u -u l L 2 (0,T ;V ) , classical error estimates are of the form [START_REF] Dautray | Mathematical Analysis and Numerical Methods for Science and Technology[END_REF]

u -u l L 2 (0,T ;V ) ≤ C u -π l V u L 2 (0,T ;V ) + ∂ ∂t (u -π l V u) L 2 (0,T ;V ) + |u l 0 -π l V u 0 | . ( 10
)
However, the POD criterion (1) does not provide a direct control on the time-derivative term in the right-hand side, and our objective is to circumvent this difficulty. To that end we use the H-projection error, still in the same L 2 (0, T ; V )-norm, i.e.

u -π l H u L 2 (0,T ;V ) .

and the following result holds.

Proposition 7. For all T > 0,

u -u l L 2 (0,T ;V ) ≤ C(|π l H u 0 -u l 0 | + u -π l H u L 2 (0,T ;V ) ).
Proof. We split u -u l into two parts u -u l = p l + q l , where p l = u -π l H u and q l = π l H u -u l . Since q l ∈ V l , and using the definition of u l , q l satisfies the variational equation

d dt (q l (t), v l ) + a(q l (t), v l ) = -(f (t), v l ) + d dt (π l H u(t), v l ) + a(π l H u(t), v l ), ∀v l ∈ V l .
The projection π l H satisfies (π l H u(t), v l ) = (u(t), v l ), so that using the definition of u we get d dt q l (t), v l + a(q l (t), v l ) = -a(p l (t), v l ).

Taking v l = q l (t), we then obtain the energy estimate

1 2 d dt {|q l | 2 }(t) + q l (t) 2 a = -a(p l (t), q l (t)),
which we now integrate on [0, T ] to obtain,

T 0 q l (t) 2 a dt ≤ 1 2 |q l (0)| 2 + T 0 a(p l (t), q l (t)) dt ,
where we dropped the term in |q l (T )| 2 in the left-hand side. Hence, by Young's inequality,

T 0 q l (t) 2 a dt ≤ |q l (0)| 2 + T 0 p l (t) 2 a dt.
This directly entails

q l L 2 (0,T ;V ) ≤ C |q l (0)| 2 + p l L 2 (0,T ;V ) ,
using the properties of the scalar product a, and the triangle inequality

u(t) -u l (t) ≤ p l (t) + q l (t)
concludes the proof.

Therefore, via the introduction of π l H we avoid the time derivative appearing in the righthand side of the more standard estimate [START_REF] Flett | Differential Analysis[END_REF]. However, we now need to deal with the approximation term u -π l H u L 2 (0,T ;V ) , which is the topic of the next section.

Galerkin estimates for linear parabolic problems with V -orthogonal projectors

Note first that, since V l is finite-dimensional, we have an inverse inequality of the form

∃α l > 0, ∀v l ∈ V l , v l ≤ α l |v l |.
Hence, π l H is continuous as an endomorphism of V , as can be seen by directly writing

π l H v ≤ α l |π l H v| ≤ α l |v| ≤ Cα l v .
However, as inverse inequality constants blow up when the dimension of the V l subspace increases, we will obtain some better insight by using the V -projection as follows

π l H v ≤ π l H v -π l V v + π l V v ≤ α l |π l H v -π l V v| + v = α l |π l H (v -π l V v)| + v ≤ α l |v -π l V v| + v .
Denoting by L(V ) the space of V -endomorphisms and by L(V, H) the space of linear operators from V into H, this entails

π l H L(V ) ≤ 1 + α l Id -π l V L(V,H) ,
where Id is the identity operator and the inverse inequality constant is now multiplied by a projection error term which can be conjectured to vanish in various cases when increasing l, since V is more regular than H. Hence, we can transform any estimate with v -π l H v into an estimate with v -π l V v . Indeed, as π l H and π l V project onto the same subspace we have

v -π l H v = (Id -π l H )(v -π l V v) = (Id -(π l H -π l V ))(v -π l V v). Since (π l H -π l V )v = π l H (v -π l V v) for all v ∈ V , we remark that (π l H -π l V )v = π l H v if v ∈ (V l ) ⊥ , 0 if v ∈ V l , denoting by (V l ) ⊥ the V -orthogonal complement of V l . Hence, π l H -π l V L(V ) ≤ π l H L(V ) . Defining ρ l = π l H L(V ) , σ l = π l H -π l V L(V )
, we can then convert the projector by writing

v -π l H v ≤ (1 + σ l ) v -π l V v (11)
≤ (1 + ρ l ) v -π l V v , (12)
from which we directly infer the following estimate.

Corollary 1. For all T > 0,

u -u l L 2 (0,T ;V ) ≤ C |π l H u 0 -u l 0 | + (1 + σ l ) u -π l V u L 2 (0,T ;V ) .
However, we do not formally assert that, for a general reduction space, neither ρ l nor σ l have a bounded behavior with respect to l. This behavior is likely to be dependent on the specific types of variational problem and Galerkin reduction considered, and can be numerical assessed when no analytical treatment is at hand.

Note that the last term in the right-hand side of (1) is the quantity that is in fact minimized in the construction of POD subspaces. Furthermore, for POD reduction subspaces the sequences (σ l ) and (ρ l ) remain bounded in a large class of situations, see Section 3.5 for some theoretical insight.

Extension to a non-linear parabolic equation

We formally introduce the abstract non-linear parabolic equation

d dt (u(t), v) + a(u(t), v) = (f (t, u(t)), v), ∀v ∈ V, (13)
u(0) = u 0 . (14)
Unlike in Equation ( 5), f is some [0, T ] × V → V function. The general theory is very delicate.

Especially the solution may explode in finite time. We provide the following proposition, where we assume that f is Lipschitz-continuous in the second variable. As proven in the appendix, this guarantees, for any T > 0, the well-posedness of Equations ( 13)-( 14) in the same spaces as in the linear case.

Proposition 8. Assume u 0 ∈ H, f ∈ C([0, T ] × H; H)
, and that f is L-Lipschitz continuous in its second variable, i.e. that there exists a constant L such that

∀t ∈ [0, T ], ∀h 1 , h 2 ∈ H, |f (t, h 1 ) -f (t, h 2 )| ≤ L|h 1 -h 2 |.
Assume also that the embedding V → H is compact. Then there exists a unique solution u of Eqs. ( 13)-( 14) such that

u ∈ L 2 (0, T ; V ) ∩ C([0, T ]; H), du dt ∈ L 2 (0, T ; V ).
Let now u l be the spatial Galerkin approximation of

u in V l u l (t) ∈ V l , ( 15
) d dt (u l (t), v l ) + a(u l (t), v l ) = (f (t, u l (t)), v l ), ∀v l ∈ V l , (16)
u l (0) = u l 0 . (17)
More simply, with the Peano existence theorem, we obtain the following result.

Proposition 9. Assume u l 0 ∈ V l , f ∈ C([0, T ] × H; H)
and f is L-Lipschitz continuous in its second variable. Then there exists a unique solution u l of (15)-( 17) such that

u l ∈ C 1 ([0, T ]; V l ).
The proof of this result can also be seen as contained in that of Proposition 8, proven in the appendix.

We now show the following result for the reduction error.

Proposition 10. For all T > 0,

u -u l L 2 (0,T ;V ) ≤ C 1 (L, T ) |π l H u 0 -u l 0 | + C 2 (L) u -π l H u L 2 (0,T ;V ) , (18)
where, for all L > 0,

C 1 (L, T ) = Ce LT , C 2 (L) = C(L + 1).
In addition, we have

u -u l L 2 (0,T ;V ) ≤ C 1 (L, T ) |π l H u 0 -u l 0 | + (1 + σ l )C 2 (L) u -π l V u L 2 (0,T ;V ) . (19)
Moreover, under the condition

L < ca C 2 Ω
, we have the improved constants

C 1 (L, T ) = C 1 (L) = C ca C 2 Ω -L , C 2 (L) = C ca C 2 Ω -L ,
where C 1 is now independent of T .

Proof. We split u -u l into two parts u -u l = p l + q l , where p l = u -π l H u and q l = π l H u -u l . Using the same property of π l H as in the proof of Prop. 7, we obtain

dq l dt (t), v l + a(q l (t), v l ) = -a(p l (t), v l ) + (f (t, u(t)) -f (t, u l (t)), v l ).
Taking v l = q l (t), and integrating on [0, t], 0 ≤ t ≤ T , we obtain

1 2 |q l (t)| 2 + c a q l 2 L 2 (0,t;V ) ≤ 1 2 |q l (0)| 2 + C a p l L 2 (0,t;V ) q l L 2 (0,t;V ) +L t 0 |u(s) -u l (s)| • |q l (s)| ds ≤ 1 2 |q l (0)| 2 + C a p l L 2 (0,t;V ) q l L 2 (0,t;V ) +L t 0 (|p l (s)| + |q l (s)|) • |q l (s)| ds ≤ 1 2 |q l (0)| 2 + (C a + LC 2 Ω ) p l L 2 (0,t;V ) q l L 2 (0,t;V ) (20)
+L q l 2 L 2 (0,t;H) .

Let us first assume

L < ca C 2 Ω
. Hence, for t = T and using the continuous embedding V → H, Eq. ( 20) entails

(c a -LC 2 Ω ) q l 2 L 2 (0,T ;V ) ≤ 1 2 |q l (0)| 2 + (C a + c a ) p l L 2 (0,T ;V ) q l L 2 (0,T ;V ) .
Using Young's inequality, we can then conclude as in Prop. 7.

Let us now consider the general case for L. By Young's inequality on [START_REF] Maday | A priori convergence theory for reduced-basis approximations of single-parameter elliptic partial differential equations[END_REF],

|q l (t)| 2 + c a q l 2 L 2 (0,t;V ) ≤ |q l (0)| 2 + (C a + LC 2 Ω ) 2 c a p l 2 L 2 (0,T ;V ) + 2L q l 2 L 2 (0,t;H) . (21)
Then, we use Gronwall's inequality for t → |q l (t)| 2 , which leads to

|q l (t)| 2 ≤ e 2Lt |q l (0)| 2 + (C a + LC 2 Ω ) 2 c a p l 2 L 2 (0,T ;V ) .
Finally, re-incorporating this estimate in [START_REF] Prud'homme | A mathematical and computational framework for reliable real-time solution of parametrized partial differential equations[END_REF] gives, for t = T ,

q l 2 L 2 (0,T ;V ) ≤ e 2LT c a |q l (0)| 2 + (C a + LC 2 Ω ) 2 c a p l 2 L 2 (0,T ;V ) ,
and we conclude for [START_REF] Kunisch | Galerkin proper orthogonal decomposition methods for a general equation in fluid dynamics[END_REF] as in Prop. 7.

Of course [START_REF] Kunisch | Proper orthogonal decomposition for optimality systems[END_REF] directly follows like in Corollary 1.

Galerkin estimates for the wave-like equation

Let us now consider the wave-like equation

d 2 dt 2 (y(t), v) + a(y(t), v) = (f (t), v), ∀v ∈ V, (22)
y(0) = y 0 , dy dt (0) = ẏ0 , (23)
for which we have the classical existence and uniqueness results [8, XVIII, §5.5.2, Th. 1, §5.5.3, Th. 2].

Proposition 11. We assume f ∈ L 2 (0, T ; H), y 0 ∈ V and ẏ0 ∈ H. Then there exists a unique solution y of Eqs. ( 22)-( 23) such that

y ∈ C([0, T ]; V ), dy dt ∈ C([0, T ]; H), d 2 y dt 2 ∈ L 2 (0, T ; V ).
As in Section 3.1, we formally introduce the spatial Galerkin approximation y l of y

y l (t) ∈ V l , (24)
d 2 dt 2 (y l (t), v l ) + a(y l (t), v l ) = (f (t), v l ), ∀v l ∈ V l , (25)
y l (0) = y l 0 , dy l dt (0) = ẏl 0 . (26)
And the following holds [8, XVIII, §5.3.1, Lemma 2].

Proposition 12. We assume f ∈ L 2 (0, T ; H), y 0 ∈ V l and ẏ0 ∈ V l . Then there exists a unique solution y l of Eqs. ( 24)-( 26) such that

y l ∈ C 1 ([0, T ]; V l ), d 2 y l dt 2 ∈ L 2 (0, T ; V l
). The error estimate between the solutions of Problems ( 22)-( 23) and ( 24)-( 26) is given by the following Proposition.

Proposition 13. For all T > 0,

y -y l L 2 (0,T ;V ) + d dt (y -y l ) L 2 (0,T ;H) ≤ C √ T ( y 0 -π l H y 0 + π l H y 0 -y l 0 + |π l H ẏ0 -ẏl 0 |) + y -π l H y L 2 (0,T ;V ) + (1 + T ) d dt y -π l H y L 2 (0,T ;V ) , (27)
and

y -y l L 2 (0,T ;V ) + d dt (y -y l ) L 2 (0,T ;H) ≤ C √ T ( y 0 -π l H y 0 + π l H y 0 -y l 0 + |π l H ẏ0 -ẏl 0 |) + (1 + σ l ) y -π l V y L 2 (0,T ;V ) + (1 + T ) d dt (y -π l V y) L 2 (0,T ;V ) . (28)
Proof. We split y -y l into two parts y -y l = p l + q l , where p l = y -π l H y and q l = π l H y -y l . Since q l ∈ V l , and using the definition of y l , q l verifies the variational equation

d 2 dt 2 (q l (t), v l ) + a(q l (t), v l ) = -(f (t), v l ) + d 2 dt 2 (π l H y(t), v l ) + a(π l H y(t), v l ), ∀v l ∈ V l .
The projector π l H verifies (π l H y(t), v l ) = (y(t), v l ), so that using the definition of y we get

d 2 q l dt 2 (t), v l + a(q l (t), v l ) = -a(p l (t), v l ).
We infer the energy balance by taking v l = dq l dt (t), viz.

1 2 d dt dq l dt 2 + q l 2 a (t) = -a p l (t), dq l dt (t) .
Performing an integration by parts over time and using Young's inequality in the right-hand side, we have

dq l dt (t) 2 + (1 -η) q l (t) 2 a ≤ dq l dt (0) 2 + q l (0) 2 a + 2a(p l (0), q l (0)) + 1 η p l (t) 2 a + θ q l 2 L 2 (0,T ;a) + 1 θ dp l dt 2 L 2 (0,T ;a)
.

By integration on [0, T ] again and taking η = 1 4 , θ = 1 4T , we get

dq l dt 2 L 2 (0,T ;H) + q l 2 L 2 (0,T ;a) ≤ C T p l (0) 2 a + q l (0) 2 a + dq l dt (0) 2 + p l 2 L 2 (0,T ;a) + T 2 dp l dt 2 L 2 (0,T ;a)
.

Using the properties of the scalar product a, we get back to the • norm, and the triangular inequality y(t) -y l (t) ≤ p l (t) + q l (t) ends the proof for [START_REF] Veroy | Reduced-basis approximation of the viscous Burgers equation: rigorous a posteriori error bounds[END_REF], whence (28) directly follows.

Boundedness of (σ l )

As already mentioned, we need some characterization of the behavior of the sequences (ρ l ) l≥1 and (σ l ) l≥1 in order for the above estimations to be meaningful. To provide some insight into this issue we give some examples of reduction subspaces for which these sequences can be proven to be bounded. Let us start by showing this boundedness when the Galerkin subspace is given by finite element discretization procedures. To fix the ideas we consider a standard P 1 discretization, but this result can be extended with ease to most other finite element procedures. Proposition 14. Let Ω be an open convex polyhedral subset of R 2 or R 3 . Let H = L 2 (Ω) and V = H 1 0 (Ω). Let (T h ) h>0 be a quasi-uniform family of triangulations of Ω, and V h the P 1 -Lagrange finite element subspace of V built on T h , with π h,H the H-orthogonal projector onto V h . Then (π h,H ) h>0 is bounded in L(V ), i.e.

∀h > 0, ∀v ∈ V, π h,H v ≤ C v . (29)
Proof. Let us introduce a family of Clément interpolation operators (C h ) h>0 associated with (T h ) h>0 and uniformly bounded from V to V h [START_REF] Clément | Approximation by finite element functions using local regularization[END_REF]. Since (T h ) h>0 is quasi-uniform, an inverse inequality holds [5, Th. 3.2.6], so that

π h,H v ≤ π h,H v -C h v + C h v ≤ Ch -1 |π h,H v -C h v| + C h v .
Remark that, by the characterization of an orthogonal projector,

|π h,H v -C h v| ≤ |v -π h,H v| + |v -C h v| ≤ 2|v -C h v|. ( 30
)
Now we use the property

∀h > 0, ∀v ∈ V, |v -C h v| ≤ Ch v , (31)
and the boundedness of (C h ) h>0 in L(V ) [START_REF] Clément | Approximation by finite element functions using local regularization[END_REF]Th. 2]. This shows our result.

We now consider spectral analysis, namely, taking Galerkin subspaces provided by the eigenmodes of the bilinear form a. We thus assume the embedding V → H to be compact, which is satisfied when Ω is bounded, H = L 2 (Ω) and V = H 1 0 (Ω). Then there exists a Hilbertian basis of H, (w i ), characterized by

a(w i , v) = ω 2 i (w i , v), 0 < ω 1 ≤ ω 2 ≤ . . . , ω i -→ i→∞ +∞.
Introducing wi = 1 ω i w i , ( wi ) i≥1 is a Hilbertian basis of V for the scalar product associated with a.

Proposition 15. Assuming V l = Span(w 1 , . . . , w l ), the sequences (ρ l ) and (σ l ) are bounded.

Proof. We remark a(v, wi ) wi = (v, w i )w i .

Summing this identity from 1 to l directly entails that π l H = π l a , where π l a is the a-orthogonal projector of V onto V l . Moreover

c 1/2 a π l a v ≤ π l a v a ≤ v a ≤ C 1/2 a v ,
which leads to

π l a L(V ) ≤ C a c a 1/2
. and the property π l H = π l a concludes the proof.

As a third example, we will consider the case of the POD subspaces arising from the analysis of the homogeneous wave-like equation. The following result is very straightforward to establish by decomposing the solution on the eigenmodes. We also refer to [START_REF] Feeny | On the physical interpretation of proper orthogonal modes in vibrations[END_REF] for related discussions.

ϕ i (T ) -→ T →∞ w σ(i) ,
where σ describes a certain reordering determined by the initial conditions. Therefore

ρ(l, T ) = π l H (T ) L(V ) -→ T →∞ C.

Numerical validations

In this section, we provide some numerical validations of the above error estimates for some examples of one-dimensional problems. As in the rest of the paper, we only consider the case of self-reduction, i.e. when the reduction space we use is the POD space generated from the trajectory of the reference solution u itself. In particular, we aim at assessing whether or not the sequences (σ l ) and (ρ l ) are bounded in several examples. Of course, since the reference solution is needed to compute the POD space, it is mostly a theoretical study on synthetic data. However, this is an important first step before tackling the practical situation of parametric variations, when a unique POD space is used to reduce a family of solutions. This issue will be addressed in forthcoming papers.

Discretization and corresponding reduction for parabolic problem

Here, we consider the reduction of ( 13)-( 14) with the one-dimensional non-linear equation

∂ t u -∂ 2 xx u = f (t, u) in (0, T ) × (0, 1), u(t, 0) = u(t, 1) = 0, u(0, x) = u 0 (x) in (0, 1),
where now f is simply a [0, T ] × R → R function. In the sequel, H = L 2 (0, 1) and V = H 1 0 (0, 1). We keep the notations (•, •) and a(•, •) for their respective scalar products.

Semi-discrete solution and reduced form

Let u h be the P 1 approximation of u on the regular mesh (

x i ) N h i=1 x i = ih, 1 ≤ i ≤ N h , h = 1 N h + 1 ,
associated with the basis of shape functions (e i ) N h i=1 . The discrete solution u h is defined by

d dt (u h (t), e i ) + a(u h (t), e i ) = (f (t, u h (t)), e i ), 1 ≤ i ≤ N h , (32)
u h (0, x) = u h,0 (x). ( 33
)
This discrete solution is the reference solution with which the reduced solutions will be compared. However, the POD basis (ϕ 1 , . . . , ϕ l ) itself will be constructed based on the fully-discrete solution u n h described below. Nevertheless, we emphasize that we do not consider time discretization issues in this paper, hence in our numerical trials we choose the time step "sufficiently small" for the discrete solution to be converged in time.

The corresponding reduced form u

l h of u h satisfies d dt (u l h (t), ϕ k ) + a(u l h (t), ϕ k ) = (f (t, u l h (t)), ϕ k ), 1 ≤ k ≤ l, (34)
u l h (0, x) = u l h,0 (x). ( 35
)
As before, we have local existence and uniqueness of the solutions u h and u l h in the classical sense.

We also similarly introduce the L 2 (0, 1)-orthogonal and H 1 0 (0, 1)-orthogonal projectors π l L 2

and π l H 1 0 from V h onto V l , and the corresponding sequences

ρ l = π l L 2 L(H 1 0 ) , σ l = π l L 2 -π l H 1 0 L(H 1 0 ) , that still verify σ l ≤ ρ l .
We can then directly adapt Proposition 10.

Proposition 17. Assume u h,0 ∈ V h , u l h,0 ∈ V l , f ∈ C([0, T ] × R; R)
, and f is Lipschitz continuous in its second variable. Then there exists unique classical and global solutions u h and u l h of Equations (32)-( 33) and (34)-( 35), respectively. Moreover, for all T > 0,

u h -u l h L 2 (0,T ;H 1 0 ) ≤ C π l L 2 u h,0 -u l h,0 L 2 (Ω) + u h -π l L 2 u h L 2 (0,T ;H 1 0 ) ≤ C π l L 2 u h,0 -u l h,0 L 2 (Ω) +(1 + σ l ) u h -π l H 1 0 u h L 2 (0,T ;H 1 0 ) . (36)
Note that -if we assume that the POD basis is constructed in a continuous-time discretespace setting -the last term in this error estimate directly corresponds to the POD remainder, recall (4), hence it is perfectly controlled in the POD construction itself.

Full discretization

We use the classical θ-method as a time discretization scheme. In order to compute the reference solution u h , we need the non-reduced mass matrix M and stiffness matrix K M = [(e j , e i )] 1≤i,j≤N h , K = [a(e j , e i )] 1≤i,j≤N h , and the reaction term application

F : R N h → R N h of coefficient (F (t, β)) i = 1 0 f t, N h k=1 β k e k (x) e i (x) dx. Then the vector U h (t) ∈ R N h concatenating the coordinates of u h (x, t) in (e i (x)) N h i=1 satisfies M Uh (t) + KU h (t) = F (t, U h (t)), U h (0) = U h,0 .
Next we apply a semi-implicit time scheme by the θ-method

1 ∆t M U n+1 h -U n h + K θU n+1 h + (1 -θ)U n h = θF (t n+1 , U n+1 h ) + (1 -θ)F (t n , U n h ),
leading to a non-linear problem in U n+1 h once U n h is known, which we can solve for using a Newton algorithm.

For the reduced solutions u l h , we follow exactly the same path (spatial discretization U l h (t), full discretization (U l,n h )), except that we substitute the POD basis (ϕ i ) l i=1 for the finite element basis (e i ) N h i=1 . This gives the reduced mass and stiffness matrices M l and K l , and the reduced reaction term F l . We emphasize that although these reduced matrices are of limited size, they are full. We call Φ l the matrix

Φ l = [ϕ 1 , . . . , ϕ l ] ∈ R N h ×l ,
where vectors ϕ i are expressed as column vectors of coordinates in (e i ) N h i=1 . Then we obtain the following relations between reduced and non-reduced operators

M l = (Φ l ) T M Φ l , K l = (Φ l ) T KΦ l , F l (t, β l ) = (Φ l ) T F (t, Φ l β l ), d β l F l (t, β l ) = (Φ l ) T d β F (t, Φ l β l )Φ l ,
where d β l F l and d β F denote the differential quantities needed in the Newton algorithm computations.

Sharpness indicators for the new estimates

Here, we define the quantities that we need to check to ensure the sharpness of the new estimates. Note first that the POD eigenvalues λ i typically decrease exponentially, hence the maximum POD rank to be considered is set as

λ lmax+1 ≤ 10 -12 λ 1 ≤ λ lmax .
in order to preserve sufficient K-orthogonality of the POD basis (ϕ i ) l i=1 when we perform the diagonalization of the covariance matrix. Indeed, since the covariance matrix is ill-conditioned, this orthogonality tends to rapidly deteriorate with l and we should preserve

(Φ l ) T KΦ l -Id l ≤ ε tol .

Summary of the estimation chain

In Prop. 17, we mainly handle three error terms:

• the reduction error R(l) R(l) = u h -u l h L 2 (0,T ;H 1 0 ) ; • the L 2 -projection error Q(l) Q(l) = u h -π l L 2 u h L 2 (0,T ;H 1 0 ) ;
• and the H 1 0 -projection error P (l)

P (l) = u h -π l H 1 0 u h L 2 (0,T ;H 1 0 ) , (37)
that coincides, in this situation of self-reduction, with the POD remainder ε(l)

ε(l) = i>l λ i 1/2
.

Note that we do not need the reduced solution u l h to compute Q(l) nor P (l). Moreover, we point out that these quantities -except for P (l) which can be obtained as a by-product of the covariance computation -are auxiliary quantities only computed to evaluate the reduction performance and accuracy.

If we prescribe the initial condition

u l h,0 = π l L 2 u h,0 ,
then the first term in the right-hand side of Eq. ( 36) vanishes. Thus we summarize the estimation chain by

R(l) ≤ CQ(l) ≤ C(1 + σ l )P (l).
Let us introduce the following sharpness indicator

S Gal (l) = R(l) Q(l) ,
which is clearly bounded under the assumptions of Prop. 17, but can be considered in a more general framework. By contrast, note that for the second inequality that only relies on [START_REF] Gugercin | A survey of model reduction by balanced truncation and some new results[END_REF], the bound Q(l) (1 + σ l )P (l) ≤ 1 always holds. Finally, we aim at numerically verifying, in various cases, that:

• the maximum POD rank l max is reasonably limited compared to the number of degrees of freedom of the system l N h , for the POD subspace to accurately approximate the solution, recall (37);

• the quantity max 1≤l≤lmax ρ l , that is an upper bound of max 1≤l≤lmax σ l , remains small;

• the indicator S Gal (l), 1 ≤ l ≤ l max , remains numerically bounded, especially in cases of strong non-linearities.

Computation of the (ρ l ) and (σ l ) sequences

In order to compute ρ l , it is useful to manipulate the L 2 -orthonormal basis (ψ k ) k that results from a Gram-Schmidt L 2 -orthonormalization on the H 1 0 -orthonormal POD basis (ϕ k ). Thus

V l = Span(ψ 1 , . . . , ψ l ),
with ψ i independent of the POD rank l. Denoting by Ψ l the matrix

Ψ l = [ψ 1 , . . . , ψ l ] ,
where the ψ i elements are expressed as column vectors of coordinates in (e i ) N h i=1 , we notice that π l L 2 has the following matrix in (e

i ) N h i=1 Π l L 2 = Ψ l (Ψ l ) T M ∈ R N h ×N h .
We use the definition of ρ l ρ 2 l = sup

v∈V h \{0} 1 0 [π l L 2 v] (x) 2 dx 1 0 v (x) 2 dx = sup β∈R N h \{0} β T (Π l L 2 ) T KΠ l L 2 β β T Kβ .
Then ρ l is the solution of the "largest K-eigenvalue" problem

ρ l = sup ω ≥ 0 | ∃β ∈ R N h \ {0}, (π l L 2 ) T Kπ l L 2 β = ω 2 Kβ .
Similarly, let Φ l be the matrix

Φ l = [ϕ 1 , . . . , ϕ l ] ,
and Π l L 2 be the matrix of the truncated projector

(π l L 2 -π l H 1 0 ) Π l L 2 = Ψ l (Ψ l ) T M -Φ l (Φ l ) T K.
Then σ l is the solution of the problem

σ l = sup ω ≥ 0 | ∃β ∈ R N h \ {0}, ( Π l L 2 ) T K Π l L 2 β = ω 2 Kβ .
These properties will allow the numerical evaluation of the sequences.

Numerical experiments and validation for the parabolic problems

We present three numerical cases of POD reduction on the generic discrete parabolic equation ( 34)-( 35). The corresponding parameters are gathered in Table 1. In all these cases, we take θ = 2 3 in the θ-method for time discretization, and N h = 100 for the spatial discretization.

Case

Case A Case B Case C nb. timesteps 10 3 10 3 800 ∆t 10 -4 10 -4 10 -5 u 0 (x)

1 [ 1 3 , 2 3 ] (x) 27 4 x 2 (1 -x) 27 4 x 2 (1 -x) f (t, u) 9u 10u 2 100u 2
Table 1: Cases of study for the reduction of parabolic equations

Case A: Lipschitz continuous reaction term

Case A satisfies the assumptions of Prop. 17 since f is linear with respect to u. We display the corresponding results in Figs. 1 and2. In these figures the POD rank l varies from 1 to l max . Figure 1, as an indication, shows the shape in space and time of the non-reduced solution u h , and a numerical comparison between the indicators P (l), Q(l) and R(l).

Figure 2 displays the sequence of POD constants ρ l and their truncated versions σ l , together with the sharpness indicator S Gal (l).

In this simple case, all our verifications are successful, namely, • the POD remainder decreases at an exponential rate, and l max = 10;

• the POD constants ρ l are of magnitude O(1) and remain bounded with l. Also, the improvement provided by σ l is limited;

• as expected with Prop. 17, the sharpness indicator is bounded and of small value, viz.

S Gal ∈ [0.3, 0.8].

Cases B and C: super-linear reaction term

By contrast, the other two cases B and C are beyond the assumptions of Prop. 17 because we consider super-linear reaction terms. Moreover, while case B remains bounded, case C appears to explode in finite time, which is why we reduced the time range in this case while keeping a similar number of time steps, see Table 1. Even though our above error estimates do not hold in these cases, we can still compute the same numerical error quantities for illustrative purposes. Case B is reported on in Figs. 3 and4, and case C in Figs. 5 and6.

In fact, the maximum POD rank as well as the behavior and magnitude of the indicators ρ l , σ l , and S Gal reveal no significant difference compared to case A. Furthermore, the reduction error still decreases as fast as the POD remainder.

Numerical assessment of the wave equation reduction

Considering now the 1D homogeneous wave equation,

∂ 2 tt y -c 2 ∂ 2 xx y = 0 in (0, T ) × (0 , 1), y(t
, 0) = y(t, 1) = 0, y(0, x) = y 0 (x) in (0, 1), ∂ t y(0, x) = ẏ0 (x) in (0, 1), we report on the numerical values obtained for the various error terms. We discretize in space with finite elements on a regular mesh, and in time with a Newmark scheme according to the classical parameters β = 1 4 and γ = 1 2 , see e.g. [START_REF] Raviart | Introduction à l'Analyse Numérique des Equations aux Dérivées Partielles[END_REF]. We take a regular cutoff function for y 0 (x), and ẏ0 (x) = 0. The corresponding results are shown in Figs. 7 and8.

We verify that the POD basis is very close to a set of H 1 0 (0, 1)-eigenmodes ( wi ) of the Dirichlet Laplacian as substantiated in Section 3.5. This also explains why σ l is much lower than ρ l , since the L 2 and H 1 0 projectors onto eigenspaces coincide. Note that the estimate of Prop. 16 contains the first-order time derivative ∂ ∂t (u-π l H 1 0 u) which is not controlled by the POD construction. Nevertheless, we observe from Figure 7 that the POD reduction is very effective, and indeed converges nearly-exponentially with the POD-rank.

Reduction of a complex system: a biomechanical heart model

In this section, we test our reduction technique and estimates with the 3D continuum mechanics model of a beating heart. The electromechanical heart model and its discretization are described in [START_REF] Sainte-Marie | Modeling and estimation of the cardiac electromechanical activity[END_REF], and a validation of the model by confrontation with clinical data is given in [START_REF] Chabiniok | Validation of a biomechanical heart model using animal data with acute myocardial infarction[END_REF].

A model describing the three-dimensional electromechanical behavior of the heart requires several important ingredients, namely,

• a constitutive law accounting for both the active and passive aspects in the behavior of the muscle fibres, • a representation of the electrical activation -the input in the constitutive law -that can be obtained from modeling approaches of various types and complexities,

• a geometrical (or "anatomical") description of the myocardium incorporating the fibre directions,

• a model of the blood circulation inside and outside of the heart cavities, and also a model describing the opening and closure of the valves that separate the cavities from each other and from the external circulation.

Let y be the displacement field, the deformation gradient is defined by F = I + ∇ y and the Green-Lagrange strain tensor is given by e = 1 2 (F T • F -I) and J = det F . Based on the above modeling ingredients, the second Piola-Kirchhoff stress tensor Σ contains the active cardiac fibre law, a viscous stress component and a hyperelastic potential accounting for passive effects, these components being combined by means of a rheological model of Hill-Maxwell type.

Using a total Lagrangian formulation and denoting by Ω H the reference domain corresponding to cardiac tissue, while the part of the boundary corresponding to ventricular endocardium is denoted by Γ, the principle of virtual work then gives

Ω H ρ ÿ • v dΩ + Ω H Σ : d y e • v dΩ + Γ P 0 ν • F -1 • vJ dΓ = 0 ∀v ∈ V,
where V denotes a suitable space of displacement test functions, ρ the mass per unit volume, d y e the differential of the Green-Lagrange strain tensor with respect to the displacement, while P 0 is a prescribed intraventricular pressure.

For the simulations presented hereafter an idealized left ventricle embedded with active fibers has been considered. The discretization is performed with P 1 -Lagrange finite elements in space (with about 1000 degrees of freedom), and a Newmark scheme in time [START_REF] Bathe | Finite Element Procedures[END_REF]. We show some snapshots of the solution for the full finite element model in Fig. 9.

Although we do not have a theoretical estimate for the reduction error in this complex nonlinear case, the three error terms appearing in the linear estimation chain still feature excellent decreasing rate and correlation, see Figs. 10 and11. Analyzing our indicators reveals that their magnitude may slightly differ from the linear one-dimensional case, but again shows the effectiveness of the POD reduction, namely,

• l max = 36 ;

• ρ l and σ l are almost identical, and numerically bounded;

• the sharpness indicator established for the linear case is still bounded, and more precisely

S Gal ∈ [0.6, 1.0].
In Figure 10, we also display the evolution in time of the relative residual e l defined as

e l (t) = y(t) -y l (t) L ∞ (Ω) y C([0,T ];L ∞ (Ω)) . ( 38
)
We observe an excellent behavior of e l , which roughly decreases by an order of magnitude for each addition of 10 modes in the POD basis.

Conclusion

We have proposed Galerkin estimates for the Proper Orthogonal Decomposition reduction of some classical PDEs. The numerical implementation of the reduction and some verifications were presented in this article. We have also demonstrated reduced simulations of a complex three-dimensional electromechanical model of the heart, where the validity of similar Galerkin estimates is numerically verified, although no formal proof can be given in this case.

A special emphasis was placed on the derivation of POD-reduction error estimates in convenient norms and which can be controlled in the construction of the POD basis.

The present study can be extended in many directions. Firstly, as far as POD reduction of PDEs is concerned, one of the major difficulties lies in achieving stability of the POD basis with respect to e.g. parameter variations, initial and boundary conditions, and so on. This subject needs to be further investigated. Secondly, filtering and estimation techniques for inverse modeling are extremely costly from a computational standpoint. This justifies -or even often requires -the use of POD-based reduced models and/or reduced filters and hence, the derivation of error estimates for such problems is crucial.

A Existence and uniqueness of solutions of variational equations with a Lipschitz continuous reaction term

Although some results pertaining to this type of problem exist in the literature, for the sake of completeness we provide the sketch of a self-contained proof for the specific result that we need in our case, namely, Proposition 8. Since we assume the embedding V → H to be compact, we can use the Hilbertian bases (w i ) and ( wi ) of H and V -respectively -made up by the eigenvectors, as already introduced in Section 3.5. Let W k , k ≥ 1, denote the subspace W k = Span(w 1 , . . . , w k ) = Span( w1 , . . . , wk ), and P k the orthogonal projector from H onto W k , i.e.

P k h = k i=1 (h, w i )w i , ∀h ∈ H.
It coincides with the orthogonal projector from (V, a) onto W k defined by

P a k v = k i=1 a(v, wi ) wi = P k v, ∀v ∈ V.
We will show the existence result by a Galerkin approach using the sequence of eigenspaces.

Proposition 18. For all k ≥ 1, there exists a unique global solution

u k ∈ C 1 (R + ; W k ) such that d dt (u k (t), v k ) + a(u k (t), v k ) = (f (t, u k (t)), v k ), ∀v k ∈ W k , (39)
u k (0) = P k u 0 .

Moreover, for all T > 0, (u k ) is bounded in C([0, T ]; H).

Proof. Let us first prove uniqueness. Consider two solutions u

1 k , u 2 k , then d dt |u 1 k -u 2 k | 2 (t) ≤ 2L|u 1 k -u 2 k | 2 (t),
and since u 1 k (0) = u 2 k (0), we infer u 1 k = u 2 k . We now tackle the global existence. By the Peano existence theorem (see e.g. [10, 2.4.4]), we have local existence. By uniqueness, the maximum time of existence T k ∈ (0, ∞] is well-defined. We test the equation with

v k = u k (t), i.e. 1 2 d dt |u k (t)| 2 + u k (t) 2 a = (f (t, u k (t)) -f (t, 0), u k (t)) + (f (t, 0), u k (t)). By Young's inequality (here |(f (t, 0), u k (t))| ≤ L|u k (t)| 2 + 1/(4L)|f (t, 0)| 2 ), d dt |u k (t)| 2 ≤ 4L|u k (t)| 2 + 1 2L |f (t, 0)| 2 .
Then, by Gronwall's lemma, for all T > 0, lim

t→T - |u k (t)| ≤ C(T ) < ∞, with C(T ) = |u 0 | + 1 2L T 0 e 4L(T -s) |f (t, 0)| 2 dt.
Finally, on the one hand we deduce global existence, i.e. T k = ∞ (e.g. [10, 2.4.3, 2.4.4]), and on the other hand we get the boundedness in C([0, T ]; H).

In order to show that (u k ) is a Cauchy sequence in the Banach spaces C([0, T ]; H) and L 2 (0, T ; V ), let us consider the decomposition

u k+p -u k = P k (u k+p -u k ) + (Id -P k )u k+p . ( 40
)
For the first term in the right-hand side, we get the following estimates.

Lemma 1. For all 0 ≤ t 0 ≤ t 1 and all k, p ≥ 1, We conclude using the Lipschitz character of f and the boundedness of (u k ) in C([0, T ]; H).

max P k (u k+p -u k ) C([t 0 ,t 1 ];H) , √ 2c a P k (u k+p -u k ) L 2 (t 0 ,t 1 ;V ) ≤ |(u k+p -u k )(t 0 )| + 2L(t 1 -t 0 ) u k+p -u k C([t 0 ,t 1 ];H) . Proof. For all v k ∈ W k , d dt (u k+p (t) -u k (t), v k ) + a(u k+p (t) -u k (t), v k ) = (f (t, u k+p (t)) -f (t, u k (t)), v k ). (41
We are ready to show the first convergence result.

Proposition 19. For all T > 0, (u k ) converges in C([0, T ]; H) to some limit u.

Proof. Let τ = 1 8L and j ≥ 1. Lemmas 1 and 2 lead to (43)

We prove by induction that the statement P(j) : (u k ) is a Cauchy sequence in C([0, jτ ]; H) holds for all j ≥ 1. We easily show P(1). Assume now that P(j -1) holds for some j ≥ 2. Let u be the limit of (u k ) in C([0, (j -1)τ ]; H). Then we decompose in (43) |(Id -P k )u k+p ((j -1)τ )| ≤ |(u k+p -u)((j -1)τ )| + |(Id -P k )u((j -1)τ )|, which proves that (u k ) is a Cauchy sequence in C([(j -1)τ, jτ ]; H), and hence that P(j) holds.

Remark that we directly obtain u(0) = u 0 .

Next, we get an estimate for the second term in (40) in the L 2 (0, T ; V ) norm. where the sequence g k,p is defined as g k,p = (Id -P k )f (•, u k+p ) L 2 (0,T ;H) , and verifies ∀ε > 0, ∃k 0 , ∀k ≥ k 0 , ∀p ≥ 0, g k,p ≤ ε.

Proof. We consider again (42) that we rewrite as Moreover, by the Parseval theorem, ḡ2 k,n and ḡ2 k are the remainders of some positive converging series, so that in particular ḡk -→ k→∞ 0, and (ḡ k,n ) k is a decreasing sequence for each n. Finally for ε > 0, there exists K such that ḡK ≤ ε 2 , and n 0 such that for all n ≥ n 0 , |ḡ K,n -ḡK | ≤ ε 2 . We conclude by taking k 0 = max(K, n 0 ). This entails the second convergence result.

Proposition 20. For all T > 0, (u k ) converges in L 2 (0, T ; V ) and its limit is u. Using (46), (u k ) is also a Cauchy sequence in L 2 (0, T ; V ). Let ũ be its limit. Since L 2 (0, T ; V ) and C([0, T ]; H) are both continuously embedded in L 2 (0, T ; H), then ũ = u.

We can finally conclude.

Proof of Proposition 8. Using the previous convergence results we can reinterpret the limit u as satisfying Equation [START_REF] Hinze | Error estimates for abstract linear-quadratic optimal control problems using proper orthogonal decomposition[END_REF] in the distribution sense. Given the regularity of u, we have that -a(u(t), v) + (f (t, u(t)), v) is in L 2 (0, T ) for any v ∈ V , hence it directly follows that du dt ∈ L 2 (0, T ; V ). We finally prove the uniqueness as in Proposition 18.

  for a.e. s ∈ [0, T ]. We have the following proposition. Proposition 4. Cov and Cov share the same non-zero eigenvalues, with identical multiplicities. Moreover Cov is compact.

Proposition 16 .

 16 Let y be the solution of the homogeneous wave equation, namely, (22)-(23) with f = 0. Denoting by (ϕ i (T )) l i=1 the POD basis constructed with y over [0, T ], for T ∈ [0, ∞),

Figure 1 :

 1 Figure 1: Case A. Left: full solution. Right: relative errors of H 1 0 -projection, L 2 -projection and reduction.

Figure 2 :

 2 Figure 2: Case A: Left: POD sequences (ρ l ) and (σ l ). Right: sharpness indicator (S Gal (l)).

Figure 3 :

 3 Figure 3: Case B. Left: full solution. Right: relative errors of H 1 0 -projection, L 2 -projection and reduction.

Figure 4 :

 4 Figure 4: Case B: Left: POD sequences (ρ l ) and (σ l ). Right: sharpness indicator (S Gal (l)).

Figure 5 :Figure 6 :

 56 Figure 5: Case C. Left: full solution. Right: relative errors of H 1 0 -projection, L 2 -projection and reduction.

Figure 7 : 5 Figure 8 :

 758 Figure 7: Homogeneous wave equation. Left: full solution. Right: relative errors of H 1 0projection, L 2 -projection and reduction.

Figure 9 :Figure 10 :Figure 11 :

 91011 Figure 9: Model of a ventricle: snapshots of the displacement field at the beginning (left) and 40% (right) of the first cardiac cycle for the full model.

1 2 u

 12 k+p -u k C([(j-1)τ,jτ ];H) ≤ |(u k+p -u k )((j -1)τ )| + |(Id -P k )u k+p ((j -1)τ )| + C ω k+1.

Lemma 3 .

 3 For all T > 0 and all k, p ≥ 1,(Id -P k )u k+p L 2 (0,T ;V ) ≤ C(|(Id -P k )u 0 | + g k,p ),(45)

2 L 2

 22 -P k )u k+p | 2 (t) + (Id -P k )u k+p (t) 2 a = (Id -P k )f (t, u k+p (t)), (Id -P k )u k+p (t) .Then by integration,c a (Id -P k )u k+p (0,T ;V ) ≤ 1 2 |(Id -P k )u 0 | 2 + g k,p (Id -P k )u k+p L 2 (0,T ;V ) ,and by Young's inequality, we obtain the estimate (45). Now, by strong convergence in C([0, T ]; H) and continuity of f ,|(Id -P k )f (s, u n (s))| -→ n→∞ |(Id -P k )f (s, u(s))|. Also, |(Id -P k )f (s, u n (s))| ≤ LC + |f (s, 0)| ≤ C,so that by the dominated convergence theorem, ḡk,n = (Id -P k )f (•, u n ) L 2 (0,T ;H) -→ n→∞ ḡk = (Id -P k )f (•, u) L 2 (0,T ;H) .

Proof.

  By the decomposition (40) and Lemma 1,u k+p -u k L 2 (0,T ;V ) ≤ C |(Id -P k )u 0 | + u k+p -u k C([0,T ];H) + g k,p .

)

  Testing this equation with v k = P k (u k+p -u k ), using orthogonality properties of P k , and finally integrating on [t 0 , t], t ∈ [t 0 , t 1 ], we have1 2 |P k (u k+p -u k )(t)| 2 + c a P k (u k+p -u k )(t) 2Using the diagonalisation of a, we obtain the following estimate for the second term in (40) in the C([t 0 , t 1 ]; H)-norm.Lemma 2. For all 0 ≤ t 0 ≤ t 1 and all k, p ≥ 1(Id -P k )u k+p C([t 0 ,t 1 ];H) ≤ |(Id -P k )u k+p (t 0 )| + C ω k+1.Proof. Applying now (39) for u k+p with the test function v l = (Id -P k )u k+p (t) yieldsId -P k )u k+p | 2 (t) + (Id -P k )u k+p (t) 2 a = f (t, u k+p (t)), (Id -P k )u k+p (t) . Id -P k )u k+p (t)2 a ≥ ω 2 k+1 |(Id -P k )u k+p (t)| 2 , so that by Young's inequality on the right-hand side of (42), we infer f (t, u k+p (t)), (Id -P k )u k+p (t) ≤ 1 4ω 2 k+1 |f (t, u k+p (t))| 2 + (Id -P k )u k+p (t) 2 a .

	Note that	L 2 (t 0 ,t;V ) 1 2 |((42) ≤ |P 1 2 d dt (

k (u k+p -u k )(t 0 )| 2 + L u k+p -u k 2 L 2 (t 0 ,t;H) .