Analysis and numeric of mixed approach for frictional contact problem in electro-elasticity

Abstract: This work handle a mathematical model describing the process of contact between a piezoelectric body and rigid foundation. The behavior of the material is modeled with a electro-elastic constitutive law. The contact is formulated by Signorini conditions and Coulomb friction. A new decoupled mixed variational formulation is stated. Existence and uniqueness of the solution are proved using elements of the saddle point theory and a fixed point technique. To show the efficiency of our approach, we present a decomposition iterative method and its convergence is proved and some numerical tests are presented.


Introduction
T he numerical study of the piezoelectric [1,2] contact problems [3][4][5][6][7] presents a great challenges, because the non coercivity and non differentiability of some terms. The linear term coupling the mechanical field and electric potential is non coercive and non symmetric. The term corresponding to the friction is convex and non differentiable, in the variational formulation in one hand. In the other hand, the non-linear coupling mechanical field and frictional contact, it can be shown in the norm of the tangential component of the mechanical field present in the frictional function.
To overcome these difficulties, the authors develop some methods like finite elements method [8,9], penalty method and fixed point method [4]. The most and sufficient method for this type of problem is the ones based on convex duality [10] and the introduction of Lagrange multipliers [11][12][13]. Primal-dual active sets strategy, which is an equivalent to infinite dimensional semismooth Newton method, is applied in [14] while in [15] the author propose a numerical approximation and based Uzawa block relaxation method. Alternating directions method of multipliers (ADMM) is based in [16].
In mechanic of structures one is always interested in the determination of the stress tensor σ more than the mechanical displacement itself. Methods have been developed for calculate an approximation of σ from u and the drawback is not easy to build an approximation space of tensors satisfying the equilibrium relations and required regularity. A mixed variational formulation have been developed to handle this difficulty. Concerning the piezoelectric contact problems, mixed formulation were developed in [17][18][19].
In this paper, we introduce a mixed variational approach based Lagrange multipliers which describes the static frictional contact between a piezoelectric body and non-conductive foundation. The standard mixed variational formulation of contact problems is formally in the following form [17,18,20]: where a(·, ·) is symmetric, coercive and the term b(·, ·) coupling the normal and tangential Lagrange multipliers. This coupling is difficult to handle numerically and for model describing the contact problem with friction it is important to use the decoupling form to identify the slid and slip on the contact zone. Moreover, when the model (1) describe problem with electroelastic body the bi-linear form a(·, ·) become not symmetric. The idea consists in decoupling the contact from the friction by introducing two convex sets, one is reserved to the contact multiplier and the second is stated for the friction. This approach leads to decoupled inequalities in mixed variational problem. Since a(·, ·) is not symmetric we follow the standards techniques and steps based on [17,18] but with different fixed point map (Step 3 in the proof) and hence more analysis is needed to get the existence and uniqueness result. The resulting problem is an system by blocks of two unknowns (displacement and potential), regarding its structure and the form of blocks we use Gauss elimination technique and then we obtain a Schur complement in the matrix corresponding to displacement subproblem. It is well known (see [16,21]) that this technique allows a suitable preconditioner for conjugate gradient method employed to solve the resulting symmetric and positive definite system.
To prove the efficiency of this approach, we state suitable numerical fixed point scheme. The convergence is proved basing abstract perturbed problem and fixed point process, Banach fixed point is hence applied. For details concerning the mathematical tools we refer to [22][23][24][25].
The paper is structured as: In Section 2, we present the model of equilibrium process of the elastic piezoelectric body in frictional contact with a non-conductive foundation,we introduce the functional spaces for various quantities, list the assumptions on given data and derive the weak formulation of the problem. In Section, 2.2 we state and prove our main existence and uniqueness result, Theorem 1. The proof of these theorem are carried out in several steps and are based on an abstract result in the study of elliptic variational inequalities and Banach fixed point technique. The successive iterative method is detailed followed by the convergence result in Section 4. In Section 5, we conclude with finite element discretization and we give some numerical experiments by simple example .

Problem setting
The piezoelectric body occupies in its reference configuration (initial configuration) the domain Ω ⊂ R d , d = 2, 3. We suppose that Ω is bounded with a smooth (enough) boundary ∂Ω = Γ. We denote by n be the outer normal to Γ and summation over repeated indices is implied and the index that follows a comma represents the partial derivative with respect to the corresponding component of the variable. The indices take values in {1, · · · , d} the summation convention over repeated indices is used.
Below we use S d to denote the space of second order symmetric tensors on R d while "." and · will denote the inner product and the Euclidean norm on S d and R d , that is We also use the notations u n and u τ for the normal and tangential displacement, that is u n = u.n and u τ = u − u n n. Similarly we denote by σ n and σ τ the normal and tangential stress tensor given by σ n = σn.n, σ τ = σn − σ n n. We introduce the following functional spaces on Ω; We recall the well known Green's formula where Div σ = (σ ij,j ) and for more details of this formula see [24]. In addition we shall use the following notations. u is the displacement field , ε(u) = (ε ij (u)) to denote the strain tensor, given by ε ij (u) = 1 2 (u i,j + u j,i ) and σ = (σ ij ) being the stress tensor. Lets ϕ denote the electric potential, E(ϕ) = (E i (ϕ)) is the electric field, which is defined by E i (ϕ) = −ϕ ,i and D = (D i ) is the electric displacement field.
The equilibrium equations are given by where the constitutive relations for the piezoelectric material are: where A = (a ijkl ) is a (fourth-order) elasticity tensor, B = (b ijk ) is the (third-order) piezoelectric tensor, B * is the transpose of B and β = (β ij ) is the electric permitivity and div(D) = D i,i (see [26]). To give the mechanical and electrical boundary conditions, we subdivide Γ into three disjoints measurable parts Γ 1 , Γ 2 , Γ 3 such that meas(Γ 1 ) > 0. The body is assumed to be clamped on Γ 1 and surfaces traction of density f 2 act on Γ 2 , on Γ 3 the body can reaches a frictional contact with the so called foundation (insulating foundation). A second partition of Γ, that is Γ = Γ 3 ∪ Γ a ∪ Γ b . Surface electric charge of density q 2 acts on Γ b , and the electric potential vanishes on Γ b . We use the same symbol v for the trace of v on Γ.
The contact and the Coulomb friction conditions: where F is the friction coefficient and g is the gap between the body and the rigid foundation. The electric contact condition is; To resume, we consider the following problem: Problem 1. Find the displacement field u : Ω −→ R d and the electric potential field ϕ : Ω −→ R such that (2)-(12) hold.
To study of Problem 1 we will assume, under Einstein summation convention, that: Let us introduce the following Hilbert spaces: If u and ϕ are regular functions which satisfy (2)-(10), then we find: Let us introduce the functional spaceṼ = V × W, which is the Hilbert space endowed with the inner product: Moreover, by Riesz's representation theorem, we definef ∈Ṽ by: Using the previous tools, we find: where ·, · H * Γ ,H Γ denotes the duality product between H * Γ and H Γ , K n the set of the normal component of admissible displacement, i.e. K n = {v n , v n ≤ g} .
It is straightforward that M T,N are two closed convex sets of H * Γ and 0 H * Γ ∈ M N,T . We introduce two dual Lagrange multipliers λ N and λ T ∈ M as follows: We define two bi-linear and continuous forms b 1 and b 2 for all v ∈ V, δ 1 , δ 2 ∈ H * Γ as follows: σ n u n dΓ, and by definition of M N we have: Also, taking into account the definition of M T , λ T , and the assumption (10), we have: Keeping in mind that the Sobolev trace operator is linear and continuous, it is clear that there exists In addition, using the properties of the Sobolev trace operator it can be shown that there exists α i > 0 such that: The following weak formulation of Problem 1 is then obtained :

Main results
In this section we present our main results.
where b(·, ·) = b 1 (·, ·) + b 2 (·, ·) :Ṽ × M N × M T −→ R. Now we introduce a numerical scheme to get numerically the solution of the Problem 2. The scheme is an fixed point iterative and is stated in the following Algorithm 1.

Algorithm 1 Decomposition method for Problem 2
Initialization: λ N and λ T Iteration: ≥ 1. Compute successivelyũ +1 , λ +1 N and λ +1 T using Algorithm Step 1: Findũ +1 ∈Ṽ such that Step 2: Step 3: Proposition 1. Let (u , λ ) be the solution generated by the Algorithm 1, then The proof of the main results will be presented in the next section.

Proof of the main result
Let X and Y two Hilbert spaces endowed with the inner product (·, ·) X and (·, ·) Y respectively and let us consider two bi-linear forms as follows: a(·, ·) : X × X −→ R, generally non symmetric, such that there exists α > 0 such that Now, let M = M N × M T ⊂ Y × Y be closed and convex set that contain 0 Y×Y , we consider the following problem: Problem 3. For given f ∈ X, find u ∈ X and λ = (λ N , λ T ) ∈ M such that: We have the following result; Theorem 2. Let f ∈ X and assume that (30)-(34) hold. Then, there exists a unique solution (u, λ) of Problem AWF Moreover, if (u 1 , λ) and (u 2 , γ) are two solutions of the Problem AWF for given data functions f 1 ∈ X and f 2 ∈ X respectively, then, there exists SC > 0 such that: Proof. We consider the symmetric a 0 (., .) and anti-symmetric c(., .) part of a(., .) respectively, defined by For given 0 < t < 1, let us present the following bi-linear form For all t ∈ [0, 1], we note that Let us consider the following auxiliary perturbed problem: Problem 4. (Auxiliary perturbed problem) For given f ∈ X, find u ∈ X and λ ∈ M, such that The rest of the proof will be treated by steps. Step 1. If t = 0, the Problem 4 has unique solution. Indeed, if l = 0 the problem is equivalent to the saddle point problem: find u ∈ X and λ ∈ M such that L(., .) has at least one solution, see [10], in fact, from L(v, 0) = 1 2 a 0 (v, v) − ( f , v) X and the coercivity of a 0 (., .) we have lim Indeed, let δ 0 be an element of M and let u δ 0 ∈ X be the solution of the equation which is equivalent that u δ 0 is the solution of the following minimization problem Substituting v = u δ 0 in (44), we get Additionally, using the inf-sup property of the form b(., .); we deduce that there exists a constant C > 0 such that From (45) we deduce (38), which implies the existence of solution of Problem 4.
To show the uniqueness of the solution, let us assumes that (u 1 , λ) and (u 2 , γ) are two solutions of the problem By subtracting these two equations, we find and by coercivity of a 0 , we have u 1 = u 2 . Moreover and by inf-sup property of b(., .), we have and finally λ = γ.
Step 2. Assume now that f ∈ X, there exists a unique solution (u, λ) ∈ X × M of the Problem 4, when we have two solutions (u 1 , λ) and (u 2 , γ) of Problem 4 corresponding to two given data f 1 ∈ X × X and f 2 ∈ X × X respectively, then In fact Since a t is coercive and, hence In addition, b(v, λ − γ) = ( f 1 − f 2 , v) X + a r (u 2 − u 1 , v) and by inf-sup property of b(., .) we have Hence (47) and (48) lead to Step 3. Let τ ∈ [0, 1]. Assume for given f , g ∈ X there exists a unique solution of Problem 4 with t = τ, (u, λ) ∈ X × M. Then for given f ∈ X there exists a unique solution (u, Indeed, given f ∈ X, we define the operator T : is the solution of the following problem: Problem 5. For given f ∈ X, find u ∈ X and λ ∈ M, such that We will show that T is a contraction. To this end, we consider two pairs (w 1 , ξ) and (w 2 , By the same argument in (48) and by definition of F s In addition, and hence, which implies that T is a contraction and by Banach fixed theorem we conclude that T has unique fixed point. Let (ū,λ) be the unique fixed point of T, using the definition of the operator T, we deduce that We substitute F s in the first equation, we find (ū,λ) be the unique fixed point of T, using the definition of the operator T, we deduce that which gives the existence of solution. In order to justify the uniqueness, let us assume that the problem with l = s ∈ [τ, τ + t 0 ] has two solutions (u 1 , λ) and (u 2 , γ), we have hence, by coercivity of a s , we get u 1 = u 2 and λ = γ.
Step 4. Using Step 3, a finite number of times, we deduce that the Problem 4 admits a unique solution (u, λ) for t = 1.
Step 5. In order to get (44), let us consider the data By subtracting this two equations, we find and by (34), we have Using (56), we can write where c 1 , c 2 are strictly positive constants. By combining this inequality and (57), we deduce that The constants c 1 and c 2 are chosen such that m a − c Finally, combining (57) and (60), we have (47).
Proof of Theorem 1. We consider X =Ṽ, Y = H * Γ and M N × M N given by (18). The subset M N × M N is a non-empty, closed, convex of H * Γ × H * Γ and 0 H * Γ ∈ M. By using (12) and (15) we deduce that there exists M a = M a (A, E , β) > 0 and m a = m a (A, β) > 0 such that the bilinear form a(., .) satisfies By the conditions (21) and (22) we deduce that the bilinear form b(., .) satisfies (32). Using inf-sup property (34) and Theorem 2 we find the result of Theorem 1.

Proof of Proposition 1.
To prove the convergence result 1 of the Algorithm 1, lets reconsider the following perturbed problem, for l ∈ [τ; τ + t 0 ] ⊂ [0, 1], where this time: Given f , g ∈ X, we define the mapping T : is the solution of the following fixed point problem: Problem 6. (Fixed point problem) For given f , g ∈ X, find u ∈ X and λ ∈ M, such that It is straightforward that By the conditions (61)-(62) and (66), the operator T is a contraction. This implies that there exists (ū,λ) such that T(ū,λ) = (ū,λ). Hence we have T(u , λ ) = (u +1 , λ +1 ) and the following scheme converges; Problem 7. For given f ∈ X, find u +1 ∈ X and λ +1 ∈ M, such that The convergence of the iterative fixed point scheme in the Algorithm 1 follows directly.

Discretization and numerics
The problem is now how to identify the multipliers λ N and λ T in the convex sets M 1,2 . One manner to do this, is the use of Projection maps. To this end, we consider the finite dimensional spaces V h ⊂ V, K h = K ∩ V h and W h ⊂ W approximating the spaces V and W, respectively, in which h > 0 denotes the spatial discretization parameter. Let us define: Let us denote also X * h n ⊂ X * n L 2 (Γ 3 ) and X * h T ⊂ X * T L 2 (Γ 3 ; R d−1 ) the finite discretizations of X * n and X * T respectively, such that the following discrete Babuska-Brezzi inf-sup conditions hold; with α independent of h. We consider the following discrete approximation of Problem 2: where r > 0 and and where P M is the projection over M. For more details we refer to [27].

Matrix formulation
In this section, we adopt the same technical discretization as in the work [27]. Let a j (j = 1, . . . , n c ) be a contact node (i.e. a j are the nodes forming Γ 3 ). The displacement vector at a j is denoted by u j , i.e. u j = u(a j ). Denoting n j and t j , the unit outward normal vector to Γ 3 and the unit tangential vector to Γ 3 , respectively. Let us introduce the linear mappings • N : R 2d → R n c , such that (Nu) j = u j n j , j = 1, . . . , n c .
The finite element discretization leads to the following matrices and vectors: • A, (2d) × (2d) the elastic matrix (symmetric and positive definite) ; • B, d × d electric potential stiffness matrix (symmetric positive definite); • E, d × (2d) coupling matrix; • M n and M τ normal and tangential mass matrices (n c × n c ); • f (the external forces in R 2d ), q (the external charges in R d ).
• λ N , λ T the vectors associated to λ h N and λ h T respectively.
With the above notations, we can solve the Problem 8 with Coulomb friction using fixed point procedure on the friction threshold (see [27,28] for more details). The Fixed point on the friction threshold procedure is in Uzawa Algorithm 2, Algorithm 2 Uzawa algorithm (UA) for Problem 8 Initialization: λ N and λ T Iteration: ≥ 1. Compute successively U +1 , λ +1 N and λ +1 T as follow Step 1: U +1 such that AU +1 = b (70) Step 2: λ +1 N such that λ +1 Step 3: λ +1 T such that B(0, −F λ +1 N ) denote the ball of center 0 and radius −F λ N > 0 and x + denote the non negative part of x i.e. x + = max(0, x).

Remark 1.
In practice, the Algorithm 2 is solved using Tresca friction with slip bound S and an fixed point iteration is used to compute the problem with Coulomb friction, i.e. S +1 = −F λ N .
To compute the solution of the system (70), we proceed by the following elimination technique: where S c is the Schur complement given by S c = A − E B −1 E. Since the matrices A and B −1 are symmetric and positive definite, the suitable method for solving the subsystem (74) is the Conjugate Gradient method (CG). We take profit from the Schur complement S c to obtain a convenient preconditioner, as discussed in [16] and reference therein the (CG)-preconditioner is P = A. The preconditioned Conjugate Gradient method for solving the system (74) is stated in the Algorithm 3. Once u is computed one can compute ϕ by the explicit formula (73). The solution of the problem with Coulomb friction is stated in the algorithm 4
Step 3: Update the slip bound by S +1 = −F λ N . displacements and the electric potential are prescribed, i.e., u = 0 and ϕ = 0 on Γ D . On Γ N , non-homogeneous Neumann boundary conditions are prescribed f 0 = σ(u).n = −2. On {1} × (0, 1) the homogeneous Neumann boundary condition is applied (σn = 0 and Dn = 0). For seek of simplicity, the normalized gap between Γ 3 and the foundation is g(x) = 0 and the friction coefficient is F = 0.6. The mesh is generated by using Matlab function "kmg.m" built by the author in [29].
The deformed configuration is showed in the Figure 2 and the the contour plot of the electric potential distribution are showed in the Figure 4. The Figure 3 show the Lagrange multipliers λ T,N in the contact zone Γ 3 corresponding to the problem with Coulomb friction condition. The Figures 5 and 6 show the Lagrange multipliers λ T,N with different choices of loads acting on Γ 2 , it is clear that the slide occurs when the load is important (large enough). The performance of the algorithm is presented in the Table 1, as showed the number of iterations is independent of mesh refinement and the time of execution is significant.

Conclusion
We have investigated numerical analysis of a model describing the process of contact between a piezoelectric body and non-conductive foundation. The behavior of the material is modeled with a electro-elastic constitutive law. The contact is formulated by Signorini conditions and Coulomb friction. In coming works, more general problem with non linear constitutive equation and non monotone friction will be treated. This is may be handled by fixed point iteration and hemivariational inequalities [30].