G Legrain

C Geuzaine

J F Remacle

N Moës

P Cresta

J Gaudin

Numerical simulation of CAD thin structures using the eXtended Finite Element Method and Level Sets

An efficient approach is proposed in order to predict the mechanical response of complex industrial parts. As these structures are usually composed by massive and thin parts, different models have to be mixed together (plate, shells, solid). The transition between these different kinematic assumptions can be problematic and non-linear models cannot be employed depending on the plate model that is considered. Moreover, Finite Element analysis usually implies tedious meshing steps in the case of large and complex assemblies, with tedious and manual meshing steps. The idealisation and simplification of these structures into a mix of 2D and 3D Finite Elements usually takes therefore significantly more time than the analysis itself. The objective of the present contribution is to explore a calculation process that enables a simple automation of the meshing steps. Even though potentially computationally more expensive, the meshing automation may lead to drastic time reduction for the CAD to mesh process and a much tighter link between CAD and calculated assembly. Finally, easier and faster design explorations would be allowed. This strategy relies on the use a nonconforming quadratic approximation that is defined on a sufficiently fine mesh. The eXtended Finite Element Method is used in order to alleviate meshing issues. The mesh and Level-Set function are built from the CAD input, by means of an automated approach. The strategy is verified against analytical solutions and real aerospace substructures.

Introduction

The assessment of the mechanical behaviour of complex industrial parts is still an issue nowadays. Usually, these parts are composed of both thin and massive areas that are assembled together along junctions. The mechanical behaviour of these massive and thin areas is different, so that different models are usually considered. The usual path for such analysis begin with the CAD representation of the structure: both thin and thick areas are represented as three-dimensional volumes which are described by their outer boundary (BREP). This volumetric representation is converted into plates or shells by means of the extraction of the mid-plane of the part. Extracting this mid-plane cannot always be automatized for complex shapes, leading to some tedious and time consuming issues (difficulties to obtain continuous mid-planes with part assembly). However, the constitutive equations must also be rewritten for plate and shells models, which can be problematic for non-linear problems. Moreover, a special care as to be taken on the junction between these 2D zones and the 3D massive parts. Usually, the kinematic link between the two areas induces stress concentrations. This is why in some highly secured applications, solid modelizations have to be considered for certification purpose [START_REF] Lee | Automatic adaptive fe analysis of thin-walled structures using 3d solid elements[END_REF]. In addition, such a kinematic hypothesis decreases the accuracy of the model near junctions. In order to overcome these issues, either solid elements can be used, or special "solid shell" elements can be considered. Using low order solid elements requires a large number of degrees or freedom in order to be able to represent the mechanical fields in the parts. This is also due to the bad behaviour of such elements when their aspect ratio increases.

The use of high order finite elements, as advocated in [START_REF] Duster | pq-adaptive solid finite elements for three-dimensional plates and shells[END_REF] allows to get rid of these issues. Using solid shell elements allows one to represent plate-like behaviour by means of displacement-only finite elements. Usually, these elements involve incompatible modes [START_REF] Wilson | Incompatible displacement models, volume numerical and computer models in structural mechanics, chapter numerical and computer models in structural mechanics[END_REF][START_REF] Rl | A nonconforming element for stress analysis[END_REF] and selective integration [START_REF] Zienkiewicz | Reduced integration technique in general analysis of plates and shells[END_REF]6,7] in order to be computationally efficient. In addition, the constitutive equation has still to be rewritten for some of these elements. Note also that equivalent mid-planes have to be defined as the integration is done along this preferential direction. Alternatively, the transition between thin and solid parts can be represented by means of the Arlequin method [8].This method allows to couple models involving different physics, and relies on the partition of the energy in the transition zone between solid and shell.

Albeit these modelization issue, the use of state-of-the-art industrial simulation tools such as FEM usually imply tedious meshing steps in the case of large and complex assemblies, with complex and manual meshing steps. The idealisation and simplification of these structures into a mix of 2D and 3D Finite Elements usually takes significantly more time than the analysis itself. This is a major drawback if many complex designs are to be explored: the meshing parametrization may have to be done separately from the CAD parametrization or with complex coding of case-dependent meshing rules. This is why an alternative and pragmatic path is proposed here for the simulation of thin structures. This strategy relies on two main ingredients: (i) the eXtended Finite Element Method (X-FEM) [9] and (ii) the Level-Set method [START_REF] Osher | Fronts propagations with curvature dependent speed: Algorithms based on Hamilton-Jacobi formulations[END_REF]. The eXtended Finite Element Method specificities are used to explore a calculation process that enables a simple automation of the meshing steps.

Even though potentially computationally more expensive, the meshing automation allowed by the X-FEM may lead to drastic time reduction of the CAD to mesh process and a much tighter link between CAD and calculated assembly. This process may therefore allow easier and faster design explorations in an industrial context.

The X-FEM belongs to the class of the partition of unity finite element methods [START_REF] Melenk | The partition of unity finite element method : Basic theory and applications[END_REF]. It was first proposed in the context of fracture mechanics as a mean to avoid remeshing issues. Since this first work, the X-FEM has been applied to a large number of applications (see [START_REF] Fries | The extended/generalized finite element method: An overview of the method and its applications[END_REF][START_REF] Yazid | A state-of-the-art review of the X-FEM for computational fracture mechanics[END_REF] for reviews on this subject), all trying to avoid or limit remeshing. In particular, the X-FEM allows to solve mechanical problems on meshes that are independent of the geometry (like advocated in fictitious domain approaches [START_REF] Saulèv | On the solution of some boundary value problems on high performance computers by fictitious domain method[END_REF]). Usually, part or discontinuity surfaces are represented by means of the Level-Set method [START_REF] Osher | Fronts propagations with curvature dependent speed: Algorithms based on Hamilton-Jacobi formulations[END_REF]. The Level-Set method allows one to define implicitly surfaces as the iso-zero of a scalar function that is usually taken as the signed distance function to the surface. The Level-Set function is evaluated at some grid points, then interpolated on the domain. The iso-zero of the Level-Set is usually located inside elements, which leads to the necessity to use non-conforming methods such as the X-FEM for the simulation over Level-Set domains.

The strategy proposed here consists in translating the explicit CAD representation of the domain into a Level-Set and a mesh that can be used in order to simulate the behaviour of mechanical parts, using the X-FEM. This method has to be tolerant with damaged CAD files and accurate. The accuracy is obtained thanks to the use of quadratic elements that have been reported to represent a good compromise between accuracy and computational cost for modelling thin structures with solid elements [START_REF] Lee | Automatic adaptive 3-d finite element refinement using different-order tetrahedral elements[END_REF][START_REF] Lee | Automatic adaptive refinement finite element procedure for 3d stress analysis[END_REF].

The paper is organised as follows: First, the strategy for efficient mesh and Level-Set creation from a CAD file is presented and illustrated on some real geometries. Then, the X-FEM is recalled.

In a third section, the use of solid elements is assessed in the case of thin and moderately thin structures, both for FEM and X-FEM. Finally, numerical examples are proposed in order to illustrate the use of the proposed approach for real structures, before concluding.

Mesh and Level-Set creation

Presentation of the algorithm

In the eXtended Finite Element Method context, the geometrical features are taken into account using the Level-Set approach. In this case, the geometry is represented implicitly as the iso-zero of a scalar function called Level-Set function φ(x). This function is usually chosen as the signed distance to the curve (resp. surface) of the object so that its gradient has unit norm. Level-Set primitives can be defined for simple geometries such as planes, cylinders, spheres together with boolean operations. This allows to define geometrical objects of quite complex shape, but not in the case of spline and nurbs-based geometries which are common in CAD parts. An approach has been proposed recently by Moumnassi et al. [START_REF] Moumnassi | Finite element analysis on implicitly defined domains: An accurate representation based on arbitrary parametric surfaces[END_REF] in order to construct the Level-Set corresponding to such surfaces with a given geometrical accuracy. This approach gives the opportunity to translate the CAD representation of parts into a Level-Set field and its mesh support. However it is not clear how does the approach behaves when faced with damaged CAD with holes, superimposition between surfaces, or in the case of parts involving numerous surfaces (in term of processing computational time). In this paper, a pragmatic alternative is considered: it is proposed to rely directly on the knowledge of the explicit surface of the part, then to compute the signed distance to this surface. Albeit being less accurate than the aforementioned strategy, the proposed approach is computationally efficient and tolerant to damaged CAD geometries. This strategy is based on the opensource software Gmsh [START_REF] Geuzaine | Gmsh: A 3-d finite element mesh generator with built-in preand post-processing facilities[END_REF] that can load CAD parts thanks to the integration of the OpenCASCADE [START_REF]Open cascade technology, an environment for 3d modeling[END_REF] opensource CAD engine. The input of the process is the definition of the geometry as a CAD file format such as IGES or STEP. Then, the procedure is defined in three steps illustrated in figure 1: 1. Immerse the CAD part in a regular mesh based on the bounding-box of the geometry; 2. Keep only the elements that contain some matter; 3. Compute the distance from the nodes of this minimal mesh to the surface of the part. As the approach is dedicated to plate-like structures, the algorithm takes advantage of the characteristic features of such structures, so that the input of the procedure relies on a small set of parameters:

1. The size of the elements; 2. The maximum thickness of the part;

3. The number of sub-grids that have to be generated. This last parameter is linked to the fact that contrary to Moumnassi et al. [START_REF] Moumnassi | Finite element analysis on implicitly defined domains: An accurate representation based on arbitrary parametric surfaces[END_REF], the different surfaces of the model are not decomposed, which means that the edges and angles of the geometry cannot be represented on coarse meshes by the Level-Set (which tends to smooth sharp geometrical features). This is illustrated in a simple 2D case in figure 2(a-b). In order to overcome this issue, the finite element mesh is refined near sharp features such as edges or angles, see figure 2(c). A multilevel approach is considered, and the elements containing such features are recursively refined. This approach differs from an octree [START_REF] Frey | Mesh Generation: Application to Finite Elements[END_REF], as no database is created to navigate through the element hierarchy (it avoids the overhead due to the classical octree database). The construction of the Level-Set mesh is now detailed, the objective is to make it as efficient as possible. The general description of the algorithm is as follows (see figure 3 for an illustration of the different steps): a) Create the STL discretization of the part (to obtain, from a prescribed geometrical accuracy a discretization of the surface of the part); b)

• Create a cloud of points on the STL discretization. This cloud is such that the minimal distance between two points is less than the elements size.

(b) (a) (c)
• Create a cloud of points on the edges of the model. This cloud is such that the minimal distance between two points is less than the elements size. e) Activate level 0 cells sufficiently close to the surface (distance smaller than the thickness of the part).

f ) Level 0 mesh is built.

g) Activate level > 0 cells containing at least one point of the edge cloud (i.e. crossed by one edge), or sufficiently close (distance smaller than the thickness of the part divided by 2 level-1 ). Note that invalid cells may be introduced at this point (i.e. cells that are not completely filled by their sons). These cells are depicted in red in Figure 3 g).

h) Cancel invalid cells (see Algorithm 2).

i) Compute the distance from all active cells to the the STL faces of the model. j) Cancel cells out of the domain.

k) Find hanging nodes.

This general description is detailed in algorithm 1, and step g) is covered with more details in algorithm 2. This methods allows to build very efficiently the mesh and the Level-Set that define implicitly the part.

Illustration

The strategy presented above is now illustrated on some CAD parts. The idealized CAD parts In particular, we will study the influence of the number of levels and the size of level 0 elements. Note that in all these examples, we will consider level 0 elements whose size is close to the thickness of the part (typically close to one half on the thickness): The use of such mesh density will be justified in section 4.

Perforated plate

The geometry of the plate is 30 × 90 × 2mm. Two types of level 0 meshes are considered: (i)

L x = L y = L z = 1.2mm and L x = L y = L z = 0.
6mm. These meshes are relatively coarse with respect to the thickness of the plate. The influence of the number of sub-grids levels is studied in figures 5 et 6 for mesh L x = L y = L z = 1.2mm and 7 and 8 for mesh

L x = L y = L z = 0.6mm.
In these figures, one can see the number of nodes, elements and the computational time required to create the mesh and the level set (on a Quad-Core AMD Opteron Processor 2376, with only one core used). The influence of the refinement near the edges of the model is clearly visible in figures 6 and 8. One can also see that the processing time for a very accurate geometrical accuracy (292100 elements) took only 10 seconds. The error in the volume of the part is given in table 1: it can be seen that the initial volume error is decreased by about 25% when increasing the number of sub-grids. it enables to obtain a good geometrical approximation with a local addition of elements. Looking at the error on the volume of the domain (Table 2), one can see that it remains very small (less than 1%) in all cases although the geometrical description of the sharp features can be quite poor.

2 ep a) b) c) d) e) f) g) h) i) j)
Concerning the computational time, level sets on extremely fine meshes (typically 1.5 millions of elements) can be generated in less that 2 minutes.

Stiffened panel

Finally, consider the stiffened panel depicted in figure 4(c). The geometry of the part is also close to structures that can be found in the aerospace industry. The size of the bounding-box of the part is 38.7 × 600 × 1000mm, and the maximum and minimum thickness are respectively 6 and 2.94mm. Level 0 mesh is composed of

L x = L y = L z = 1.
5mm elements (about one half of the minimum thickness). By imposing this size for the elements, one can see in figure 11 that only one level of grid ensures a good geometrical approximation (see figure 12). In this example, it is also possible to keep the same geometrical accuracy by taking into account the particular geometry of the part: it is extruded along z direction. Thus, it is appealing to choose bigger elements in this direction: a L x = L y = 1.5mm, L z = 4mm level 0 grid is considered. The resulting mesh is depicted in figure 13, and the corresponding Level-Set in figure 14. The size of the mesh is drastically reduced without sacrificing the geometrical description. Note however that finding preferential directions can be difficult in practice, as these directions are usually associated to local axes of the part: these axes are usually not aligned with the global axis of the computational domain. The processing time is longer in this example (6140s), because of the size of the point cloud (61, 244, 302 points): this is linked to the thickness that is extremely small with respect to the size of the structure.

(a) (b) (c)
To conclude this section, we have proposed a systematic and pragmatic approach for the construction of the implicit representation of complex CAD structures. The computational efficiency of the approach was demonstrated. One can note that such a strategy can easily be parallelized in order to further improve its performances (the faces of the model and the nodes of the cloud can be scattered on different processes). Moreover, note that the strategy can also be applied to massive parts, although a loss of efficiency may occur due to the introduction of a large number of cells out of the part (as the "thickness" parameter is large in this case). Algorithm 1: Computation of the Level-Set.

Step i) are referring to Figure 3.

Input: CAD definition, Element size (Level 0) L x , L y , L z , maximal thickness R max , number of levels Levels Deducted parameters: resolution = min(L x , L y , L z ) and R tube = 2 max(L x , L y , L z ).

Step a) Create the STL discretization of the part.

Step b) foreach Face f do

Based on the STL discretization of the geometry, create a cloud of points on the part's surface, such as the maximum distance between every points is smaller than resolution. end if Levels > 1 then

foreach Edge e do Create a cloud of points on the part's edges, such as the maximum distance between every points is smaller than resolution. Step i) Computation of the level set : foreach CAD Surface f do foreach STL face f stl do Compute the distance from all active nodes to face f stl end end

Step j) Cancel cells out of the domain (positive Level-Set for all nodes).

Step k) Find orphan (hanging) nodes Output: Mesh containing the CAD + level set + hanging nodes data Algorithm 2: Elimination of the invalid cells.

Input

L x = L y = L z = 1.2mm 0 10.38 L x = L y = L z = 1.2mm 1 7.48 L x = L y = L z = 1.2mm 2 6.92 L x = L y = L z = 0.6mm 0 2.69 L x = L y = L z = 0.6mm 1 2.1 L x = L y = L z = 0.6mm 2 2.04
L x = L y = L z = 6mm 0 0.73 L x = L y = L z = 6mm 1 0.689 L x = L y = L z = 6mm 2 0.526

The eXtended Finite Element Method

The X-FEM is an extension of the finite element method (FEM) that was developed from the need to improve the FEM approach for problems with complex and evolving geometries. In contrast to classical finite elements, the X-FEM does not require the mesh to conform the geometry [START_REF] Sukumar | Modeling Holes and Inclusions by Level Sets in the Extended Finite Element Method[END_REF][START_REF] Belytschko | Structured extended finite element methods of solids defined by implicit surfaces[END_REF].

Instead, one can consider a regular mesh that defines the domain of interest without describing geometrical features such as cracks, inclusions or free surfaces. Then, the approximation is enriched by means of suitable functions to take into account these features. More precisely, the X-FEM approximation of the displacement field, u, over an element is given by:

u(x)| Ωe = n α=1 N α   u α + ne β=1 a α β ϕ β (x)   (1)
where the approximation can be divided into a classical one that depends only on the shape functions N α (x) and classical degrees of freedom (dofs in the following) u α , and an enriched one that depends on enrichment functions ϕ β (x) and enriched dofs a α β . Those functions allow to take the non-conforming geometrical features into account in the approximation. The additional degrees of freedom are only added to the nodes whose support is split by the interface, which means that typically only a few of them are added. In the case of free surfaces, as in this contribution, no enrichment has to be considered: only the matter part of the domain is considered when integrating the weak form, and no enriched dofs are introduced [START_REF] Sukumar | Modeling Holes and Inclusions by Level Sets in the Extended Finite Element Method[END_REF]. Because of the non-conforming geometry, the integration process is conducted using a modified Gauss quadrature scheme, as described in [9]. The number of integration points over each subdomain is chosen so that the integration is 'exact'. In the present contribution, quadratic approximations will be considered in order to model the mechanical fields that occurs in plate structures with a minimum number of elements. As presented in the previous section, the geometry on the structure is defined implicitly on nested meshes. This type of mesh produces so-called hanging nodes that have to be treated properly in order to ensure the conformity of the approximation. Some strategies have been proposed in order to overcome this issue [START_REF] Tabarraei | Extended finite element method on polygonal and quadtree meshes[END_REF][START_REF] Legrain | On the use of the eXtended Finite Element Method with Quatree/Octree meshes[END_REF][START_REF] Fries | Hanging nodes and XFEM[END_REF]: In the present contribution, the approach proposed by Legrain et al. [START_REF] Legrain | On the use of the eXtended Finite Element Method with Quatree/Octree meshes[END_REF] will be considered. It is now necessary to characterize the approximation that is used in the approach. Numerous contributions studied the simulation of thin structures with finite elements. If low order elements are used, their bad approximation of the physical behaviour of such structures leads to locking phenomena. This phenomenon can appear for both 3D and 2D modelling of plates and shells. In both cases dedicated formulations has to be taken into account. If the plate is idealized by its mid-plane, a 2D model can be considered. However, mixed formulations [START_REF] Zienkiewicz | A robust triangular plate bending element of the Reissner-Mindlin type[END_REF][START_REF] Bathe | Four-node plate bending element based on Mindlin/Reissner plate theory and a mixed interpolation[END_REF][START_REF] Bathe | Formulation of general shell elements -The use of mixed interpolation of tensorial components[END_REF][START_REF] Bathe | A simple and effective element for analysis of general shell structures[END_REF][START_REF] Batoz | Evaluation of a new thin plate quadrilateral element[END_REF] or reduced integration approaches [START_REF] Zienkiewicz | Reduced integration technique in general analysis of plates and shells[END_REF]6,7] must be used in order to avoid locking issues. Mixed formulations are based on the introduction of rotational degrees of freedom in the approximation. In contrast, under-integrated elements allows to use a displacement based formulation. Unfortunately, selective reduced integration causes zero energy modes that have to be eliminated [START_REF] Flanagan | A uniform strain hexahedron and quadrilateral with orthogonal hourglass control[END_REF]. Moreover, the material law have to be modified for these models, which can be difficult for non-linear constitutive laws. In addition, the transition between such bidimensional areas to massive one can lead to stress concentrations if a rigid link is prescribed between the two models. Richer approaches have thus to be considered, but they lead to a computational overhead, especially for complex structures. An alternative consists in the use of three-dimensional elements.

It is known that displacement-based low order elements are not able to represent the behaviour of slender structures, as locking occurs for both thick and thin plates. This is why so called solid-shell elements have been proposed in the literature. Two main approaches arise: the first one is based on the use of low order finite elements whose behaviour has to be enhanced in order to represent the kinematics of shells [START_REF] Simo | A class of mixed assumed strain methods and the method of incompatible modes[END_REF][START_REF] Simo | Geometrically non-linear enhanced strain mixed methods and the method of incompatible modes[END_REF][START_REF] Simo | Improved versions of assumed enhanced strain trilinear elements for 3d finite deformation problems[END_REF]. The numerical efficiency of those elements is usually achieved using selective integration approaches [START_REF] Zienkiewicz | Reduced integration technique in general analysis of plates and shells[END_REF][START_REF] Hughes | Generalization of selective integration procedures to anisotropic and nonlinear media[END_REF][START_REF] Prathap | Reduced integration and the shear flexible beam element[END_REF][START_REF] Flanagan | A uniform strain hexahedron and quadrilateral with orthogonal hourglass control[END_REF], which can also lead to hourglass modes. The last possibility was advocated recently by Duster et al. [START_REF] Duster | pq-adaptive solid finite elements for three-dimensional plates and shells[END_REF]: The authors proposed to use very large conforming elements, but with a high polynomial order. This approach was shown not to be sensitive to locking issues if sufficiently high order polynomial basis functions are used (typically order 6 to 8). The authors applied this strategy to both Reisner-Mindlin plates (using 2D p-FEM) [START_REF] Rank | High order finite elements for shells[END_REF] and three-dimensional curved structures (using 3D p-FEM) [START_REF] Duster | pq-adaptive solid finite elements for three-dimensional plates and shells[END_REF].

An intermediate possibility for the use of solid elements was proposed by Lee and Xu [START_REF] Lee | A new automatic adaptive 3d solid mesh generation scheme for thinwalled structures[END_REF][START_REF] Lee | Automatic adaptive fe analysis of thin-walled structures using 3d solid elements[END_REF] for applications in offshore structures where both thin and thick plates coexist. Performing a cost study for quadratic and cubic approximations, the authors deduced that quadratic elements were a good compromise between computational cost and accuracy [START_REF] Lee | Automatic adaptive 3-d finite element refinement using different-order tetrahedral elements[END_REF][START_REF] Lee | Automatic adaptive refinement finite element procedure for 3d stress analysis[END_REF]. Moreover, they showed that it was more efficient to refine the computational mesh first in the plane of the plate then in the thickness. This is the strategy that we will pursue here. A non-conforming 3D quadratic mesh will be considered for both approximation and geometrical representation. In order to be able to describe implicitly the typical geometry of thin structures with Level-Sets, the size of the mesh must be close to the minimal thickness of the part. This geometrical requirement guided the geometrical examples presented in the last section. The objective of the next examples is to assess this strategy on some model problems.

Verifications

Thin structure

Consider a clamped square plate shown in figure 15. The length of the plate is 6mm and its thickness h is set to 0.023438mm so that the Kircchhoff theory can be considered. If the external load is denoted f , the deflection at the centre of the plate can be approximated by the following expression:

u c = 0.00126 f L 4 D ( 2
)
Where D is called 'plate stiffness': First, consider the case where the geometry is exact (the Level-Set is taken analytically), so that only the approximation induces errors. In addition, the mesh is aligned with the plate. A reference shell computation (18726 dofs) involving 60 × 60 abaqus S4R51 elements is considered for comparison with the analytical solution and to define a reference in-plane shear stress field.

D = Eh 3 12(1 -ν 2 ) ( 3
The discrepancy between this shell solution and the analytical one is 0.433% (u c = 1.391 10 6 mm in this case). The amplitude of the reference in-plane shear stress is ±3356MPa (see figure 16). The results are summarized in table 3: It can be seen that the mesh with a unit aspect ratio induces a 0.31% error on the deflection and 0.02% for the in-plane shear stress. Increasing the aspect ratio to 8 increases moderately the errors, while with a 16 aspect ratio begins to exhibit large errors. If tetrahedral elements are considered, then the influence of the aspect ratio is much more pronounced. This illustrates the fact that tetrahedral elements are more sensitive to distortions than hexahedral elements.

For X-FEM and using similar elements (see for example figure 19), the error on the deflection is consistent with the finite element case (2.16%), while the error on the shear stress is lower than the finite element counterpart (0.13%). The difference stems from the fact that two half elements are used in the thickness of the plate, so that the approximation space is not exactly the same as the finite element one. However, this validates the fact that using the X-FEM do not modify the performance of the underlying finite element approach. Note also that, as explained in [START_REF] Lee | Automatic adaptive fe analysis of thin-walled structures using 3d solid elements[END_REF], the large increase of the error in the case of a 16 aspect ratio with tetrahedrons is due to an insufficient discretization in the plane of the plate. This small verification shows that thin structures can be accurately simulated with X-FEM using low order solid elements, provided that a sufficient number of elements are generated in the plane. This minimum number of elements is obtained naturally within the proposed strategy, as

Thick structure

The same tests can be applied on a thick plate whose dimensions are 6 × 6 × 0.3mm (aspect ratio=20). The results obtained with different discretizations are presented in table 4 and are compared to a reference solution which was computed using 120 × 120 abaqus S8R6 2 elements (deflection 6.957 10 2 mm, σ M ax 12 = 19.993MPa). In these examples, variations of both aspect ratio and element size were considered. From these results, one can see that a good accuracy can be obtained using solid elements of size comparable to the thickness of the plate (less than one percent of discrepancy with respect to the reference solution). Obviously, using smaller elements improves the accuracy. Using flattened elements increases the error above to about two percent, but this also stems from the fact that the number of dofs has dropped. As in the thin case, tetrahedral elements 16 aspect ratio 16 aspect ratio elements) In the X-FEM case (see figure 20 for a typical mesh), the conclusions are similar than in the case of thin plates (see Table 6): the use of the X-FEM do not modify the behaviour of the elements.

Nevertheless, some discrepancies can appear, due to the fact that the approximation space obtained after cutting the elements is different from the one obtained with a conforming mesh.

Application of the proposed approach

Now, the complete strategy proposed above is applied with the example of the clamped plate.

Both thick and thin cases are considered. In the two cases, a CAD file representing the plate is defined. Then, this file is taken as an input of the mesh-generation routine. The boundary conditions are imposed by means of a penalty method. Note that Lagrange multiplier method [START_REF] Moës | Imposing Dirichlet boundary conditions in the extended finite element method[END_REF][START_REF] Béchet | A stable lagrange multiplier space for the stiff interface conditions within the extended finite element method[END_REF] or Nitsche's method [START_REF] Nitsche | Über ein Variationprinzip zur lösung von Dirichlet-Problem bei Verwendung von Teilräumen, die keinen Randbedingungen unterworfen sind[END_REF] could be equally used. As the process tends to smooth geometrical angles, an approximation is introduced in the boundary conditions. In practice, the elements of the iso-zero are selected for the boundary conditions when the distance of their centroid to a given plane 19: Clamped plate subjected to a uniform pressure: typical X-FEM mesh (aspect ratio: 16).

Aspect ratio Type of elements Deflection error (%) In-plane shear stress error (%) is lower than a prescribed threshold. The higher this threshold, the higher the number of elements whose displacement is prescribed, and thus the stiffer the behaviour of the plate in this particular case. This is illustrated in figure 21 where the actual boundaries of the plate are compared with the CAD one. The deflection results are given in table 5. One can see the influence of the selection of the boundaries: smaller boundaries results in a higher deflection at the centre of the plate. Depending on this parameter, the errors in deflection and shear stress vary respectively from 1.28% to 2.16% and 0.86% to 1.1%. Moreover, it can be seen that a good accuracy can be obtained in this case, with elements whose aspect ratio is moderate. The case of the thick plate (thickness 0.3mm) is now considered: The base mesh is composed of a 0.15 × 0.15 × 0.15mm quadratic elements level 0 mesh.

In order to represent accurately the boundaries of the domain, 1, 2 and 3 levels are added: Figure 22 depicts the influence of adding multiple levels to the initial mesh. The deflection and shear stress errors are given in table 6. One can see that although the deflection error is decreasing, the shear error do not decrease monotonously upon mesh refinement. This is due to the reference that is used: Aspect ratio Elts/thickness Type of elements Deflec. error (%) In-plane shear stress error (%) it was obtained using shell elements that may not represent the three dimensional behaviour of the plate. An overkill FEM solution is computed in order to build a reference (0.06 × 0.06 × 0.06mm quadratic elements). The figures inside the brackets show the errors (in both deflection and stress) with respect to this new reference: It can be seen now that the error converges to zero when the number of levels increase, which is consistent.

Finally, we still consider the same example, but now the plate is rotated by an angle of 60 • in the xy plane in order to study the influence of the orientation of the part with respect to the mesh. 1, 2, 3 and 4 levels are considered, level 0 mesh being composed of 0.15 × 0.15 × 0.15mm quadratic elements. The corresponding meshes and boundaries are depicted in figure 23, and the improvement on the geometry is clearly visible. Deflection errors are presented in table 7 (stress are not depicted, as the plate is not oriented along x and y axis). Like in the previous example, the improvement in term of error is clear when increasing the number of mesh levels. It is also clear that the rotation of the plate almost did not modify the error level (with a monotonous decrease of the error with respect to the overkill FEM solution).

Numerical examples

The proposed approach now illustrated on some idealized CAD parts provided by EADS IW.

These structures are representative of substructures that are commonly found in the aerospace industry. Of increasing complexity, the objective here is to illustrate the potential of the approach. The results are presented in figures 25 , 26 and 27 for the proposed cases. One can see that despite of the fact that the level 0 mesh in figures 26 and 27 is twice coarser than the mesh in figure 25, the mechanical fields (both displacement and Von-Mises norm of the stress) are similar and extremely smooth. Adding mesh levels allows similar or better geometrical accuracy depending on the number of levels involved, but with a smaller amount of dofs than homogeneous refinement.

Quantitative results

In order to deal with quantitative results, the error in the strain energy is computed from the previous results. The reference energy is obtained by means of an overkill conforming finite element solution with a fourth order polynomial approximation (594 060 dofs). The results are summarized in Table 8. In particular, it can be seen that the use of multiple mesh levels greatly improves the accuracy and the efficiency of the approach: the coarse mesh L x = L y = L z = 1.2mm with two levels gives a lower error and fewer dofs than the fine mesh (L x = L y = L z = 0.6mm) with only one level.

The convergence of the approach is finally assessed by means of homogeneous mesh refinement.

Meshes of size L x = L y = L z = 0.8, 0.4, 0.3 and 0.25mm are considered, and the error in the strain energy is computed with respect to the reference solution. The results are presented in figure 28. It can be seen that the convergence rate is close to the optimal one (i.e. 2 as a quadratic approximation is used). Theoretically, the convergence rate should be bounded by O(h 3/2 ) (see number of points depend on both mesh and geometry. The advantage of such an approach stems from the fact that the burden of creating the sub-elements is avoided. A better alternative would consist in the use of this incompatible integration strategy together with the exact CAD description of the geometry (like with the Finite Cell [START_REF] Duster | pq-adaptive solid finite elements for three-dimensional plates and shells[END_REF]).

(c) (d)

Glued junction

The case of the glued junction is now considered. Here, only one part of the junction is presented in order to run on a personal computer.

Conclusion

In this contribution, an integrated approach was proposed in order to simulate the behaviour of aerospace CAD-based thin structures. The proposed strategy relies on the combined use of the eXtended Finite Element method and the Level-Set approach. As the method is non-conforming, there is no need to explicitly mesh the volume of the domain of interest: the computation is done on a regular grid of elements which also supports the geometrical description (thanks to the Level-Set approach). An algorithm was proposed in order to convert complex CAD parts into a Level-Set and a mesh support. The use of a nested set of grids makes possible the description of sharp edges with Level-Sets. Although less accurate than the algorithm proposed by Moumnassi et al. [START_REF] Moumnassi | Finite element analysis on implicitly defined domains: An accurate representation based on arbitrary parametric surfaces[END_REF], our implementation is tolerant to damaged CAD geometries (with overlap or holes).

Note however that the parallelization of the algorithm would improve its efficiency. The use of solid elements was also proposed, as it avoids the construction of mid-plane surfaces, and also makes natural the modelization of junctions between thin and massive parts. The approach was applied to some academic examples, then to representative CAD parts. It was shown that a very good accuracy could be obtained. The efficiency of the approach could be improved by means of the parallelization of the Level-Set computation. Current work on this topic focus now on the use of higher order approximations together with the separation between geometry and mechanical approximation using the approaches developed in [START_REF] Legrain | High order eXtended Finite Element Method and Level Sets: Uncoupling geometry and approximation[END_REF][START_REF] Legrain | High order X-FEM and levelsets for complex microstructures: Uncoupling geometry and approximation[END_REF] (see the illustrative example shown below), or [START_REF] Legrain | A NURBS Enhanced eXtended Finite Element Approach for Unfitted CAD Analysis[END_REF]. It could also be applied to geometries defined by means of surface scanners as a tool for reverse engineering.

Outlook

The improvement of the approach proposed in the conclusion is illustrated here. In order to save degrees of freedom with a fixed geometrical accuracy, an approach derived from [START_REF] Legrain | High order X-FEM and levelsets for complex microstructures: Uncoupling geometry and approximation[END_REF] is considered.

The geometry is represented on the fine grid depicted in figure 33, and the approximation is defined on the mesh shown in figure 34. The problem is subjected to the same BCs as in section 5.1. A P6 approximation is considered, leading to 106 953 dofs (which is far less than in section 5.1). The resulting stress and displacement fields are given in figure 35 (note that the stress field was not smoothed at all). The corresponding error in the strain energy is 0.007% which couldn't be attained using the approach advocated in this paper, unless a million of dofs are used.

Figure 1 :

 1 Figure 1: CAD to Level-Set process.

Figure 2 :

 2 Figure 2: Illustration of the smoothing of the sharp geometrical features: (a) Mesh and initial geometry (grey); (b) mesh and Level-Set geometry (grey); (c) Multi-Level mesh and Level-Set geometry (grey).

c)

  Compute the bounding-box of the STL discretization and enlarge it. d) Activate level 0 cells containing at least one point of the surface cloud (i.e. crossed by the interface).

(see figure 4 )

 4 were provided by EADS IW, and are representative of substructures that are commonly found in the aerospace industry. The first one (figure 4(a)) represents a plate containing a hole at its centre, the second one (figure 4(b)) a glued junction and the third one (figure 4(c)) a stiffened panel.

Figure 3 :Figure 4 :

 34 Figure 3: Illustration of the algorithm for constructing the level set: a) Initial geometry; b) face points cloud;c) enlargement of the bounding-box; d) activation of the cells containing "surface points"; e) activations of cells sufficiently close to the surface; f) level 0 mesh; g) activation of level > 0 cells (in red, invalid son cells); h) Level-Set mesh; i) computation of the Level-Set at the nodes of the mesh; j) final mesh + Level-Set (after cancelling cells out of the domain)

Figure 5 :

 5 Figure 5: Perforated plate: mesh, Lx = Ly = Lz = 1.2mm. (a), (b) and (c): 0, 1, 2 sub-grids levels. (a) 6200 nodes 23500 elements (0.561s); (b) 13000 nodes 49700 elements (0.885s); (c) 27630 nodes 97768 elements (1.553s)

Figure 6 :

 6 Figure 6: Perforated plate: iso-zero, Lx = Ly = Lz = 1.2mm. From left to right 0, 1, 2 sub-grids levels. (a) 6200 nodes 23500 elements (0.561s); (b) 13000 nodes 49700 elements (0.885s); (c) 27630 nodes 97768 elements (1.553s)

Figure 7 :

 7 Figure 7: Perforated plate: mesh, Lx = Ly = Lz = 0.6mm. From left to right 0, 1, 2 sub-grids levels. (a) 38400 nodes 179000 elements (6.747s); (b) 53300 nodes 227200 elements (7.739s); (c) 76500 nodes 292100 elements (10.095s).

Figure 8 :Figure 9 :

 89 Figure 8: Perforated plate: iso-zero, Lx = Ly = Lz = 0.6mm. From left to right 0, 1, 2 sub-grids levels. (a) 38400 nodes 179000 elements (6.747s); (b) 53300 nodes 227200 elements (7.739s); (c) 76500 nodes 292100 elements (10.095s).

Figure 10 :Figure 11 :

 1011 Figure 10: Glued junction: iso-zero, Lx = Ly = Lz = 6mm. From left to right 0, 1, 2 sub-grids levels. (a) 112000 nodes 538000 elements (99s); (b) 187500 nodes 851000 elements (111s); (c) 399000 nodes 1 745 980 elements (160s)

Figure 12 :Figure 13 :

 1213 Figure 12: Stiffened panel: level set (Lx = Ly = Lz = 1.5mm), 4 492 189 nodes 8 400 768 elements, 6140s

  end end Step c) Compute the bounding-box of the cloud of nodes Enlarge the bounding-box Compute N x , N y et N z (number of elements along x, y and z) thanks to the knowledge of L x , L y , L z and the bounding-box Create a regular N x × N y × N z grid of cubic elements Step d), e), f ) foreach Point p on the faces point cloud do Activate level 0 cells containing p or sufficiently close to p (max distance ±R max along x, y et z) end Step g) for Level CurrentLevel=1 to Levels do foreach Point p on the edges point cloud do Activate level CurrentLevel cells containing p or sufficiently close to p (max distance ±R tube /2 Levels-CurrentLevel-1 along x, y et z) end end Step h) Cancel invalid son cells (see Algorithm 2)

Figure 14 :

 14 Figure 14: Stiffened panel: level set (Lx = Ly = 1.5mm, Lz = 4mm), 688000 nodes 3 141 000 elements, 2460s

Figure 15 :

 15 Figure 15: Clamped plate subjected to a uniform pressure

Figure 16 :

 16 Figure 16: Clamped plate subjected to a uniform pressure: In-plane shear stress.

Figure 17 :

 17 Figure 17: Clamped plate subjected to a uniform pressure: solid meshes

Figure 18 :

 18 Figure 18: Clamped plate subjected to a uniform pressure: solid meshes

Figure 20 :

 20 Figure 20: Clamped plate subjected to uniform pressure: Typical X-FEM mesh (aspect ratio: 1)

Figure 21 :

 21 Figure 21: Influence of the threshold used to select the boundaries. (a) Large threshold; (b) Small threshold

Figure 23 :

 23 Figure 23: Influence of the number of levels on the representation of the boundaries (rotated case).

5. 1 .Figure 24 :

 124 Figure 24: Geometry and boundary conditions.

Figure 25 :

 25 Figure 25: Drilled plate in traction and bending (Lx = Ly = Lz = 0.6mm elements, one mesh level) : (a), (b), (c) : displacement ; (d) stress (Von Mises norm)

Figure 27 :

 27 Figure 27: Drilled plate in traction and bending (Lx = Ly = Lz = 1.2mm elements, three mesh levels): (a), (b): displacement ; (c), (d): stress (Von Mises norm)

Figure 29 :

 29 Figure 29: Boundary conditions for the junction

Figure 30 :Figure 31 :

 3031 Figure 30: Glued junction: results for a two levels mesh (Level 0: Lx = Ly = Lz = 6mm)

Figure 33 :

 33 Figure 33: Geometrical grid and iso-surface

Figure 34 :

 34 Figure 34: Approximation mesh and iso-surface

Figure 35 :

 35 Figure 35: Displacement and stress fields. Note that the stress field was not smoothed at all.

  [6] T.J.R. Hughes, M. Cohen, and M. Haroun. Reduced and selective integration techniques in finite element analysis of plates. Nuclear Engineering Design, 46:203-222, 1978.[7] T.J.R. Hughes. Generalization of selective integration procedures to anisotropic and nonlinear media. International Journal for Numerical Methods in Engineering, 15:1413-1418, 1980.[8] H. Ben-Dhia and G. Rateau. The arlequin method as a flexible engineering design tool. International Journal for Numerical Methods in Engineering, 62:1442-1462, 2005.[9] N. Moës, J. Dolbow, and T. Belytschko. A finite element method for crack growth without remeshing. International Journal for Numerical Methods in Engineering, 46:131-150, 1999.

:

  Hierarchy of cartesian grids for CurrentLevel=1 to Levels do foreach Active cell c in Grid of level CurrentLevel do Get sons of c in Grid of level CurrentLevel+1 if At least one son is active and not All son active then (Only consider partially filled cells)

	Deactivate active sons;
	end
	end
	end
	Output: Cartesian grids without invalid sons
	Level 0 elements	Levels Volume error (%)

Table 1 :

 1 Volume error, perforated plate

	Level 0 elements	Levels Volume error (%)

Table 2 :

 2 Volume error, glued junction

Table 3 :

 3 Influence of mesh density and aspect ratio on the error (thin plate)

Table 4 :

 4 Errors depending on mesh density and aspect ratio (thick plate)

Table 5 :

 5 Influence of mesh density and aspect ratio (thin plate) on the error, using the proposed approach. The two rows show the influence of the selection of the boundaries.

	Asp. ratio Elts/thick. Lev. Deflec. err. (%) [/FEM] Shear err. (%) [/FEM]
	1	2	1	8.857 [8.3568]	0.5247 [0.1481]
	1	2	2	2.230 [1.6929]	1.2409 [0.9531]
	1	2	3	0.0985 [0.45021]	1.2079 [0.9202]

Table 6 :

 6 Influence of mesh density and aspect ratio (thick plate) on the error, using the proposed approach. Influence of the number of sub-grids. Figure 22: Influence of the number of levels on the representation of the boundaries.

			Y Z	X	Y Z	X
		level		2 levels
				Y Z	X
			3 levels
	Aspect ratio Elts/thickness Levels Deflec. error (%) Deflec. error /FEM(%)
	1	2	1	5.741	5.2233
	1	2	2	2.644	2.1096
	1	2		1.249	0.9744
	1	2	4	1.249	0.7069

Table 7 :

 7 Influence of mesh density and aspect ratio (thick plate) on the error, using the proposed approach. Influence of the number of sub-grids for the rotated case.

4-node thin shell, reduced integration with hourglass control, using five degrees of freedom per node (two rotations and three displacements)

8 nodes doubly curved thick shell with reduced integration, using six degrees of freedom per node (three rotations and three displacements).

[43, p 119] and [START_REF] Legrain | Stability of incompressible formulations enriched with X-FEM[END_REF]). However, it was shown in [START_REF] Dréau | Studied X-FEM enrichment to handle material interfaces with higher order finite element[END_REF] that optimal convergence could nevertheless be obtained in practice, depending on the problem (which seems to be the case here). The influence of the integration is now assessed: rather than partitioning the elements along the iso-zero of the Level-Set, a classical integration rule is considered. In order to be able to converge, the number of integration points is increased in the elements that are cut by the interface: a 72 points degenerated quadrangle integration rule is considered in this case: the convergence of this approach is depicted in figure 28: it can be seen that a loss of accuracy can occur depending on the mesh density and most likely on the geometry of the problem. Moreover, increasing the number of integration points on such fine meshes has a huge impact on the computational time. Note also that if the number of integration points is insufficient, then the stiffness matrix can become singular. This minimum Mesh Levels dofs Error (%) Strain energy error 2 X-FEM -Quadratic X-FEM -Quadratic incomp. integ.

Figure 28: h-convergence of the proposed approach. It can be seen that optimal rate is obtained when the integration is done properly. On the contrary, incompatible integration can lead to a loss of convergence if the number of integration points is not severely increased.

A two level mesh with L x = L y = L z = 6mm for level 0 elements is considered. The mesh is composed of 249 816 quadratic tetrahedra (1 051 404 dofs). Thanks to the use of solid elements, the junction between the different parts is taken into account without any hypothesis on the geometry of the different sub-parts. The displacement and stress fields are displayed in figure 30. Stress concentrations are captured by the solution at the centre of the junction. These informations could not be obtained using a plate model.

Stiffened plate

Finally, the stiffened plate presented in section 2.2 is considered. Like in the previous example, only one part of the plate is considered in order to allow its treatment on a personal computer. The