Jérémie Lasry
email: j_lasry@insa-toulouse.fr

Julien Pommier

Yves Renard

Michel Salaün

eXtended finite element methods for thin cracked plates with Kirchhoff-Love theory

Keywords: structures, fracture, eXtended finite element method, plates, Kirchhoff-Love, rate of convergence

A modelization of cracked plates under bending loads in the XFEM framework is addressed. The Kirchhoff-Love model is considered. It is well suited for very thin plates commonly used for instance in aircraft structures. Reduced HCT and FVS elements are used for the numerical discretization. Two kinds of strategies are proposed for the enrichment around the crack tip with, for both of them, an enrichment area of fixed size (i.e. independant of the mesh size parameter). In the first one, each degree of freedom inside this area is enriched with the nonsmooth functions that describe the asymptotic displacement near the crack tip. The second strategy consists in introducing these functions in the finite element basis with a single degree of freedom for each one. An integral matching is then used in order to ensure the C 1 continuity of the solution at the interface between the enriched and the non-enriched areas. Finally, numerical convergence results for these strategies are presented and discussed.

INTRODUCTION

This paper deals with an adaptation of XFEM (eXtended Finite Element Method) to the computation of thin plates, which present a through-thickness crack. Let us recall that XFEM is a way to introduce the discontinuity across the crack and the asymptotic displacement into the finite element space. It has been initially developed for plane elasticity problems (see [START_REF] Moës | A finite element method for crack growth without remeshing[END_REF][START_REF] Moës | X-fem: Nouvelles frontières pour les eléments finis[END_REF]) and is now the subject of a wide literature (see for instance [START_REF] Moës | Non-planar 3D crack growth by the extended finite element and level sets, Part I: mechanical model[END_REF][START_REF] Laborde | High order extended finite element method for cracked domains[END_REF][START_REF] Béchet | Improved implementation and robustness study of the X-FEM for stress analysis around cracks[END_REF][START_REF] Sukumar | Partition of unity enrichment for bimaterial interface cracks[END_REF][START_REF] Areias | Analysis of three-dimensional crack initiation and propagation using the extended finite element method[END_REF][START_REF] Bordas | Enriched finite elements and level sets for damage tolerance assessment of complex structures[END_REF]).

Nowadays, as far as we know, there are few previous works devoted to the adaptation of XFEM to plate or shell models [START_REF] Areias | Analysis of fracture in thin shells by overlapping paired elements[END_REF][START_REF] Areias | Analysis of finite strain anisotropic elastoplastic fracture in thin plates and shells[END_REF][START_REF] Wyart | A sustructured FE-shell/XFE-3D method for crack analysis in thin-walled structures[END_REF][START_REF] Dolbow | Modeling fracture in Mindlin-Reissner plates with the extended finite element method[END_REF]. In [START_REF] Areias | Analysis of fracture in thin shells by overlapping paired elements[END_REF][START_REF] Areias | Analysis of finite strain anisotropic elastoplastic fracture in thin plates and shells[END_REF], shell models are used: since the near tip asymptotic displacement in this model is unknown, no singular enrichment is formulated in these references. In [START_REF] Wyart | A sustructured FE-shell/XFE-3D method for crack analysis in thin-walled structures[END_REF], which deals with cracked shells, the cracked part of the domain is modelized by a three-dimensional XFEM formulation. It is matched with the rest of the domain, formulated with a classical finite element shell model. However, this is not the way we choose to work. In this paper, a plate model is kept on the entire domain, and we consider singular enrichment.

In [START_REF] Dolbow | Modeling fracture in Mindlin-Reissner plates with the extended finite element method[END_REF], the plate model used is the Mindlin-Reissner one. However, in this reference, an important locking effect for thin plates has been detected despite the use of some classical lockingfree elements. This suggests that this locking effect is due to the XFEM enrichment. In order to avoid locking, some special treatments are mandatory, for example selective reduced integration in the case of the QUAD 4 element (see [START_REF] Mcneal | A simple quadrilateral shell element[END_REF]). But these treatments are not straightforwardly adaptable to an enrichment by nonsmooth functions, which are not polynomials.

The most popular models for the small deformations of plates are the so-called Mindlin-Reissner and Kirchhoff-Love ones. The Mindlin-Reissner model is convenient for moderately thin plates. It takes into account the transverse shear strain but finite element approximations are subject to the so-called shear locking phenomenon for very thin plates. The Kirchhoff-Love model, which is not subject to the shear locking phenomenon, provides a realistic description of the displacement of thin plates, and especially very thin ones, since it is the limit model when the thickness vanishes (see [START_REF] Ciarlet | A justification of two-dimensional linear plate model[END_REF]). It has already been used for the purpose of fracture mechanics (for instance, see [START_REF] Hui | A theory for the fracture of thin plates subjected to bending and twisting moments[END_REF]). Moreover, for through-thickness cracks, the limit of the energy release rate of the three-dimensional model can be expressed with the Kirchhoff-Love model solution (see [START_REF] Ph | Une Théorie Asymptotique des Plaques Minces en Élasticité Linéaire[END_REF][START_REF] Ph | Sur une interprétation mathématique de l'intégrale de Rice en théorie de la rupture fragile[END_REF]).

The Kirchhoff-Love model corresponds to a fourth-order partial differential equation. Consequently, a conformal finite element method needs the use of C 1 (continuously differentiable) finite elements. To avoid too costly elements we consider the reduced Hsieh-Clough-Tocher triangle (reduced HCT) and Fraeijs de Veubeke-Sanders quadrilateral (reduced FVS) (see [START_REF] Ciarlet | The Finite Element Method for Elliptic Problems[END_REF]). They lead to a satisfactory theoretical accuracy, with a reasonable computational cost. Let us also remark that, for this plate model, the crack tip asymptotic bending displacement is well-known for an isotropic plate: it corresponds to the bilaplacian problem one. Now, concerning the specific XFEM enrichment, a 'jump function' (or Heaviside function) is used in order to represent the discontinuity due to the crack. Following the ideas already presented in [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF], we propose two strategies for the crack tip enrichment. In both of them, an enrichment area of fixed size is defined, centered on the crack tip. In the first strategy, each node contained in the enrichment area has all its degrees of freedom (dof) enriched with the crack tip asymptotic bending. For the second strategy, the asymptotic bending displacements are introduced in a global way. Then a matching condition is needed in order to ensure the continuity of the displacement and its derivatives across the interface between the enrichment area and the remaining part of the domain.

Along this paper, it will be underlined that our XFEM formulation presents a sense of optimality, both about accuracy and computational cost, since they are nearly equal to those of a classical finite element method on a regular non-cracked problem. This paper is organized as follows. Section 2 describes the model problem. Section 3 is devoted to some aspects of the finite element discretization of the Kirchhoff-Love model. In Section 4, the enrichment strategies are presented which are evaluated on two test problems in Section 5.

THE MODEL PROBLEM

Notations and variational formulation

Let us consider a thin plate, i.e. a plane structure for which one dimension, called the thickness, is very small compared with the others. For this kind of structures, starting from a priori hypotheses on the expression of the displacement fields, a two-dimensional problem is usually derived from the three-dimensional elasticity formulation by means of integration along the thickness. Then, the unknown variables are set down on the mid-plane of the plate. In all the following, this mid-plane will be denoted by . It is an open subset of R 2 . So, in a three-dimensional cartesian referential, the plate is the set

{(x 1 , x 2 , x 3 ) ∈ R 3 , (x 1 , x 2 ) ∈ and x 3 ∈]-ε ε[}.
The x 3 coordinate corresponds to the transverse direction. All the mid-plane points have their third coordinate equal to 0 and the thickness is 2ε (see Figure 1). Finally, we assume that the plate has a through-thickness crack (see Figure 1) and that the material is isotropic of Young's modulus E and the Poisson ratio .

In plate theory, it is usual to consider the following approximation of the three-dimensional displacements for (

x 1 , x 2 , x 3 ) ∈ ×]-ε, ε[ u 1 (x 1 , x 2 , x 3 ) = u 1 (x 1 , x 2 )+ x 3 1 (x 1 , x 2 ), u 2 (x 1 , x 2 , x 3 ) = u 2 (x 1 , x 2 )+ x 3 2 (x 1 , x 2 ), u 3 (x 1 , x 2 , x 3 ) = u 3 (x 1 , x 2 ). (1)
In these expressions, u 1 and u 2 are the membrane displacements of the mid-plane points while u 3 is the deflection, 1 and 2 are the section rotations. In the case of an homogeneous isotropic material, the variational plate model splits into two independent problems: the first, called the membrane problem, deals only with the membrane displacements, while the second, called the bending problem, concerns the deflection and the rotations. The membrane problem corresponds to the classical plane elasticity problem and has been already treated in many references (see for instance [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF][START_REF] Béchet | Improved implementation and robustness study of the X-FEM for stress analysis around cracks[END_REF]). In this paper, we shall only address the bending problem.

For reasons mentioned in the introduction, we consider the Kirchhoff-Love model, which can be seen as a particular case of (1), as it is obtained by introducing the Kirchhoff-Love assumptions, which read as

= -∇u 3 ⇔ 1 = -* 1 u 3 , 2 = -* 2 u 3 , (2)
where, the notation * stands for the partial derivative with respect to x (for = 1, 2). A first consequence of this relation is that the transverse shear strain is identically zero, which avoids the shear locking problem. A second consequence of ( 2) is that the section rotation only depends on the transverse displacement. It means that this displacement is the only unknown function for the bending problem. For convenience, it will be denoted by u all along the following of this paper. So, in the Kirchhoff-Love framework and for a pure bending problem, the three-dimensional displacement reads as

u 1 (x 1 , x 2 , x 3 ) = -x 3 * 1 u(x 1 , x 2 ), u 2 (x 1 , x 2 , x 3 ) = -x 3 * 2 u(x 1 , x 2 ), u 3 (x 1 , x 2 , x 3 ) = u(x 1 , x 2 ).
For the sake of simplicity, we assume the plate to be clamped on its boundary, except along the crack faces where a traction-free condition is considered. The plate is subject to a volume force f of coordinates ( f 1 , f 2 , f 3 ), and two surface forces, say g + and g -, applied on the top and bottom surfaces. The variational formulation (or virtual work formulation) of the Kirchhoff-Love model reads as

     Find u ∈ H 2 0 ( ) such that for any v ∈ H 2 0 ( ) 2Eε 3 3(1-2 ) [(1-)* 2 u + u ]* 2 v dx = [Fv + M * v] dx. (3)
where

F = g + 3 + g - 3 +
ε -ε f 3 dx 3 is the resulting transverse loading and M = ε(g + -g -)+ ε -ε x 3 f dx 3 is the resulting moment loading. Moreover, stands for the Kronecker's symbol and the summation convention over repeated indices is adopted, Greek indices varying in {1, 2}. Finally, H 2 0 ( ) is the classical Sobolev space, the definition of which can be found in [START_REF] Adams | Sobolev Spaces[END_REF] for instance.

Asymptotic displacement near the crack tip and fracture modes

In the Kirchhoff-Love plate model, there are two fracture modes. Applying a symmetric bending leads to the first fracture mode, while applying an anti-symmetric bending or a transverse shear, leads to the second one (see Figure 2). Left: a symmetric bending leads to mode I. Right: a shear bending leads to mode II.

To characterize them, let us recall that the governing equation related to the bending variational problem (3) reads

2Eε 3 3 (1-2 ) 2 u = F +* M (4)
on the mid-plane . It is a bilaplacian problem for which the asymptotic behavior near the crack-tip is well-known (see [START_REF] Grisvard | Singularities in Boundary Value Problems[END_REF]). So, close to the crack tip, the displacement may be written as u = u r +u s , where u r stands for the regular part of the transverse displacement and belongs to H 3 ( ). The nonsmooth part u s belongs to H 5/2-( ), for any >0, and reads

u s (r, ) = Ar 3/2 K 1 +7 3( -1) cos 3 2 +cos 2 + K 2 3 +5 3( -1) sin 3 2 +sin 2 (5)
in polar coordinates relatively to the crack tip and with A = √ 2

Remark 1

The singular modes of the Kirchhoff-Love theory, which are presented in Figure 2, are different from the ones of the well-known three-dimensional elasticity theory. Some interesting details on this subject are given in [START_REF] Zucchini | Crack tip stress fields for thin, cracked plates in bending, shear and twisting: a comparison of plate theory and three-dimensional elasticity theory solutions[END_REF].

FINITE ELEMENT APPROXIMATION OF THE KIRCHHOFF-LOVE MODEL

Choice of the finite element discretization

Let us introduce now the finite element discretization of the variational formulation (3). As usual, a finite-dimensional approximation space V h of the continuous space of solutions H 2 0 ( ) is built. In order to have a conformal method, i.e. V h ⊂ H 2 0 ( ), the functions of V h must be continuously differentiable, which needs the use of C 1 finite elements. Among the available elements having this regularity (see [START_REF] Ciarlet | The Finite Element Method for Elliptic Problems[END_REF]), the reduced HCT triangles and FVS quadrangles are of particular interest. For the HCT (resp. FVS) element, the triangle (resp. quadrangle) is divided into three (resp. four) sub-triangles (see Figure 3). The basis functions of these elements are P 3 polynomials on each sub-triangle and matched C 1 across each internal edge. In addition, to decrease the number of degrees of freedom (dof), the normal derivative is assumed to vary linearly along the external edges of the elements (this assumption does not hold on the internal edges). Finally, both for triangles and quadrangles, there are only three dof on each node: the value of the function and its first derivatives. So, these elements have the two following properties:

1. With the reduction, the computational cost is limited to three dof for each node of the mesh, like a classical Mindlin-Reissner element for which each node has also three dof (the deflection and the two section rotations). 2. For regular problems (see [START_REF] Ciarlet | Basic error estimates for elliptic problems[END_REF]), the theoretical error is in O(h) for the energy norm and O(h 2 ) for the displacement norm, where h stands for the mesh parameter.

Then, the reduced HCT or FVS elements and the standard Mindlin elements have the same features as far as numerical cost and accuracy are concerned.

Remark 2

One can also notice that reduced HCT/FVS elements have the same dof. Consequently, it is possible to use meshes made of triangles and quadrangles at the same time.

The shape functions of the reduced HCT and FVS are presented in Figures 4 and5 one corresponding to the x 2 -derivative dof can be obtained by symmetry). Note that the expression of the shape functions is rather complex because it is a polynomial of degree three in each subtriangle. The implementation used in Getfem++ [START_REF] Pommier | Getfem++, An Open Source Generic C ++ Library for Finite Element Methods[END_REF] uses the expression of the shape functions of the (not reduced) HCT and FVS elements on the reference element (complete expressions are freely available in the Getfem++ sources). The reduction operation is performed during the transformation from the reference element to the real one.

Integral matching

As mentioned in the introduction, one of the strategies, which will be presented further, needs to match the enriched subdomain with the non-enriched one. So we shall have to carry out calculations in the case of a domain splitted into two subdomains, including an integral matching at the subdomains boundary.

Some theoretical results exist on the convergence of the mortar methods for fourth-order problems. For instance, in [START_REF] Marcinkowski | Mortar methods for some second and fourth order elliptic equations[END_REF][START_REF] Marcinkowski | A mortar finite element method for fourth order problems in two dimensions with lagrange multipliers[END_REF], it is shown that the use of mortar methods with reduced HCT/FVS elements can achieve an error in O(h)/O(h 2 ) for the energy/displacement norms. Nevertheless, as far as we know, no numerical experiment checking the accuracy of such a method has been published. Moreover, we do not use exactly the method described in [START_REF] Marcinkowski | Mortar methods for some second and fourth order elliptic equations[END_REF][START_REF] Marcinkowski | A mortar finite element method for fourth order problems in two dimensions with lagrange multipliers[END_REF]. So, before introducing our enrichment strategy on a cracked domain, it appeared meaningful to test our integral matching conditions first.

Let us consider the problem (4) without crack. The loading F is such that the exact solution is u(x 1 , x 2 ) = sin(10(x 1 + x 2 )). Let us observe that it also implies non-homogeneous Dirichlet boundary conditions. In a first step, we compute the solution to this problem on structured and non-structured meshes with the reduced HCT/FVS elements. This enables to check that the errors in energy and displacement norms are in accordance with what is expected. It also establishes a reference error level which will be compared with the mortar method one.

Then, we consider that the domain is splitted into two subdomains, say 1 and 2 , with matching triangulations. It means that, on both sides of the boundary between the two subdomains, the nodes are located at the same place, see Figure 6. It is identical to the so-called 'mortar method with matching triangulations', in [START_REF] Marcinkowski | Mortar methods for some second and fourth order elliptic equations[END_REF]. The matching conditions, which involve the values of the unknowns, say u i on each part i are introduced in order to ensure the continuity of the function and its normal derivative. The following relations were chosen at this aim

u 1 d = u 2 d ∀ ∈ , * n u 1 d = -* n u 2 d ∀ ∈ M, (6)
where and M are appropriate multiplier spaces. Let us notice that the change of sign in front of the normal derivative * n u i (second equation of ( 6)) is due to the fact that the outside normal vector has an opposite sign whether it refers to 1 or 2 . In our numerical experiments, the multipliers spaces and M are spaces of piecewise polynomials of degree 2 and 1, respectively. The numerical results for the mortar method, on structured, non-structured, triangular and quadrangular meshes, are given in Figure 7, where they are compared with the classical method. It shows that our integral matching does not affect the accuracy: the convergence curves are nearly identical in all cases.

XFEM ENRICHMENTS

Basics

As in previous papers [START_REF] Moës | A finite element method for crack growth without remeshing[END_REF][START_REF] Laborde | High order extended finite element method for cracked domains[END_REF][START_REF] Béchet | Improved implementation and robustness study of the X-FEM for stress analysis around cracks[END_REF], the displacement discontinuity across the crack is represented using a Heaviside-like function, which is multiplied by the finite element shape functions. An enrichment area of fixed size is defined and some nonsmooth functions are added inside this area. Then, the global unknown function may be written as:

u h = i a i i + j∈J b j H j + k∈K c k f k , (7)
where function H is equal to +1 on one side of the crack and -1 on the other and i stands for the reduced HCT/FVS shape functions. The set J denotes the dof whose associated shape functions support is completely crossed by the crack (see Figure 8). The set K and functions f k are connected with the crack tip enrichment; their exact meaning and expressions will be detailed further, after next section which deals with a numerical problem that we met and the way we chose to overcome it.

An incompatibility between H-enrichment and FVS element

A kind of numerical incompatibility appears between the Heaviside-like enrichment and the FVS element, reduced or not. If a straight crack crosses two adjacent edges of a quadrangle, then a non-invertible linear system may result. The reason is there exists a function, in the FVS functions space, that is completely zero on two adjacent sub-triangles of the quadrangle and not on the two others (see Figure 9). Consequently, if the quadrangle is cut along those two sub-triangles (see Figure 10, left), the regular dof associated with the opposite node are linear combinations of the H -enriched ones. This is not a proof that the resulting matrix will be singular, but is at least an insight. In our numerical applications, with not too coarse meshes, the linear system was always singular. An example of such a mesh is given in Figure 11. In addition, this problem is independant of the distance between nodes and crack: whether they are far or not, the non-invertibility remains.

Let us also remark that it is not particularly due to the fact that the mesh is structured or not, even if, in our tests, the problem arises only with non-structured meshes. By the way, in the case of a straight crack and structured meshes, the quadrangles had always three sub-triangles crossed by the crack (see Figure 10, right), which does not cause any non-invertibility.

To solve this problem, several ideas were studied, like moving the crack or the nodes. However, they are not fully convenient. Moving the crack in order to avoid to cross two sub-triangles is quite complicated and, moreover, it may introduce a significant error. If we choose to move the nodes, it is mandatory to move almost all the nodes belonging to cracked quadrangles and to line up them with the crack. This method can be considered as a remeshing of a part of the domain. It is numerically expensive and not in good agreement with the philosophy of XFEM. So we tried something else. Starting from the fact that this incompatibility problem does not occur with the HCT triangle, the idea is then to divide each quadrangle crossed by the crack into two triangles (see an example in Figure 12). It is an easy and quite cheap method, and we emphasize on the fact that it is not a remeshing: the nodes of the 'new' mesh are exactly the same as those of the 'old' one. To conclude this section, let us notice that it may certainly exist other possibilities to overcome this non-invertibility. This could be an area for further research.

Remark 3

Let us specify that we also tested the possibility to remove the degree of freedom in the square, Figure 10. We found that it leads to a significant loss in accuracy, which prevent to obtain the optimal rate of convergence.

First enrichment strategy: adding dof on each node

This strategy is directly inspired from [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF][START_REF] Béchet | Improved implementation and robustness study of the X-FEM for stress analysis around cracks[END_REF] for plane elasticity problems. Some optimal a priori error estimates are given in (Nicaise et al.; submitted). The main idea is that every node included in a fixed area will be enriched with the nonsmooth functions (see Figure 13). This enrichment is obtained by multiplying the nonsmooth functions with the finite element shape ones of the reduced HCT/FVS element (both nodal and derivatives values), denoted i . Contrary to [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF], the shape functions are the same both in enriched and non-enriched terms. Then, the unknown reads

u h = i a i i + j∈J b j H j + k∈K 4 l=1 c kl F l k . ( 8
)
where the nonsmooth functions are linked to the expression of the bilaplacian asymptotic solution ( 5) and given by F 1 (r, ) = r 3/2 sin 3 2 ; F 2 (r, ) =r 3/2 sin 2 ;

F 3 (r, ) = r 3/2 cos 3 2 ; F 4 (r, ) =r 3/2 cos 2 , (9)
in polar coordinates relatively to the crack tip. Moreover, K denotes the set of dof associated with the nodes that are inside the enrichment zone: they are indicated by a small circle in Figure 13. At each circled node, both reduced HCT and FVS elements have three dof, named for instance

i 1 , i 2 , i 3 .
For each of these dof, we add four nonsmooth enrichment functions: it means that for i 1 we introduce

F 1 i 1 , F 2 i 1 , F 3 i 1 , F 4 i 1 .
The procedure is the same for i 2 and i 3 . So it makes 12 additional dof due to enrichment functions, on each node inside the enrichment area.

Let us also remark that some nodes are enriched both with H and nonsmooth functions.

Two dofs enrichment strategies

From [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF], we know that a large number of additional dof in the enrichment area generally leads to a very high condition number for the associated linear system. A first attempt in decreasing the number of nonsmooth enrichment dofs, and then decreasing the condition number, is obtained by making the following modification in the nonsmooth enrichment.

The asymptotic displacements near the crack ( 5) are the function of the two unknown stress intensity factors K 1 and K 2 . But, in the XFEM enrichment [START_REF] Bordas | Enriched finite elements and level sets for damage tolerance assessment of complex structures[END_REF], these functions are developed on four dof. So an idea is to develop the asymptotic displacements on only two dof. It means that the four F l functions are replaced by the two following ones, say

G l G 1 (r, ) = r 3/2 +7 3( -1) cos 3 2 +cos 2 ; G 2 (r, ) = r 3/2 3 +5 3( -1) sin 3 2 +sin 2 . ( 10
)
So the nonsmooth enrichment term 8) of the approximate solution u h . To conclude, let us just remark that this idea was tested in several references (for instance, in plane elasticity, see [START_REF] Liu | XFEM for direct evaluation of mixed mode SIFs in homogeneous and bi-materials[END_REF]).

Second enrichment strategy: global nonsmooth functions and integral matching

The second strategy is inspired by the so-called 'XFEM dof gathering with pointwise matching', introduced in [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF]. The idea of such a method is to reduce the number of dof that are introduced in the finite element basis. For this, the support of the nonsmooth added functions is the whole enrichment area (see Figure 14) and they are not multiplied by the finite element basis functions. So, instead of 12 additional dof per node inside the enrichment area, there are only four nonsmooth dof for the whole system: they are associated with the functions given by [START_REF] Areias | Analysis of fracture in thin shells by overlapping paired elements[END_REF]. We introduce the boundary of the enrichment area, say . It cuts into two sub-domains: the enrichment area, say 1 , and the rest of , say 2 (see Figure 14). The unknowns defined on each domain i are denoted by u h i and their expressions read

u h 1 = i∈N 1 a i i + j∈J 1 b j H j + 4 i=1 c i F i ; u h 2 = i∈N 2 a i i + j∈J 2 b j H j , (11)
where N 1 and N 2 are the sets of dof that are located in 1 and 2 (N 1 ∩ N 2 is not empty: it corresponds to the set of nodes that are on the boundary ). In a same way, J i are the sets of degrees of freedom of J that are located in i (J 1 ∩ J 2 is not empty for the same reason). Naturally, a matching condition is needed at the interface between the enrichment area and the rest of the domain, in order to insure the continuity of the function and its derivatives. Instead of the pointwise matching introduced in [4], we use the integral matching given by ( 6). Note that a mathematical analysis of the optimality of such a method in the framework of plane elasticity is provided in [START_REF] Chahine | Étude mathématique et numérique de méthodes d'éléments finis étendues pour le calcul en domaines fissurés[END_REF] (Chahine et al.; submitted).

Remark 4

As in the previous section, it is also possible to develop the nonsmooth displacement on the two dofs given by [START_REF] Areias | Analysis of finite strain anisotropic elastoplastic fracture in thin plates and shells[END_REF]. This possibility will be considered further.

NUMERICAL EXPERIMENTS

The numerical experiments presented in this section were performed with the open-source finite element library Getfem++ [START_REF] Pommier | Getfem++, An Open Source Generic C ++ Library for Finite Element Methods[END_REF] developed by our team.

Test problems

We test the previous methods on two benchmark problems of a plate with a straight through crack. The plate is the unit square and the material is assumed to be isotropic. In the experiments, the values of E and were taken at 14.98 kPa and 0.3. Moreover, half the thickness of ε is 0.045. However, let us recall that, for the Kirchhoff-Love model, the thickness is irrelevant with respect to numerical accuracy.

First problem: The exact solution is the pure mode II (see Figure 15) which reads

u ex 1 (r, ) = 1 10 r 3/2 3 +5 3( -1) sin 3 2 +sin 2 , (12)
in polar coordinates. It is quite easy to check whether this function is (2Eε 3 /3(1-2 )) 2 u = 0, and satisfies free-edge boundary conditions on the crack faces.

Second problem: The exact solution is the sum of the pure modes I and II, noted G 1 and G 2 on Equation [START_REF] Areias | Analysis of finite strain anisotropic elastoplastic fracture in thin plates and shells[END_REF], and of a regular part. The expression reads:

u ex 2 = G 1 + G 2 +10[x 4 1 -6 x 2 1 x 2 2 +(2 -1)x 4 2 +2]. (13)
Figure 15. Graph of the first solution, with the Von Mises stress criteria.

The three terms G 1 , G 2 and the regular part remaining are all such that (2Eε 3 /3(1-2 )) 2 u = 0, and satisfy the free-edge boundary conditions on the crack faces. The crack makes a /6 angle with the line {x 2 = 0}.

For both these solutions ( 12)-( 13), there is no volume nor surface forces: the loading is due to non-homogeneous boundary conditions, given by u ex 1 and u ex 2 .

Convergence Curves

The classical and eXtended Finite Element Methods (FEM and XFEM) were applied to the two previous problems. In particular, three XFEM strategies were investigated: they are called XFEM 1 (first strategy with four nonsmooth functions, see Section 4.3), XFEM 1 bis (first strategy with two nonsmooth functions, see Section 4.4) and XFEM 2 (second strategy, see Section 4.5). The L 2 (displacement) and H 2 (energy) errors were measured, on triangular, quadrangular, structured and non-structured meshes. Let us notice that the coarsest meshes had 11 elements on each edge of the domain, and the most refined ones 151. The enrichment area diameter was 0.30.

First test-case: The convergence curves are given in Figures 16 and17. They show that the errors decrease when the mesh is refined. More precisely, we can first observe that the three XFEM are more accurate than classical FEM. Among them, the most accurate is the second enrichment strategy (XFEM 2) and the less one is the first strategy with two nonsmooth functions (XFEM 1 bis). So it exhibits that the use of four nonsmooth functions instead of the two asymptotics nonsmooth displacements leads to more accurate numerical results. A kind of superconvergence behavior also appears for XFEM 1, especially with quadrangular meshes: thanks to it, its accuracy reaches XFEM 2 for the thinnest meshes. Moreover, it can be noticed that, for non-structured meshes and particularly for triangular ones, the displacement norm error is quite 'oscillatory'. However, we can see that the error globally decreases. So, even though this curve is less smooth than the others, the numerical results remain perfectible, but satisfactory.

Then, the rates of convergence were examined. Let us recall that, denoting u and u h the analytical and numerical solutions, the so-called rate of convergence is the number such that u h -u V = O(h ). It also depends on the choice of space V . For example, on regular problems, with reduced HCT/FVS elements, we can expect u h -u H 2 = O(h) and u h -u L 2 = O(h 2 ) (see [START_REF] Ciarlet | Basic error estimates for elliptic problems[END_REF]). So, the results we obtained are summarized in Table I. First, we observe that, for classical FEM, due to nonsmooth displacement at the crack tip, the rates of convergence remain limited around 0.5 in energy norm and 1 in displacement norm. Second, as expected, the XFEM rates of convergence are close to be optimal, i.e. they are the same as those of a regular problem. An exception is XFEM 1 that appears to be superconvergent. Second test-case: Only the two main strategies, named XFEM 1 and XFEM 2, were tested. We did not perform a numerical comparison with a classical finite element method. However, we compared them with XFEM without singular enrichment, which means the use of Equation ( 7) without the last term (it gives similar results compared with classical finite element method). The convergence curves are given in Figures 18 and19, and the rates of convergence in Table II. These results lead to the same conclusions than the first test case:

• XFEM 2 is more accurate than XFEM 1 (however, roughly equal on non-structured meshes).

• The poor rate of convergence due to the crack is improved by the use of singular enrichment.

The optimal rates are reached by both enriched strategies.

Condition number

The application of a finite element method leads to solve a linear system which can be more or less well-conditioned. This section deals with this point in XFEM strategies framework. The numerical results of the first test-case are given in Figures 20 and21, the second test-case in Figure 22.

The classical FEM is taken as a reference. It can be observed that its condition number grows at a rate close to 4 which is the expected rate due to the fact that we solve a fourth-order boundary value problem. As expected, since XFEM uses additional nonsmooth shape functions, its condition number is higher. More precisely, XFEM 1 produces the most ill-conditioned linear systems. It is close to 10 16 for the thinnest meshes. So this phenomenon could limit the capabilities of such a method as, beyond this value, a matrix is close to be non-invertible at the machine precision. Moreover, the rates of growth are the highest: between 6.5 and 8.4. XFEM 2 and XFEM 1 bis produce better conditioned linear systems since their rate of growth remains close to 4, like the classical FEM.

Let us recall that the high condition number of XFEM 1 was already mentioned in [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF][START_REF] Béchet | Improved implementation and robustness study of the X-FEM for stress analysis around cracks[END_REF]. So, to avoid this drawback, alternative enrichment strategies were developed. In [START_REF] Béchet | Improved implementation and robustness study of the X-FEM for stress analysis around cracks[END_REF], the authors introduced a numerical procedure that allows, in a sense, to orthogonalize the enrichment functions. Some key features of our 'Integral Matching' strategy (XFEM 2) are directly inspired from [START_REF] Laborde | High order extended finite element method for cracked domains[END_REF]: the use of global enrichment functions instead of enriching each node decreases the number of dofs and also the condition number. Hence, it is not surprising to observe that the condition number of XFEM 2 is significantly lower than XFEM 1. However, what can be observed here is that the development of the nonsmooth displacement on two dofs, instead of four, is also a very efficient way to reduce the condition number, even if it leads to a slight loss of accuracy.

Remark 5

We also tested the nonsmooth displacement developed on two dofs for the Integral Matching enrichment strategy. The nonsmooth enrichment term 4 k=1 c k F k of ( 11) is simply replaced by k=1 c k G k , where functions G k are expressed in [START_REF] Areias | Analysis of finite strain anisotropic elastoplastic fracture in thin plates and shells[END_REF]. With this modification, the errors remain unchanged: their numerical values are exactly the same on about 3 significant figures. Hence, the curves are the same as those of XFEM 2 in Figures 16 and17. It is why we do not present any graphic here. The condition number is just slightly decreased compared with XFEM 2: it is divided by a factor that varies between 1 and 3 for most of the meshes.

CONCLUDING REMARKS

We presented an adaptation of XFEM to the Kirchhoff-Love plate model. Several XFEM formulations have been investigated and we stress the fact that the features of these methods make possible implementation in an industrial software:

• for all formulations, it was shown that the accuracy was nearly equal to the one of a classical finite element method on a regular non-cracked problem; • the computational cost is nearly equal to the one of the Mindlin-Reissner model, which is the reference model for plate problems in industrial softwares; • our formulation is not subject to numerical locking due to the enrichment, unlike the Mindlin-Reissner model and the strategy presented in [START_REF] Dolbow | Modeling fracture in Mindlin-Reissner plates with the extended finite element method[END_REF].

The numerical experiments suggest that the use of the integral matching formulation is the best choice in terms of accuracy (especially on coarse meshes) and in terms of condition number of the linear system.

Finally, the computation of stress intensity factors, which is one of the key applications of the numerical methods on cracked domains, will be presented in a forthcoming paper.

Figure 1 .

 1 Figure 1. Cracked thin plate (the thickness is oversized for the sake of clarity).

Figure 2 .

 2 Figure 2. Fracture modes for Kirchhoff-Love bending model ( C is the cracked part of the boundary).Left: a symmetric bending leads to mode I. Right: a shear bending leads to mode II.

Figure 3 .

 3 Figure 3. HCT triangle and FVS quadrangle. Location of degrees of freedom and sub-triangles.

Figure 4 .

 4 Figure 4. Shape functions of the reduced HCT. Left: nodal value. Right: x 1 -derivative value.

  . The shape functions corresponding to the Lagrange dof, and to the x 1 -derivative dof are shown (the remaining

Figure 5 .

 5 Figure 5. Shape functions of the reduced FVS. Left: nodal value. Right: x 1 -derivative value.

1 Figure 6 .Figure 7 .

 167 Figure 6. Matching decomposition of the domain .

Figure 8 .Figure 9 .

 89 Figure 8. Set of nodes to be enriched along the crack.

Figure 10 .

 10 Figure 10. Left: only two sub-triangles are crossed, which leads to a singular matrix. Right: three sub-triangles are crossed, which leads to an invertible matrix.

Figure 11 .

 11 Figure 11. Left: example of mesh leading to a non-invertible linear system, with the crack (in bold line). Right: zoom on the crack.

Figure 12 .

 12 Figure 12. Left: entire modified mesh (the straight line is the crack). Right: zoom on the crack.

Figure 13 .

 13 Figure 13. Set of nodes enriched by the nonsmooth functions.

k∈K 4 l=1 c kl F l k becomes k∈K 2 l=1

 42 c kl G l k in the expression (

1 Figure 14 .

 114 Figure 14. Set of elements which represents the support of the nonsmooth functions (set 1 ).

Figure 20 .

 20 Figure 20. First test-case, condition numbers of FEM and XFEM first strategy (with 2 or 4 nonsmooth dofs).

Figure 21 .

 21 Figure 21. First test-case, condition numbers of FEM and all XFEM.

Figure 22 .

 22 Figure 22. Second test-case, condition numbers of all XFEM (with/without singular enrichment).

2

 2

Table I .

 I First test-case, convergence curves for FEM and XFEM on structured meshes. Top: triangular meshes. Bottom: quadrangular meshes. Rates of convergence for FEM and XFEM.

			Errors of FEM and XFEM (Triang. Struct. meshes)
			Structured meshes	Non-structured meshes
			H 2 norm	L 2 norm	H 2 norm	L 2 norm
	FEM Tri.		0.51	0.95	0.50	0.94
	FEM Quad.		0.51	0.94	Not tested
	XFEM 1 bis Tri.		0.95	1.77	1.10	2.11
	XFEM 1 bis Quad. XFEM 1 Tri. XFEM 1 Quad. XFEM 2 Tri. XFEM 2 Quad.	Relative Error Relative Error	0.94 1.19 1.44 1.02 1.05	1.73 2.36 3.23 2.09 2.11	1.35 1.27 1.72 1.04 1.06	2.61 2.67 3.65 1.94 1.95
	Mesh parameter Errors of FEM and XFEM (Quadr. Struct. meshes) Mesh parameter Mesh parameter Figure 16. Mesh parameter Relative Error Relative Error Figure 17. First test-case, convergence curves for FEM and XFEM on non-structured meshes.
			Top: triangular meshes. Bottom: quadrangular meshes.

Table II .

 II Figure 18. Second test-case, convergence curves for FEM and XFEM on structured meshes.Top: triangular meshes. Bottom: quadrangular meshes. Rates of convergence for FEM and XFEM.

			Test 2: errors of XFEM (Tri. Struct. meshes)
			Structured meshes	Non-structured meshes
			H 2 norm	L 2 norm	H 2 norm	L 2 norm
	No. singul. Tri.		0.60	1.01	Not tested
	No. singul. Quad.	0.58	1.01	Not tested
	XFEM 1 Tri.		1.05	2.35	1.06	2.07
	XFEM 1 Quad. XFEM 2 Tri. Relative Error XFEM 2 Quad.		1.11 0.99 0.99	2.42 2.13 1.65	1.09 1.02 1.08	2.44 1.94 2.30
			Condition Number of FEM, XFEM 1 and XFEM 1 bis
		10 16
		10 14
		10 12		Mesh parameter
		10 10	Test 2: errors of XFEM (Quad. Struct. meshes)
		10 8
		10 6
	Relative Error	Relative Error		Mesh parameter
				Mesh parameter Mesh parameter
	Figure 19. Second test-case, convergence curves for XFEM on non-structured meshes.
			Top: triangular meshes. Bottom: quadrangular meshes.

((1-2 )/Eε(3+ )). The scalar coefficients K 1 and K 2 are the so-called 'Stress Intensity Factors'. They are widely used in fracture mechanics for crack propagation.

ACKNOWLEDGEMENTS

This work was supported by a joint contract with Airbus (no 001847) and 'Centre National de la Recherche Scientifique' (no 060131), and also supported by 'Agence Nationale de la Recherche' (project ANR-05-JCJC-0182-01). The authors wish to thank Marc Balzano (Airbus France) and Patrick Laborde (University of Toulouse, France) for their confidence and support all along this work.

Contract/grant sponsor: Airbus; contract/grant number: 001847 Contract/grant sponsor: Centre National de la Recherche Scientifique; contract/grant number: 060131 Contract/grant sponsor: Agence Nationale de la Recherche; contract/grant number: ANR-05-JCJC-0182-01


Yves Gourinat
email: yves.gourinat@isae.fr

Uwe Apel

Franck Delbart

The Aerocrew Mission : Training Space Session at Ny Aalesund Arctic Base

The Aerocrew mission was realized at Ny-Aalesund Arctic Base in December 2007, in the frame of the International Polar Year, and in cooperation with the French Polar Institute. The team made an original training experience, constituting an integrated space crew, including physicians, aerospace trainers and engineers, implied in a seminar with 4 sessions, dealing with the training capabilities of Arctic Bases. This kind of base constitutes a pertinent and affordable facility for aerospace teams, and the specific aerospace crew training techniques, are fruitful for the polar scientists (glaciologists, geologists, specialists of the atmosphere).

The sessions, given by professionals of aerospace, robotics and medicine, covered the training methods for crews, robotics for outdoor and indoor activities, engineering of embedded systems, and the arrangement of crafts. The experience has shown the efficiency of a transverse visiting multidisciplinary team for training, and possible synergies with the resident scientists. In addition, the sessions were enriched by demonstrations such as mini-robot for observation, micro-helicopter, and also the comparison between EVA Russian glove and Polar Suits.

2) Scope and General Context

The previous works developed with Swedish and Norwegian specialists aimed at using the polar bases as testbed for space crew management -risk, isolation, scientific team in exotic situation -and spacecraft internal arrangement. Concerning the risk, a specific list had been published (ref. [START_REF] Gourinat | Sweden Specific Polar Experiments for Long Range Flight Ergonomics[END_REF]), including the classical concepts mentioned in fundamental studies realized on arctic and space crews (ref. [START_REF] Stuster | Lessons from Polar and Space Exploration[END_REF]), but pointing in addition some other ones such as : hearth paradigm, training in relation with telework, and even mourning in space. Concerning the internal arrangement, the basic study was published by Louise Jagger in 1998 (ref. [START_REF] Gourinat | Spacecraft Architectural Design with Minimal Artificial Weightiness Concept[END_REF]), concluding a study financed by Matra Marconi Space and dealing with the Minimal Artificial Weightiness Concept. During the following studies (Reine Kurth, Francesca Molendini), it appeared clearly that the internal arrangement of arctic (an Antarctic) bases were in strong synergy with spacecraft for long range flight, including in particular the perception of light in polar night (and day). Consequently, the scientific cooperation with northern bases and centers went on during the following years, and it appeared that the International Polar Year could be an excellent opportunity to make a concrete mission with realistic resources, and in a specific scientific frame.

The goal defined for this mission was not only the investigation of classical direct spinoff from arctic bases in relation with space training, but also the inverse one considering that space techniques -with focus on training and robotics -are unique for the users of polar bases. In addition, the crew aspect has been emphasized, considering that any technique is finally useful for human purpose (again, this is valid in space and in base…).

For this purpose, the crew for the AEROCREW Mission was chosen in function of coherence and complementarity with an actual planetary crew, taking into account the training priority.

It gathered eleven members (four women & seven men) with three physicians (one surgeon, one aerospace & sport specialist, one ear specialist), one specialist in life sciences, one specialist in physical sciences, one satellite engineer, one researcher in robotics and command systems, one aircraft crew trainer, one researcher in aerospace embedded systems, and two researchers and professors in mechanics and space techniques.

3) Training Sessions

Four technical sessions were realized during the mission, each session dealing with a specific thematics, with a formal presentation by the specialist, and a practical application of proposed protocols. The first one was precisely devoted to Astronaut Training Methodology, and possible applications for the preparation and follow-up of arctic crews .This session was proposed by Delphine Gourdou, former trainer for astronauts and currently trainer for aircraft crews (which also represents an interesting spinoff).

The different procedures were presented, with emphasis on the necessity of rigorous and clear specification for any hugh level technical training for crews, and also the importance of initial consideration of crew communication protocols for defined tasks. The different approaches of written documents (in particular in US an Russian context) was highlighted.

A practical exercise was then made about Extra Vehicular Activities, with time limitated action and standard orders given by the base to the crew member in EVA. This is directly applyable for external worker in polar medium, with the safety and performance stakes. This exercise was made in the context of robotic assistance in EVA, which was developed in other sessions.

It is also necessary to precise that the crew human context was comparable, for training, to the context for space crews : coexistence specialized high level scientists, with major mission scopes, in an exotic context, and transverse crews whose mission is rather technical and organisational. Actually, this matricial aspect of crew activity must be reflected by the organisation of training strategy and communication method.

In fact, the three other sessions were devoted to scientific presentations for all the other crew members : command and autonomy in robotics, system engineering and specification, space mechanisms and long range flights. These sessions aimed at sharing the technical stakes in a international and multicultural team, which is the common case for a polar base.

4) Man-Machine Interface and Robotics

Beside the sessions devoted to training and crew management aspects, a specific session conducted by Marc Bertin -researcher in command theory at Paris Sorbonne -dealed with robotics and commandability of systems, and had a direct impact on human-machine complementarity for planetary exploration, in the arctic context. The demos showed two micro-vehicles, and demonstrated in particular ability for transmission from the boarded camera to vritual-reality glasses. The stake is here, on one hand the possibilities open for pioneering robots (in new medium), and on the other hand, the autonomy in repetitive tasks for servitude robots in the craft itself. This kind of small and adaptative robot will certainly constitute the standard facility both for ice researcher, in order to search, find and examine the convenient specimens out of the base, and for planetar future explorations. At the end of the session, a micro-helicopter (6" long) was tested in external flight, with difficult control due to temperature, wind and polar night ; the pertinence of this kind of micro flying vehicle is clear for observation and access to sharp places.

These techniques are also commonly envisaged for terrestrial applications such as the observation of animals, vegetation, and any mission for survey where the environment could be hostile.

As the human capability cannot be replaced by robot in all situations, the Orlan Russian Extra Vehicular spacesuit glove (given by Space City, Toulouse) was tested in the snow, as it could be used by ice specialists or by planetar geologists. As already observed in space, it appeared that the thermal protection on silicon finger tips is relatively weak.

The test of space suits and garments in arctic zone, with enventual ice and ground drilling in situation is an excellent testbed not only for space planetar applications -future lunar or martian missions -but also constitutes a good opportunity to prepare the use in severe Antarctic conditions.

5) Spacecraft and Base Arrangement

The goal is here to use the arctic base as an analogue for a spacecraft or a future base on the Moon or Mars; extreme terrestrial environments are indeed privileged places to conceive and investigate habitats for space exploration. The first activity (A. Pacros) consisted in a general survey about the overall architectural design of the polar base; it will be done using a questionnaire adapted from the one developed by the MDRS Crew 43.

The second theme (Y. Gourinat) was a study of the importance of Hearth and Light in hostile environments, based on previous studies devoted to spacecraft (ref. [START_REF] Gourinat | Spacecraft Architectural Design with Minimal Artificial Weightiness Concept[END_REF]). Lectures and panel sessions have been performed in a collective approach about the hearth paradigm (fireplace, light, wardroom, crew quarters), to initiate an exchange between the scientists staying in the station and the visiting Aerocrew.

The third aspect was called "Art-Habitat" (R. Binot). Exploration is often facing hostile conditions for which psychological balance is necessary. How can the technicalities of the habitat and the local environment conditions be used towards an esthetical environment and psychological balance ? The idea is here to use Art to overcome the contrast between the nature (cold temperatures, boreal auroras, very primitive world dominated but rich/diverse, rare/very localised biology), and the functional habitat and motivated scientific/technical crews. It should be an evolutionary art that depends on the competences and interests of the crew (residents and [successive] visitors), the availability of materials, and the actual base architecture. Examples could be to use video sequences, painting with light to reflect seasons, or creating a Japanese garden.

6) Conclusion and Perspectives

The mission demonstrated clearly the pertinence of aerospace training methodology for polar crews, the challenge and context being similar. It was also the occasion to demonstrate some particular situations such as a visio-conference with training session in France, which also induces analogies with space missions : the astronauts keeping in touch with their terrestrial job. In any way, the crossed cooperation between aerospace trainers and arctic teams is obviously promising and opens wide and affordable possibilities.

The Arctic bases constitute are unique to organize high level scientific seminars, not only for the exceptional environment, but also for the human context : specialized crew with competences in complementarity. The space training techniques are also very pertinent for the training of polar crews. In addition, the specific robotic applications are the occasion to test space systems, but also to help the polar researchers in particular activities such as ice specimen return, atmosphere analysis on site and in flight or bio-specimens in extreme media.

6 )

 6 Conclusion and Perspectives 1 Professor, Head of SUPAERO Mechanical Training Unit 2 Professor, Vice-Rektor 3 CNRS Research Engineer, IPEV French Polar Institute Arctic Logistics Dept Manager 1) Introduction In the frame ot the International Polar Year, the Alfred Wegener Institut and Institut Polaire Français Paul-Émile Victor have welcome the AEROCREW Mission at Ny-Ålesund, during one week in December 2007. This short mission was nevertherless original in its scope, Organized by Terre&Espace, a French Association devoted to applications of space techniques (www.terrespace.org), after discussion with the Polar Institute. The AEROCREW mission aimed at training and identify possible spinoff between aerospace crews preparation and techniques on one hand, and scientific research in polar bases on the other hand. One fundamental scope was actually the demonstration of training synergies between arctic and aerospace crews, as the scientific and human challenges are similar : high level scientific team in hostile and isolated medium, with operational scopes. The AEROCREW After the presentation of the general context of the mission, the paper describes the training sessions, and then the activities devoted to man-machine complementarity in base. The conclusion and perspectives are finally mentioned.

Acknowledgements

Gérard Jugie, Director of Institut Polaire Français Paul-Émile Victor, CNRS Director of Research Rainer Vockenroth, Alfred Wegener Institut, Head of AWIPEV Base Anne Pacros, Roger Binot, European Space Agency Delphine Gourdou, Airbus Training