NONCONFORMING MIXED FINITE ELEMENT APPROXIMATION OF A FLUID-STRUCTURE INTERACTION SPECTRAL PROBLEM

. We aim to provide a ﬁnite element analysis for the elastoacoustic vibration problem. We use a dual-mixed variational formulation for the elas- ticity system and combine the lowest order Lagrange ﬁnite element in the ﬂuid domain with the reduced symmetry element known as PEERS and introduced for linear elasticity in [1]. We show that the resulting global nonconforming scheme provides a correct spectral approximation and we prove quasi-optimal error estimates. Finally, we conﬁrm the asymptotic rates of convergence by numerical experiments.

1. Introduction. We are concerned with the computation of the free vibration modes of a coupled system consisting of an elastic structure which is in contact with an internal compressible fluid. We refer to [3,5,15,18] for the analysis of different formulations of this eigenvalue problem. Here, we follow [12,16,17] and consider a dual-mixed formulation with reduced symmetry in the solid. This leads to a symmetric saddle point problem that delivers direct finite element approximations of the stresses and that is immune to the locking phenomenon that arises in the nearly incompressible case. Recently, a Galerkin scheme based on the lowest-order Lagrange finite element in the fluid and the lowest order Arnold-Falk-Winther [2] mixed finite element in the solid has been analyzed in [17]. It was shown that such a mixed finite element Galerkin approximation is spectrally correct and provides optimal convergence error estimates for eigenvalues and eigenfunctions. Our purpose here is to show that the same order of convergence can be achieved, at a lower computational cost, when PEERS element [1] is used in the solid in association with the lowest order Lagrange finite element in the fluid. Compared with [17], the main 270 SALIM MEDDAHI AND DAVID MORA technical difficulty is related with the fact that, in this case, the Galerkin scheme is nonconforming.
The paper is organized as follows. In Section 2 we recall the mixed formulation with reduced symmetry of the fluid-structure eigenvalue problem and provide a spectral description of the corresponding solution operator. In Section 3 we introduce the mixed finite element approximation of the saddle point eigenproblem and characterize the spectrum of the discrete solution operator. In Section 4 we provide the conditions under which the numerical scheme is spectrally correct and we provide abstract convergence error estimates for the eigenfunctions and the eigenvalues. In Section 5 we establish asymptotic error estimates and finally, in Section 6, we present numerical tests and confirm that the experimental rates of convergence are in accordance with the theoretical ones.

Notations.
In all what follows we will denote the vectorial and tensorial counterparts of order n (n = 2, 3) of a given Hilbert space V by V n and V n×n respectively. We use standard notation for the Hilbertian Sobolev space H s (Ω), s ≥ 0, defined on a Lipschitz bounded domain Ω ⊂ R n and denote by · s,Ω the norms in H s (Ω), H s (Ω) n and H s (Ω) n×n .
Given two Hilbert spaces V and W and a bounded bilinear form c : V × W → R, we say that c satisfies the inf-sup condition for the pair {V, W}, whenever there exists β > 0 such that Finally, 0 stands for a generic null vector or tensor and denote by C generic constants independent of the discretization parameters, which may take different values at different places.
2. The spectral problem. We consider an elastic structure occupying a Lipschitz and polyhedral domain Ω S . We assume that the structure is fixed at ∅ = Γ D ⊂ ∂Ω S and free of stress on Γ N := ∂Ω S \ (Γ D ∪ Σ). We are interested by the simplified model in which the stress tensor σ is related to the linearized deformation tensor ε := 1 2 [∇u + (∇u) t ] through the constitutive law where λ S and µ S are Lamé coefficients, I is the identity matrix of R n×n and C : R n×n → R n×n is given by Cτ := λ S (tr τ ) I + 2µ S τ . We will also consider the rotation r := 1 2 [∇u − (∇u) t ] as a further variable. The interior fluid domain is given by a Lipschitz and polyhedral domain Ω F and the fluid-structure interface is represented by by Σ, see Figure 1. The boundary ∂Ω F of the fluid domain is the union of the interface Σ and the open boundary of the fluid Γ 0 (we don't exclude the case Γ 0 = ∅).
The spectral structural-acoustic coupled problem described in Figure 1, with natural frequencies ω, can be formulated as follows in terms of the stress tensor and the pressure (see, for instance, [5,18]): Find σ : Ω S → R n×n symmetric, r : Ω S → R n×n skew symmetric, p : Ω F → R and ω ∈ R such that, where c is the acoustic speed, g is the gravity acceleration and ρ F and ρ S represent the fluid and solid mass densities respectively. Notice that the displacement can be recovered, and also post-precessed at the discrete level, from div σ + ω 2 ρ S u = 0.
The solution operator corresponding to this eigenvalue problem is where ((σ * , p * ), r * ) ∈ Y × Q solves the source problem: for all (τ , q) ∈ Y and s ∈ Q.
Theorem 2.1. The linear operator T is well defined and bounded. Moreover, the norm of this operator remains bounded in the nearly incompressible case (i.e., when λ S → ∞). Moreover, it is clear that µ = 1 is an eigenvalue of T : Let us now rewrite the equations of problem (9)-(10) as follows: Find λ ∈ R and 0 = ((σ, p), r) ∈ Y × Q such that, where A and B are the bounded bilinear forms in Y × Q defined by To continue with the spectral description of T : Moreover, we have the direct and stable decomposition Proof. We refer to [16,Proposition A.1] for the proof of (13) and proceed as follows to obtain the splitting of a given ((σ, p), r) ∈ Y × Q according to (14). We consider the problem: The first equation shows that the linear form proves the existence of u 0 ∈ L 2 (Ω S ) n such that for all (τ , q) ∈ Y R and (s, v) ∈ Q × L 2 (Ω S ) n . The saddle point problem (15)- (16) is well-posed (see [17,Section 3]) and ((σ−σ 0 , p−p 0 ), r−r 0 ) belongs to [ker(a)×Q] ⊥ B by construction. The decomposition (14) follows then from It is clear now that the solution of the continuous eigenvalue problem (9)-(10) relies on the spectral description of To this end, we need to provide a characterization of the unique projection P : and kernel ker(a) × Q associated to the splitting (14).
In what follows,q := q − 1 |ΩF| ΩF q stands for the zero mean value component It is convenient to incorporate the divergence restriction on σ by means of a Lagrange multiplier and use a shift argument to deal properly with the affine transmission condition relating σ andp on Σ. For this purpose, given q ∈ H 1 (Ω F ), we let u ∈ H 1 (Ω S ) n and σ ∈ H(div; Ω S ) be the solution of the following linear elasticity problem: and we define the bounded linear operator E : H 1 (Ω F ) −→ W given by Eq := − σ. We notice that E provides a symmetric divergence-free extension of a given pressure field q to the solid domain. Classical regularity results for the elasticity equations in polyhedral (polygonal) domains (cf. [10,13]) ensure the existence of t S ∈ (0, 1], which depends on the geometry of Ω S and the Lamé coefficients, such that Eq ∈ H tS (Ω S ) n×n and We consider Eq := (Eq, q) ∈ Y and introduce the operator The arguments given for the wellposedness of (15)- (16) are valid for the saddle point problem (21)-(22). Moreover, it is clear from (17) The following regularity results obtained in Lemma 3.1 and Proposition 4.1 of [17] are essential for the forthcoming analysis.
As a first consequence of Lemmas 2.3 and 2.4 and the fact that follows. This permits us to announce the following spectral characterization of T .
is a sequence of finite-multiplicity eigenvalues of T which converge to 0 and the corresponding eigenspaces lie on [ker(a) × Q] ⊥ B ; moreover, the ascent of each of these eigenvalues is 1; iii) µ = 0 is an infinite-multiplicity eigenvalue of T and its associated eigenspace is For the sake of simplicity, in the forthcoming analysis we assume that T h (Ω S ) and T h (Ω F ) induce on Σ a coincident triangulation denoted Σ h . In what follows, given an integer k ≥ 0 and a subset S of R n , P k (S) denotes the space of polynomial functions defined in S of total degree ≤ k.
We consider the first order Raviart-Thomas finite element and denote by b T the usual bubble function on T ∈ T h (Ω S ). We introduce We denote by ̺ h the L 2 (Σ)-orthogonal projection onto P 0 (Σ h ) and introduce the finite element subspaces We point out that Y h ⊂ Y is not a subspace of Y. In addition, for the analysis below we will also use the space Notice that W h × Q h × U h is the lowest-order mixed finite element of the PEERS family introduced for linear elasticity by Arnold, Brezzi and Douglas (see [1]). In particular we have the inf-sup condition [1,8]: There exists β * > 0, independent of h, such that The discrete counterpart of problem (9)-(10) reads as follows: for all (τ h , q h ) ∈ Y h and s h ∈ Q h . The discrete version of the operator T is then given by for all (τ h , q h ) ∈ Y h and s h ∈ Q h . We can use the classical Babuška-Brezzi theory to prove that T h is well defined and bounded uniformly with respect to h. Indeed, we already know from [17, Lemma 2.1] that A is elliptic on the whole W × H 1 (Ω F ) (and in particular on Y h ), whereas the discrete inf-sup condition  (11)-(12) and (27)-(28) respectively. The following identity holds true, Proof. We have from (27)-(28) that for all ((τ h , q h ), s h ) ∈ Y h ×Q h . On the other hand, testing (11) with (τ , 0), (0, q) ∈ D(Ω S ) n×n × D(Ω F ) ⊂ Y yields Applying an integration by parts formula to (11) and using the last two equations we deduce that Testing now (11) with an appropriate (τ , q) ∈ Y we can show that Substituting the last identity in (31) and taking into account (30) we deduce (29).
To describe the spectrum of this operator, we proceed as in the continuous case and consider reduces to the identity, which means that here again µ h = 1 is an eigenvalue of T h with associated eigenspace ker h (a) × Q h . Let us also consider We have the following discrete analogue to Lemma 2.2.
Moreover, we have the following direct and uniformly stable decomposition Proof. Taking into account the inf-sup condition (24), the proof is similar to the one given for Lemma 2.2.
We denote by P h : Y h × Q h → Y h × Q h the unique projection with range [ker h (a) × Q h ] ⊥ B and kernel ker h (a) × Q h associated to the discrete direct splitting (32). Given ( We are now in a position to provide a characterization of the spectrum of T h and, hence, of the solutions to problem (25)-(26). Proof. See [17, Theorem 6.7] for more details. 4. Abstract convergence analysis. For the sake of brevity, we will denote in this section X := Y × Q, X := Y × Q and X h := Y h × Q h . Moreover, when no confusion can arise, we will use indistinctly x, y, etc. to denote elements in X and, analogously, x h , y h , etc. for those in X h . Finally, we will use · L(X h , X) to denote the norm of an operator restricted to the discrete subspace X h ; namely, if S : X → X, then where the consistency error is given by with T x := ((σ * , p * ), r * ).
Proof. We deduce from the well-posedness of problem (27)-(28) that the operator has a uniformly bounded inverse. It follows that (see [9]), there exists γ > 0 independent of h such that Let x * := T x ∈ X and x * h := T h x ∈ X h be the solutions of (11)- (12) and (27)-(28) respectively with data x = ((F , f ), g). The triangle inequality and (38) show that, for allỹ h ∈ X h , Denoting by A the norm of the bilinear form A and using again the triangle inequality we deduce that and the result follows from (29).
Lemma 4.1. There exists C > 0, independent of h, such that

SALIM MEDDAHI AND DAVID MORA
Proof. Given x h ∈ X h , we have that where the last equality is because both T and T h become the identity when restricted to ker h (a) × Q h . On the other hand, yields the estimate . (40) Then, if m is the multiplicity of an eigenvalue µ / ∈ {0, 1} of T , there exists exactly m eigenvalues {µ i,h , i = 1 · · · , m} of T h such that Moreover, if E(µ) is the eigenspace corresponding to µ and E h (µ) is the T h -invariant subspace of X h spanned by the eigenspaces corresponding to {µ i,h , i = 1 · · · , m} then lim h→0 δ(E(µ), E h (µ)) = 0.
Proof. We deduce from Theorem 4.1 that lim h→0 T − T h L(X h , X) = 0 and the result is a consequence of Section 2 of [11].
Proof. It is well-known that G is a projector in X with range E(µ) and we deduce from Theorem 4.2 that, for h sufficiently enough, G h is a projector in X h with range E h (µ), cf. [11]. Hence, for all for all x h ∈ E h (µ). We deduce from (40) that there exist constants C > 0 and h 0 > 0 such that, for all h < h 0 , On the other hand, we notice that Gx = x for all x ∈ E(µ). Then, for all y h ∈ X h and for h small enough, where the last inequality follows from the triangle inequality, (40), Lemma 4.1 and (41). It follows that Using (40) and the fact that T • P x = T x = µx for any x ∈ E(µ) yield .
The result follows now from Lemma 4.1 by noticing that both δ(T • P (X), X h ) and δ(T • P (X h ), X h ) are smaller than δ(T • P ( X), X h ).
Proof. Let u i,h be an eigenfunction corresponding to µ i,h such that x i,h = 1. We know from (42) that, if h is sufficiently small, .

SALIM MEDDAHI AND DAVID MORA
Then, there exists an eigenfunction x := ((σ, p), r) ∈ E(µ) satisfying , which proves that |x | is bounded from below and above by constant independent of h. Proceeding as in the proof of Lemma 3.1 we obtain that for all y h := ((τ h , q h ), s h ) ∈ X h . With the aid of (43), it is easy to show that the identity for h sufficiently small. The result follows now from the fact that A and B are continuous bilinear forms on X.

5.
Asymptotic error estimates. We begin this section by recalling some wellknown approximation properties of the finite element spaces introduced above. Given s ∈ (0, 1], let Π h : H s (Ω S ) n×n ∩ W → W h be the usual lowest-order Raviart-Thomas interpolation operator (see [9]), which is characterized by the identities for all faces (edges) F of elements T ∈ T h (Ω S ), with ν F being a unit vector normal to the face (edge) F . It is well known that Π h is a bounded linear operator and that the following commuting diagram property holds true (cf. [9]): where L h : L 2 (Ω S ) n → U h is the L 2 (Ω S ) n -orthogonal projector. In addition, it is well-known that the arguments leading to [14,Theorem 3.16] allow showing that there exists C > 0, independent of h, such that τ − Π h τ 0,ΩS ≤ Ch s τ s,ΩS + div τ 0,ΩS ∀τ ∈ H s (Ω S ) n×n ∩ H(div; Ω S ).
(45) Finally, we denote by R h : Q → Q h the orthogonal projector with respect to the L 2 (Ω S ) n×n -norm and by π h : H 1 (Ω F ) → V h the orthogonal projector with respect to the H 1 (Ω F )-norm. Then, for any s ∈ (0, 1], we have q − π h q 1,ΩF ≤ Ch s q 1+s,ΩF ∀q ∈ H 1+s (Ω F ).
The following estimate holds true.
Next, we introduce the discrete counterparts of E and E, defined for any q ∈ H 1 (Ω F ) by It is clear that E h q ∈ Y h for all q ∈ H 1 (Ω F ). Indeed, as ν is piecewise constant on Σ, where ̺ h : L 2 (Σ) n → P 0 (Σ h S ) n stands for the vectorial counterpart of ̺ h . Moreover, we have the following result.
Proof. On the one hand, we have that (T • P )( X) ⊂ T (X). On the other hand, it is straightforward that, for any x ∈ X, ∀y ∈ ker(a) × Q.
Thus, (I − T ) • P ( X) ⊂ P (X). It follows from Lemmas 2.3 and 2.4 that and the result is a consequence of Lemma 5.3 and the approximation properties (46)-(49).
Let E h be the operator defined in (50) and let Equations (51)-(52) constitute a conforming finite element discretization of the mixed problem (21)-(22) used to define P . The uniform discrete inf-sup condition of B for the pair {Y h,R , Q h × U h } is an easy consequence of (24). Moreover, [17, Lemma 2.1] guarantees the uniform ellipticity of d on W × H 1 (Ω F ) ⊃ ker h (a), whereas the fact that div(W h ) ⊂ U h implies that ker h (B) ⊂ ker h (a). Hence, as a consequence of the Babuška-Brezzi theory, problem (51)-(52), is well posed. Furthermore, thanks to the definition of E hp , the first estimate from Lemma 5.1, and the fact that π hp 1,ΩF ≤ p 1,ΩF (since π h is a projection), we can claim that the operators P h are bounded uniformly with respect to h and the following Strang-like estimate holds true: where (( σ 0 , c), ( u, r)) and (( σ h,0 , c h ), ( u h , r h )) are the solutions to (21)-(22) and (51)-(52), respectively. As a consequence, we have the following estimate.
Lemma 5.5. There exists C > 0, independent of h, such that Proof. See [17, Lemma 6.3] for more details.
Lemma 5.6. There exists C 1 > 0, independent of h, such that Moreover, if µ / ∈ {0, 1} is an eigenvalue of T . Then, there exists a constant C 2 > 0, independent of h such that Proof. We first notice that by definition which proves the first estimate of the Lemma.
On the other hand, if x := ((σ, p), r) ∈ E(µ) then T • P (x) = µx and we deduce again from (23) that there exists a constant C > 0 such that where u is the displacement field given by u = µ (µ−1)ρS div σ ∈ H 1+tS (Ω S ) n . With this regularity result at hand, we can proceed as in (54) to obtain and the second estimate of the Lemma follows.
We conclude that we have have the following asymptotic convergence for the eigenfunctions and eigenvalues of problem (1)-(7).
6. Numerical results. We use a two-dimensional benchmark test that is identical to the one carried out in [17]. The geometrical data representing an elastic container (steel) filled with a compressible liquid (water) is shown in Figure 2. The physical parameters are given by:   We use several meshes which are successive uniform refinements of the coarse initial triangulation shown in Figure 2. The refinement parameter N is the number of element layers across the thickness of the solid (N = 1 for the mesh in Figure 2).
We can distinguish between two types of vibrations corresponding to sloshing and elastoacoustic modes. We refer to [4,6] for a more detailed discussion on sloshing (or gravity) and elastoacoustic frequencies. We report the lowest computed sloshing vibration frequencies ω S h,k in Table 1 and the elastoacoustic vibration frequencies ω E h,k in Table 2. The tables also include the estimated orders of convergence, as well as more accurate values of the vibration frequencies extrapolated from the computed ones by means of a least-squares fitting. A double order of convergence can be clearly observed in all cases. We finally notice that our results are in agreement with those obtained in [17] and based on the Arnold-Falk-Winther element [2]. This happens even though the computer cost of the lowest-order PEERS element is lower then that of the Arnold-Falk-Winther (AFW) element. Indeed, the global number of unknowns for PEERS elements is ≈ 12N v while which the number of unknowns for AFW elements is ≈ 20N v , where N v represents the number of vertices in T h (Ω S ).