Thomas-David Jayet

Emmanuel Baranger

Guillaume Couegnat

Sébastien Denneulin

Feasability of a weakly intrusive Generalized Finite Element Method implementation in a commercial code: application to Ceramic Matrix Composite micro-structures

Keywords: GFEM, Handbook problem, Implementation, Composites

In this contribution is addressed the issue of modeling a SiC/SiC composite micro-structure at the fiber scale with elementary features. As the development of new grades of Ceramic Matrix Composites (CMC) for civil aviation grows, different manufacturing processes are used successively and lead to different types of micro-structures and a variable material quality. Consequently a versa-

tile model should be developed in order to compare these materials and create a tool to help engineers to predict the mechanical behavior at the fiber scale.

Here the Generalized Finite Element Method (GFEM) is proposed which consists in enriching the classical Finite Element (FE) approached displacement by numerical functions to deliver an accurate description of the fiber-scale structure while limiting the number of degrees of freedom compared to a classical finite element description. A pattern-based description of the microscale is depicted using an industrial code for an engineering purpose. Four main difficulties are highlighted (i) the choice of the enrichment functions regarding the literature (ii) their stiffness matrix computation in a commercial code (iii) the construction of the pattern-based structure and (iv) the post-processing. Two GFEM strategies are presented and demonstrate the feasibility of an enriched kinematics within a classical finite element modeler. The selection of such modeler

Introduction

Ceramic Matrix Composites are good candidates for the aeronautic or the nuclear industries due to their high refractoriness, their low density and their resistance to oxidation compared to metals. In the specific application of gas turbine engines for aviation propulsion, the research of high temperature performance is motivated by a desire of reducing fuel consumption. Indeed a weight loss is not the only motivation for integrating CMCs within hot parts: elevating the operating temperatures of the gas turbine engines enables a gain in propulsion efficiency [START_REF] Zok | Ceramic-matrix composites enable revolutionary gains in turbine engine efficiency[END_REF]. Improvements in temperature capabilities allowed an increasing work under service conditions of approximately 35 o C per decade [START_REF] Corman | Development history of GE's prepreg melt infiltrated ceramic matrix composite material and applications[END_REF] in the past.

CMC research has led to the development of various manufacturing processes that have resulted to the creation of several types of structures with different thermal and mechanical properties. For chosen components, two factors have an influence on the thermomechanical behavior: the architecture of the textile reinforcement with unidirectional plies [START_REF] Pryce | Behaviour of unidirectional and crossply ceramic matrix composites under quasi-static tensile loading[END_REF], 2D or 3D fabrics [START_REF] Yamada | Mechanical and thermal properties of 2D and 3D SiC/SiC composites[END_REF][START_REF] Boisse | Textile reinforcements: Architectures, mechanical behavior, and forming[END_REF] and the matrix deposition process. Typically four types of processes can be used for the matrix deposition: the Chemical Vapor Infiltration (CVI), the Chemical Vapor Deposition (CVD), the Melt Infiltration (MI) and the Polymer Impregnation and Pyrolisis (PIP). The chemical vapor processes are employed for a thicknesscontroled preform coating [START_REF] Yang | Preparation and mechanical properties of carbon fiber reinforced (bc x-SiC) n multilayered matrix composites[END_REF] at the cost of a high porosity rate [START_REF] Brennan | Interfacial characterization of a slurry-cast melt-infiltrated SiC/SiC ceramic-matrix composite[END_REF] and a long deposit time. They also can be used to lay an interphase coating such as Pyro-Carbon (PyC) or Boron Nitride (BN) around the fibers. MI and PIP come as a complement to fill the residual porosities (see Fig. 1) and accelerate the matrix deposition. A review of the different manufacturing processes is available in [START_REF] Dicarlo | Advances in SiC/SiC Composites for Aero-Propulsion[END_REF].

In particular the matrix deposition processes lead to different micro-structures. Furthermore many studies have been conducted at the fiber scale in order to quantify the first damages that lead to the material ruin. Degradation mechanisms like fiber debonding [START_REF] Evans | The physics and mechanics of fibre-reinforced brittle matrix composites[END_REF], fiber cracking due to oxidation [START_REF] Evans | Models of hightemperature, environmentally assisted embrittlement in ceramic-matrix composites[END_REF] and creep [START_REF] Abbé | Flexural creep behavior of a 2D SiC/SiC composite[END_REF][START_REF] Carrère | Creep behaviour of a SiC/Si-BC composite with a self-healing multilayered matrix[END_REF] have been observed. Silicon Carbide (SiC) textile and matrix (SiC/SiC) composites behavior under quasi-static indentation [START_REF] Herb | Damage assessment of thin SiC/SiC composite plates subjected to quasi-static indentation loading[END_REF] or low velocity impact [START_REF] Herb | Damage analysis of thin 3D-woven SiC/SiC composite under low velocity impact loading[END_REF] have been observed through micrographic or tomographic studies in order to analyze the damages at the yarn scale. Initial state is also of importance since the post-processing cooling implies a residual stress [START_REF] Bobet | Study of thermal residual stresses in ceramic matrix composites[END_REF][START_REF] Mei | Measurement and calculation of thermal residual stress in fiber reinforced ceramic matrix composites[END_REF] due to different values of the expansion coefficient for the different phases of the micro-structure [START_REF] Bobet | Thermal residual stresses in ceramic matrix composites-i. axisymmetrical model and finite element analysis[END_REF].

The succession of possible manufacturing processes coupled with the existence of non linear phenomena such as degradation mechanisms at the microscale raises the issue of a modeling strategy to predict CMCs behavior and compare solutions for different configurations and loads. Models for Ceramic Matrix Composites at the fiber scale have been developed for computational applications like crack initiation and propagation criteria [START_REF] Pineau | Virtual testing applied to transverse multiple cracking of tows in woven ceramic composites[END_REF][START_REF] Mazars | Multi-scale damage modeling of 3D ceramic matrix composites from in-situ X-ray tensile tests[END_REF] or environmental effects modeling [START_REF] Pailler | Micromechanics based model of fatigue/oxidation for ceramic matrix composites[END_REF]. However, the simulation of such phenomena requires the development of calculation tools whose effectiveness depends on the quality of the description of the structure. Thus an adapted strategy for the micro-structure description must be implemented upstream to facilitate the mechanical modeling. For an industrial purpose such modeling method must be implemented in a commercial code so the following developments are chosen to be made with the Abaqus 6.14 software.

The aim of the following paragraphs is to introduce classical strategies used to model composites numerically. A review of the homogenization and domain decomposition strategies is proposed with the purpose of contextualizing the chosen method which will be precised further.

Homogenization consists in averaging quantities such as stresses and strains within a heterogeneous Representative Volume Element (RVE) whose definition can be found in [START_REF] Nemat-Nasser | Overall stresses and strains in solids with microstructure[END_REF]. A satisfying RVE must be statistically representative which implies an upstream study of the volume fraction of the different constituents within the RVE and their allocation. It must have a small characteristic length compared to the macrostructure dimensions so the RVE is assumed to be homogeneous at the macroscale. For non linear investigations applied to periodic RVE, iterative methods based on an asymptotic homogenization [START_REF] Sanchez-Palencia | Homogenization method for the study of composite media[END_REF] to coupled scales (reviewed in [START_REF] Fish | Multiscale analysis of composite materials and structures[END_REF]), FE 2 [START_REF] Feyel | Multiscale FE2 elastoviscoplastic analysis of composite structures[END_REF] or Fast Fourier Transform [START_REF] Gélébart | Non-linear extension of FFT-based methods accelerated by conjugate gradients to evaluate the mechanical behavior of composite materials[END_REF] have been suggested for composite applications. However these homogenization strategies suffer from the simplifying assumption of periodicity and uncoupled global-local methods encounter significant errors in the portions with high gradients. Moreover no solution for the micro-structure description is provided and the information to be given is too rich compared to the information at the end of the operation, then homogenization has a limited interest for this study.

Domain Decomposition Methods (DDM) enable to split a global structure into non-overlapping substructures and interfaces. Consequently they allow a fine description of a large structure due to a parallelized computation of the different submodels. Then computation is made with accuracy. Studies have been conducted for laminate composites submitted to time-dependant mechanisms [START_REF] Ladeveze | Multiscale modelling and computational strategies for composites[END_REF][START_REF] Guo | A domain decomposition approach for static and dynamic analysis of composite laminated curved beam with general elastic restrains[END_REF]. A crack propagation has also been performed using a hybrid eXtended FEM and DDM formulation [START_REF] Guidault | A multiscale extended finite element method for crack propagation[END_REF]. Nonetheless DDM have limitations concerning the sorting of extracted information caused by the accuracy of the computation. Indeed the output information may be too rich to enable an impact on the selection of the manufacturing processes. Consequently the desired versatility for the composite depiction seems lost.

In addition to the Galerkin formulation, other kinematics descriptions are possible to approach the real displacement. The Superposition Finite Element Method (S-FEM) [START_REF] Fish | The s-version of the finite element method[END_REF], the Variable size h and Polynomial degree p called hpadaptivity [START_REF] Duarte | An hp adaptive method using clouds[END_REF] are some examples. The Partition of Unity Method [START_REF] Melenk | The partition of unity finite element method: basic theory and applications[END_REF] can be seen as a generalization of the hp-adaptivity and constitutes the foundation of the Generalized Finite Element Method (GFEM) [START_REF] Strouboulis | The design and analysis of the generalized finite element method[END_REF]. The key point of the GFEM is to enrich the kinematics with enrichment -also called specialfunctions in order to describe a structure and its behavior on each element patch of the mesh.

The GFEM is a method halfway between those of homogenization and those of HPC. Indeed, an upstream calculation allows to generate patterns on a macroscopic scale from a local description before storing the solutions in a Handbook.

The enrichment functions enable to keep a coupling between the micro scale and the macro scale. The problem to be solved is then broken down into patterns but the calculation is done in a monolithic way.

For that reason Section 2 of this contribution emphases the Generalized Finite Element Method as a Partition of Unity Method and underlines the difficulties of the method. Then three GFEM implementation strategies are presented in Section 3. Section 4 draws the statement of the tools available in Abaqus in order to estimate the most relevant strategy to implement. Eventually Section 5 presents a weakly intrusive implementation of the GFEM with the help of Abaqus and Python in the particular case of an elastic problem. The novelty of the presented method do not lies in the numerical enrichments generation or the choice of the patterns but in their use in a commercial software. It will be shown that the feasibility of the implemented method is easily automatable through a scripting language.

Generalized Finite Element Method

The following section draws a state of the art of the GFEM. Strategies of kinematics enrichment generation are summarized from the literature and discussed in order to select the most appropriate for an easy implementation.

Problem definition

A domain denoted Ω is considered. Boundary conditions can be applied on

∂Ω = ∂Ω F ⊕ ∂Ω u as: σ.n = F d on ∂Ω F (1)
u = u d on ∂Ω u (2)
Where n is the outward unit normal vector to ∂Ω F and σ is the Cauchy stress tensor. Such Ω domain is approached by a mesh denoted Ω h . Centered on the node number k can be defined a patch of elements Ω k . An illustration of the model problem definition is available in Fig. 2. The strain is assumed to be small and given by the symmetric gradient of the displacement u:

(u) = 1 2 ∇(u) + ∇(u) T (3)
For the heterogeneous material the Hooke's elasticity tensor at the coordinate x denoted H(x) can be broken down as the sum of a homogeneous part H 0 and a perturbation part ∆H(x) such as:

∀x ∈ Ω, H(x) = H 0 + ∆H(x) (4)
Then the stress tensor σ is linked to the strain tensor by:

σ(x) = H(x) : (u) (5)
Defining U 0 = {u * ∈ H 1 (Ω) 3 , u * = 0 on ∂Ω u } and U = {u ∈ H 1 (Ω) 3 , u = u d on ∂Ω u } enable to write the Principle of Virtual Work which reads: Find Based on the Partition of Unity (PU) hypotheses [START_REF] Melenk | The partition of unity finite element method: basic theory and applications[END_REF], the GFEM formulation suggested in [START_REF] Strouboulis | The design and analysis of the generalized finite element method[END_REF] enables to write the displacement approximation as the sum of the regular Finite Element (FE) interpolation and an enriched contribution. In a point x of the studied domain, such a displacement u h (x) is written as presented in [START_REF] Brennan | Interfacial characterization of a slurry-cast melt-infiltrated SiC/SiC ceramic-matrix composite[END_REF].

u ∈ U such that ∀u * ∈ U 0 , - Ω (u) : H(x) : (u * )dΩ + ∂Ω F d .u * d(∂Ω) = 0 (6)

Partition of unity

Ω h Ω Ω p Ω q Ω k ∂Ω F ∂Ω u u d F d k
u h (x) = N k=1 ϕ k (x)u k FE interpolation + N k=1 ϕ k (x) n (k) enr j=1 ψ (j) k (x)a (j) k enriched contribution (7)
N is the number of element patches inside the mesh representing the studied structure, ϕ k is the same "hat" function as defined in the classic Finite Element Method (FEM) for the element patch number k called Ω k and illustrated in Fig. 2, u k is the displacement vector calculated at the node number k, n (k) enr is the number of enrichment functions inside the patch number k, ψ (j) k is the enrichment function number j associated with the patch Ω k and a (j) k is the additionnal nodal unknown associated with the enrichment number j of Ω k . It has been shown in [START_REF] Melenk | The partition of unity finite element method: basic theory and applications[END_REF] that the richer the enrichment functions, the closer the approached displacement u h is to the exact displacement u. The immediate consequence is the emancipation of an explicit mesh due to the presence of features within the micro-structure, thanks to the enrichment functions ψ

(j) k . The ϕ k ψ (j)
k products ensure the continuity of the fields between two superimposed or adjacent patches.

XFEM and GFEM

A distinction between the Extended Finite Element Method (XFEM) and the GFEM is made in this contribution: despite both of the two strategies are consequences of the PU hypotheses, the first method will designate the extension of the standard FE with analytical enrichments and the second one will refer to the use of numerical special functions.

The XFEM has already been integrated into commercial softwares and is principaly used to model the position and the propagation of cracks by using level-set functions and adding a set of tip functions at the crack fronts explicited in [START_REF] Fleming | Enriched element-free galerkin methods for crack tip fields[END_REF]. These tip functions are defined from the asymptotic displacement solution and determine the crack deflection. A modified Partition of Unity using Shepard analytical functions [START_REF] Shepard | A two-dimensional interpolation function for irregularly-spaced data[END_REF] was also suggested to describe a crack in a homogeneous field [START_REF] Duarte | A generalized finite element method for the simulation of three-dimensional dynamic crack propagation[END_REF] or to treat problems with multiple junctions in polycristalline structures under thermal loads [START_REF] Aragón | Generalized finite element enrichment functions for discontinuous gradient fields[END_REF]. It can be noted that the XFEM has been used in Abaqus to model crack deflection at the interface between the matrix and the fiber interphase [START_REF] Braginsky | Simulation of crack propagation/deflection in ceramic matrix continuous fiber reinforced composites with weak interphase via the extended finite element method[END_REF]. A homemade code also made it possible to perform crack propagation in a homogenized composite structure via XFEM using the abaqus modeler [START_REF] Talebi | A computational library for multiscale modeling of material failure[END_REF]. Many studies [START_REF] Zi | New crack-tip elements for xfem and applications to cohesive cracks[END_REF][START_REF] Menouillard | Efficient explicit time stepping for the extended finite element method (X-FEM)[END_REF][START_REF] Fries | A corrected XFEM approximation without problems in blending elements[END_REF][START_REF] Bordas | On the performance of strain smoothing for quadratic and enriched finite element approximations (XFEM/GFEM/PUFEM)[END_REF][START_REF] Wu | An improved stable XFEM (Is-XFEM) with a novel enrichment function for the computational modeling of cohesive cracks[END_REF] have been lead in order to improve the implementation and the stability of the method. Nonetheless two limitations of the analytic description of a crack can be highlighted.

First the asymptotic displacements used as crack tip enrichment functions are generally valid in an infinite homogeneous domain. Second level-set functions give a rough approximation of a crack path if a coarse mesh is used (see Fig. 3), then element partitioning is required since usual Gauss quadrature is only accurate for polynomials. Moreover a multi-material domain demands a fine mesh in order to take into account the structural details which can be small within CMCs micro-structures (such as residual porosities for instance). The XFEM is able to describe a weak discontinuity as a frontier between two materials but it would demand a large number of additional nodal unknowns to describe the whole micro-structure. The GFEM seems to be more flexible due to the numerical enrichment functions, for instance a 3D model of a crack has been treated in [START_REF] Duarte | Analysis and applications of a generalized finite element method with global-local enrichment functions[END_REF] avoiding asymptotic expansion of the elasticity solution required in XFEM through the help of a local description of the crack by a geometrically conforming refined mesh. The number of structural details delivers mesh complication as noted in [START_REF] Strouboulis | The generalized finite element method[END_REF] and GFEM can potentially overcome such issue.

The difficulty of treating numerous structural details were also overcome [START_REF] Strouboulis | Generalized finite element method using mesh-based handbooks: application to problems in domains with many voids[END_REF].

Proper Orthogonal Decomposition (POD) modes were used as GFEM enrichment functions [START_REF] Aquino | Generalized finite element method using proper orthogonal decomposition[END_REF]. More recently a Proper Generalized Decomposition (PGD) enrichment generation strategy adapted to transient problems has been developped [START_REF] Canales | Vademecum-based GFEM (V-GFEM): optimal enrichment for transient problems[END_REF] and allows an off-line storage of the enrichment functions which can be seen as a Handbook problem defined in [START_REF] Strouboulis | The generalized finite element method[END_REF][START_REF] Strouboulis | Generalized finite element method using mesh-based handbooks: application to problems in domains with many voids[END_REF]. The Handbook classifies pre-calculated patterns used as numerical enrichment functions like a library of features. The main idea of this method is the richer the Handbook, the more Coarse mesh scenarios and structures can be described.

XFEM crack path

Real crack path

Intersection between mesh edges and crack

Before introducing the different implementation strategies existing in the literature, a review of the difficulties associated with the enrichment functions generation is proposed in the following section.

Difficulties associated with the GFEM and its implementation

The aim of this paper is to present the implementation of a versatile model through the Abaqus software. Such commercial code uses a classical finite element description of the kinematics so the first difficulty is the next:

1. The approached displacement (7) must be imposed on the mesh Ω h using classical kinematic description tools.

Another point of the implementation of the GFEM is to build the stiffness matrix with regard to the displacement formulation. The virtual work of the internal forces Ω (u h ) : H(x) : (u * h )dΩ with u * h ∈ U 0 provides the form of the stiffness matrix [START_REF] Dicarlo | Advances in SiC/SiC Composites for Aero-Propulsion[END_REF] and its components [START_REF] Evans | The physics and mechanics of fibre-reinforced brittle matrix composites[END_REF].

K GF EM =   K F EM K F EM/enr K T F EM/enr K enr   (8) (K F EM ) IJ = Ω I ∩Ω J H(x)∇(ϕ I )∇(ϕ J )dΩ (9a) (K F EM/enr ) Iα = Ω I ∩Ωp H(x)∇(ϕ I )∇(ϕ p ψ (k) p )dΩ (9b) (K enr ) αβ = Ωp∩Ωq H(x)∇(ϕ p ψ (k) p )∇(ϕ q ψ (l) q )dΩ (9c)
Indices α and β are explicited in [START_REF] Evans | Models of hightemperature, environmentally assisted embrittlement in ceramic-matrix composites[END_REF] and a reminder of the development used to obtain these different terms is available in Appendix 1.

α = p-1 m=1 n (m) enr + k ; β = q-1 n=1 n (n) enr + l (10)
The question of the GFEM stiffness matrix conditioning has been raised [START_REF] Babuška | Stable generalized finite element method (SGFEM)[END_REF]. If R(A) denotes a scaled condition number of a matrix A and h << 1 is the characteristic length of the mesh then ( 11) is verified and highlights the ill-conditioned GFEM stiffness matrix compared to the classical FEM stiffness matrix.

R(K

F EM ) R(K GF EM ) = O(h 2 ) ( 11
)
As a consequence convergence of an iterative solution method may be deteriorated. Modified enrichment functions have been used in order to correct the ill-conditioning of the GFEM stiffness matrix and reach the same order than R(K F EM ). This result is obtained by deleting the linear part of the enrichment functions from a patch [START_REF] Babuška | Stable generalized finite element method (SGFEM)[END_REF][START_REF] Gupta | A stable and optimally convergent generalized FEM (SGFEM) for linear elastic fracture mechanics[END_REF]. An illustration is shown in Fig. 4, the modified enrichment function ψ

(j)
k enables to obtain a condition number of the same order of magnitude as for the FEM one.

The main difficulty associated with the GFEM stiffness matrix is to evaluate the integration of the cross-products ∇(ϕ

I )∇(ϕ p ψ (k) p ) and mainly ∇(ϕ p ψ (k) p )∇(ϕ q ψ (l)
q ) between superimposed patches Ω I ∩ Ω p and Ω p ∩ Ω q respectivelly. Indeed hat functions are associated to the coarse mesh Ω h while special functions ψ an integration issue may occur due to non-conforming meshes. Then a refined grid is necessary to project fields or the developpement of an adaptive quadrature can be suggested [START_REF] Strouboulis | The generalized finite element method[END_REF]. However, such an operation requires additional calculations to be performed that require the permanent rewriting of files in a commercial code. The calculation time would potentially become too long, especially in order to deal with non-linear problems in the future. The field projection on a commercial code like Abaqus would therefore undermine the effectiveness of the method proposed in this paper. Consequently the second difficulty to overcome is the next: 2. The integration of cross-products ∇(ϕ p ψ

(k) p )∇(ϕ q ψ (l)
q ) in ( 9) must be simplified.

With the aim of responding to the two difficulties highlighted in this section a review of the different enrichment strategies is suggested in Section 3. Then an analysis of the Abaqus capacities will be presented in Section 4 in order to help to choose a relevant implementation method using this modeler and weakly intrusive languages such as Python for instance.

Enrichment functions generation strategies

This section draws the state of the art of the existing numerical enrichment strategies available in the literature in a chronological order.

Enrichment by precalculated solutions

A natural idea lies in the generation of enrichment functions coming from a combination of canonic loads illustrated in Fig. 5 and introduced by Strouboulis in [START_REF] Strouboulis | The generalized finite element method[END_REF]. These canonic loads allow these patterns to accommodate a combination of the precalculated loads in an elastic problem. Such precalculated solutions are stored in a Handbook and overlaid to a coarse mesh representing the geometry of the homogeneous part of the studied structure. A Computer-Assisted Design (CAD) representation is required so that the geometry of the Handbook problems can be extracted and meshed. Nonetheless the nonconformity between superimposed meshes implies a step of field interpolation and the difficulties in the calculus of the cross-products introduced in (9). A heavy adaptive procedure has been developped [START_REF] Strouboulis | The generalized finite element method[END_REF] to solve that non-conforming problem.

Enrichment by geometrically conforming meshes

Each macroscopic patch which contains a structural detail such as cracks, porosities, fibers or other observable phases of the composite micro-structure must be enriched. The main idea of this strategy later called Duarte's method is to impose the GFEM kinematics to the macroscopic nodes of a patch knowing the local behavior within. To simplify the cross-products for superimposed patches a refined mesh containing the structural details is added to the entire coarse global mesh in order to ensure the geometric conformity between the patches then to guarantee the same local description at the intersection between

Ω h Ω Ω p Ω q ∂Ω u u d ∂Ω F F d σ 11 σ 22 σ 12
Handbook k associated with the "hat" function ϕ k is computed following at a point x :

ψ (j) k (x) = u (j) tot (x) -u (j) hom (x) (12) ψ (j)
k (x) is the difference between the displacement u (j) tot (x) obtained for the real behavior and the displacement u

Enrichment by scale separation

In addition to the partition of unity hypotheses, scale separation hypotheses have been introduced in [START_REF] Tranquart | A computational approach for handling complex composite microstructures[END_REF] and withdrawn below:

Ω h Ω Ω p Ω q ∂Ω F ∂Ω u u d F d u (q) d u (p) d ψ (j2) p ψ (j1) q
Enrichment generation Handbook Fig. 6. Illustration of the different steps to fill the Handbook. First a simulation is performed on a macroscopic coarse mesh (at the left). Secondly a geometrically conforming local mesh is generated in the whole structure. Then patches are extracted to compute the special functions before their reinjection in the model. The first assumption can be interpreted as the hat functions ϕ k are almost constant on the support ω

(j) k , so if C (j) k is the center of ω (j)
k it can be written:

ϕ k (x) ≈ ϕ k (C (j) k ) ∀x ∈ ω (j) k (13)
The second assumption results in an empty intersection between two different enrichment supports. Since the enrichment functions are located on a support ω (j) k with a center C (j) k a modified then by setting ψ [START_REF] Brennan | Interfacial characterization of a slurry-cast melt-infiltrated SiC/SiC ceramic-matrix composite[END_REF] can be simplified and written as shown in [START_REF] Herb | Damage analysis of thin 3D-woven SiC/SiC composite under low velocity impact loading[END_REF].

(j) k = ϕ k (C (j) k )ψ (j) k the approached displacement
u h (x) = N k=1 ϕ k (x)u k + N k=1 n (k) enr j=1 ψ (j) k (x)a (j) k (14)
It can be noted that the enrichment contribution is independent from the hat functions. Consequently the additional nodes a Moreover the terms of the associated stiffness matrix can be calculated similarly to the demonstration of the GFEM stiffness matrix terms (9) and give:

(K F EM ) IJ = Ω I ∩Ω J H(x)∇(ϕ I )∇(ϕ J )dΩ (15a) (K F EM/enr ) Iγ = ω (k) p H(x)∇(ϕ I )∇( ψ (k) p )dΩ (15b) (K enr ) αβ = ω (k) p H(x)∇( ψ (k) p )∇( ψ (k) p )dΩ (15c)
The issue of cross-products highlighted in ( 9) is then avoided, more precisely the scale separation hypotheses allow the superimposition of nonconforming meshes used as support of features. The enrichment functions generation strategy using load base described ahead is brought more versatile, at the price of an approached kinematics which does not fully respect the Partition of Unity hypotheses. Nonetheless such simplification is considered sufficient through a development of a comprehension tool [START_REF] Tranquart | A computational approach for handling complex composite microstructures[END_REF]. An application in the crack propagation case has been implemented in Castem [START_REF] Friderikos | Multiscale gfem with superposition of crack enrichment functions driven by finite fracture mechanics: Theory, first computation and open problems[END_REF]. To summarize multi-scale GFEM (MS-GFEM) enrichment strategy is close to the load based one and illustrated in Fig. 7. Either a stress or strain basis can be used to generate the numerical enrichment functions. the main difference between the GFEM and the MS-GFEM is the emancipation of interactions between feature supports. It can also be noted that the description of the structure only relies on the precalculated Handbook and the enrichment functions as defined in (12) satisfy the SGFEM conditions.

Available tools in commercial codes

The following section presents useful Abaqus tools in order to prove the feasability of the different strategies concerning the GFEM implementation in-

Ω Ω p Ω q ∂Ω u u d ∂Ω F F d 11 22 12
Handbook Fig. 7. Enrichment functions patches chosen from the Handbook positioned in the global domain to take in account the structural details, regarding the MS-GFEM hypotheses. Here a canonic strain base is used to generate the enrichment functions.

troduced previously. For this section the names of Abaqus functionalities will be written in Italic.

Geometrically conforming meshes

A geometrically conforming refined mesh with a macroscopic one allows a same local description in regards to the intersection of two patches. It would result in a simplification of the cross-products within the GFEM stiffness matrix.

Nonetheless only one mesh can be associated with a given CAD structure in Abaqus. Consequently it is not possible to realize a coarse mesh and to super- Then each local node is subject to a displacement boundary condition (in orange).

Enrichment functions seen as perturbation functions

Since the numerical enrichment functions can be seen as a perturbation function between the heterogeneous and the homogeneous structure, the question arises as to how to substract two different displacement fields from two different Abaqus simulations. It is possible to extract the nodal displacements and store them in a Python list for instance. Indeed Python is able to read the Abaqus result file. For instance in the points 1 and 2 of the Fig. 10 the same macroscopic kinematics is imposed at the boundaries of the homogeneous patch and the heterogeneous patch. The difference between the two displacement field can be made with Python and gives the point 3 if reinjected in Abaqus.

Imposing the GFEM approached displacement

The GFEM approached displacement [START_REF] Brennan | Interfacial characterization of a slurry-cast melt-infiltrated SiC/SiC ceramic-matrix composite[END_REF] unknown associated with the enrichment function number j (with j ∈ 1, n (p) enr ) then it is possible to impose the GFEM (or the MS-GFEM) kinematics thanks to the Equation tool available in the Constraint module by writing (here in the case of the GFEM):

u (local) i - N (p) global k=1 ϕ (i) k u k - N (p) global k=1 ϕ (i) k n (p) enr j=1 ψ (j,i) k a (j) k = 0 (16)
Considering a patch Ω p , N enr refer respectively to the number of nodes within the local mesh, the number of nodes within the macroscopic mesh and the number of enrichment functions. ϕ

Substructuring operation with GFEM kinematics

In the GFEM approached displacement [START_REF] Brennan | Interfacial characterization of a slurry-cast melt-infiltrated SiC/SiC ceramic-matrix composite[END_REF], unknowns u k and a k the linear problem to solve for the patch illustrated in Fig. 11 can be expressed as shown in [START_REF] Bobet | Thermal residual stresses in ceramic matrix composites-i. axisymmetrical model and finite element analysis[END_REF], with f l and f b the nodal forces.

  K ll K lb K bl K bb     U l U b   =   f l f b   ( 17
)
As a matter of fact, the set of internal displacements U l is equal to zero because of the links between the macroscopic displacements and the local ones established in [START_REF] Mei | Measurement and calculation of thermal residual stress in fiber reinforced ceramic matrix composites[END_REF]. Consequently the previous system of equations [START_REF] Bobet | Thermal residual stresses in ceramic matrix composites-i. axisymmetrical model and finite element analysis[END_REF] does not

Patch connection applied to the global-local GFEM

Once the Duarte's enriched patches are generated separately they must be repositioned. Then the difficulty lies in the macroscopic nodes superimposition.

Substructure positioning in Abaqus

For the specific application to the MS-GFEM, once the enriched patches have been generated and stored in the Handbook, the question is how to link them to the macroscopic kinematics of the structure. A natural idea is to superimpose the patches over the global structure as shown in Fig. 13 but two issues can be underlined caused by the enrichment procedure:

1. The stiffness matrix of each patch takes into account the total behavior H in the finite element part K F EM in (8)

A substructuring operation forbids an interpolation between the macroscopic kinematics and the u k nodes in each patch

The first difficulty implies a stiffening of the structure because of the superposition of the homogeneous part between the patch and the global structure.

Indeed if in the decomposition of behavior H defined in (4) the homogeneous part corresponds to the silicium MI matrix behavior and ∆H is induced by the structural details then the homogeneous behavior is taken into account both in

MS-GFEM patch

Global structure the finite element part of the GFEM stiffness matrix (15a) and in the homogeneous global structure. Consequently the behavior H 0 must be removed from the finite element part of K GF EM . It is possible by inserting a macroscopic element with the same size as the patch. A negative homogeneous behavior -H 0 is adopted for that element and a master-slave nodes procedure imposes the same nodal displacement between the enriched patch and its finite element support.

Substructure positioning

An interpolation of the global kinematics must be performed on the features placed to model the desired micro-structure. It is technically possible with Abaqus to interpolate the macroscopic FE kinematics from a mesh another with the function Embedded Elements. So the patch support presented previously enables such operation and avoids the second difficulty highlighted. The patch positioning procedure is illustrated in Fig. 14.

Partial conclusion: choice of an enrichment strategy

Strouboulis defines a Handbook problem. Each precalculated pattern is superimposed over a mesh used as a support. Nonetheless the pre-treatment by adaptive procedure and the constraint imposed by the commercial code Abaqus encourages not to adopt that strategy. The Duarte's global-local GFEM approach permits an accurate solution thanks to the respect of the PUM in the GFEM kinematics provided that the local description (and then the mesh size) is fine enough. Moreover the enrichment generation can be parallelized since each patch can be processed independently and patch connections are easy to process. A limitation of this method is that each heterogeneity is taken into account by different enrichment functions due to the local description so the same structural details such as circular fibers with a constant radius are subject to numerous enrichment operations. It also can be noticed that Abaqus requires licenses for each calculus thus a parallelization is potentially greedy1 .

MS-GFEM patch

Similarly to the Strouboulis's Handbook problem, stored precalculated patterns can be used as enrichment functions and called in the global structure in order to build the desired model. The scale separation hypotheses give the advantage of the mesh independency between the features and the macroscopic support.

In terms of Abaqus tools available to write the GFEM kinematics, the enrichment functions can be generated for both the three different strategies introduced previously. Indeed for each presented method a macroscopic load is imposed to the boundaries of the features and allow the numerical calculus of the enrichment function as a perturbation function. The substructuring tool enables local information to be fed back to the scale of the macroscopic patch.

The main advantages of the MS-GFEM compared to the global-local strategy is (i) its ability to generate once for all the different enrichment functions and (ii) the possibility to call one pattern many times thanks to the non-interaction hypothesis. Then the MS-GFEM is potentially a low-cost method next to the Duarte's one.

In any case all the previous steps to build the numerical enriched patches can also be automatized with a Python script. It is then possible to process easily all the micro-structures, their behavior, their mesh parameters, the boundary conditions and the nodal displacement collection to generate the enriched patches.

The different Abaqus tools introduced in the current section prove the feasability of the GFEM implementation in Abaqus, with the Duarte and the multiscale point of view. In the following section a numerical example of the multiscale GFEM is presented. A choice not to present an application by geometrically conforming meshes is volontary made due to the high cost of the reading and writing files related to the generation of enrichment functions. Moreover the tools described in Section 4 are the same for the MS-GFEM. Eventually the enrichment functions generation algorithm is summarized in the Appendix 2 in the case of the MS-GFEM.

Application of the MS-GFEM

Reference problem

A micro-structure of CMC is considered. It is composed of 17 Hi-Nicalon TM -S fibers coated with Boron Nitride (BN). The matrix is heterogeneous and composed of a SiC CVI layer complemented with a silicium MI matrix. A tensile test is performed imposing three surface sliding boundary conditions to prevent rigid body movements and a unitary displacement along direction x as illustrated in Fig. 15. The material properties are summarized in Table 1 and correspond to the HiPerComp R SiC/SiC CMCs [START_REF] Dunn | The effect of fiber volume fraction in HiPerComp R SiC-SiC composites[END_REF][START_REF] Braginsky | Simulation of crack propagation/deflection in ceramic matrix continuous fiber reinforced composites with weak interphase via the extended finite element method[END_REF]. Geometric dimensions are deduced from micrographies in the literature and are given in Table 2. Fibers radius is assumed to be constant and given by the manufacturer [54].

Young Modulus Poisson ratio

Hi-Nicalon S fiber is drawn on Fig. 16 with the Matlab software by reading the Abaqus report file containing the strain components at the Gauss points.

E f ib = 420 GP a ν f ib = 0.2 BN interphase E BN = 62 GP a ν BN = 0.17 SiC CVI matrix E SiC = 400 GP a ν SiC = 0.3 Si MI matrix E Si = 165 GP a ν Si = 0.3

Choice of the features to store in the handbook

In order to recreate the configuration of the reference REV from patterns picked up in a Handbook it is necessary to choose such patterns. Different features can be extracted from the observation of the reference REV: single circular inclusion or two close circular inclusions as illustrated in Fig. 17 and introduced in [START_REF] Tranquart | A computational approach for handling complex composite microstructures[END_REF]. One inclusion is composed of a Hi-Nicalon TM -S fiber, the BN coating and the CVI SiC matrix. BN coating and CVI SiC matrix can interpenetrate in a multi-circular inclusions pattern. The main objective of this section is to find a criterion in order to choose an enrichment support size denoted c thereafter. In the particular case where the load is plane it is possible to draw inspiration from the work of Lekhnitskii [START_REF] Lekhnitskii | Anisotropic plates[END_REF]. For one circular inhomogeneity in an infinite plate like illustrated in Fig. 18 the stress field is given analytically according to potential functions and Dundurs's coefficients [START_REF] Dundurs | Mathematical theory of dislocations[END_REF].

Remembering that the composite inclusions are made of several phases, just one material can be retained to use the Lekhnitskii's model. For a given external stress solicitation the higher the rigidity of the inclusion, the more extensive the deformation field around it. Consequently the chosen radius of the circular inclusion can be taken equal to R inc = R + e BN + e SiC and its Young modulus can be set at E f ib = 420 GP a. If a patch with two fibers must be sized, a circular inclusion embracing the two fibers will be considered. Such choice enables to oversize the patch. Since the patches to be sized must contain as much information as possible about the structure detail that is defined inside, a strain energy criterion is chosen to conclude about the size of each patch. Indeed the knowledge of the analytical strain energy ω(x) leads to the knowledge of its minimum value reached for an infinite distance of the fiber ω ∞ . If x = 0 refers to the center of the circular inclusion and c refers to the length of the patch then then the following criterion is defined:

ω x.x = c 2 = (1 + )ω ∞ ( 18
)
Here x is the base vector from Fig. 18, x is the position within the structure and is an adjustable setting. For a tensile test along direction x component σ xx is predominant and spreads along the same axis x. Fig. 19 shows the evolution of the volumic strain energy along x-axis, for = 10% the chosen relative gap allowed here. For a patch including one fiber the theory leads to a patch size of c = 60 µm which is 10 times the radius of a Hi-Nicalon TM -S fiber. For a patch including two fibers the maximum size calculated is of c = 150 µm which is 25

times the radius of one fiber.

Numerical results

Once the Handbook has been defined, one wishes to reproduce the reference structure defined in Fig. 15. For that purpose the different enriched patterns are positioned on a macroscopic finite element mesh with the procedure described in the last point of Section 4, patterns with one or two fibers are arranged as shown in Fig. 20. Since the patterns have been prefabricated upstream it is not possible to respect the identical arrangement of the fibers on the reference structure, this will lead to a first error when comparing the reference field with the field resulting from the MS-GFEM calculation. However this is not a problem in practice since the arrangement of the fibers is random, it is simply necessary to respect the volume fraction of the different components of the CMC. Moreover by observing Fig. 20 more carefully it is possible to anticipate bad interactions between the fibers. Indeed some features are missing in order to take into account all the interactions between fibers. Eventually some patterns overlap in part which will have the effect of accumulating local enrichments.

The macroscopic structure has a larger volume than the reference REV caused by the necessity of containing entirely the features and their support.

Indeed the Embedded Element function is reliable if the nodes with interpolated displacements are located inside the master structure. To find conditions at the limits similar to those of the reference problem, the homogeneous strain is the same for the new model. The number of degrees of freedom is 708 for this model, it is about a thousand times less than for the finite element model. Then a second idea is to convert the strain fields into a color map so the color of a pixel is associated with a value of a component of the strain tensor. To realize the color map of each patch the Matlab software is used. The strain field is extracted from Abaqus in order to be recreated in Matlab. It turns out that image generation from fields is effective in Abaqus, hence its use. Fig. 21 shows a sample of the first component of the strain tensor calculated for three patches and their equivalent in terms of color map. A resolution of 512 colors is chosen and the color scale is calculated from the strain scale which explains why the colors around the fibers are different for the three patches.

After removing the homogeneous contribution in the strain tensor, the color maps can be summed in order to superimpose the enrichments. Once this step is completed, a rescaling can be made to add the homogeneous part of the strain. The result of this operation is illustrated in Fig. 22. It can be noticed that the range of values taken by the strains is wider than that of the reference calculation. This is explained by the superposition of enrichments when two features are locally overlaid as visible in Fig. 20. Nevertheless there is a correspondence between the location of the strain maxima for the reconstructed model and for the reference model. It is also possible to visualize the strain difference (gap) 11

between the reference strain field denoted

The (gap) 11

is shown in Fig. 23. The greatest visible differences correspond to the areas where the patterns overlap in Fig. 20. This result was already anticipated earlier in the article. Other gaps are calculated at the edges. These differences are caused by (i) strain jumps between the different phases (ii) the impossibility of having identical meshes between the reference structure and the feature-generated structure (iii) the slight possible offset between the fibers of the reference model and the fibers of the reconstructed model. The relative gap is also proposed in Fig. 24, it can be noticed that the maximum error calculated is greater than 400% in some very localized areas. In order not to overwrite the the color scale a thresholded colormap is suggested in Fig. 25. As illustrated in Fig. 26, far from the edges the relative absolute error is less than 10% in the matrix, the interphase or the fiber and the error is higher in the zones of superimposed features. This last example points the importance of the patterns choice to minimize the error made in the reconstruction of the structure. The white zones are defined as "Not a Number" and refer to the errors higher than 100% in absolute value.

Conclusion

This paper has shown the possibility of the numerical implementation of two GFEM strategies within the Abaqus software. A discussion has been conducted in order to select the most versatile method regarding the desire to develop a tool for CMC micro-structures description. Then an application of the MS-GFEM to a micro-structure of a CMC has been made in an elastic problem. It can be noticed that the presented Abaqus tools can also be used to build the Duarte's GFEM.

It has been shown that Abaqus possesses functionalities to (i) generate numerical enrichment functions with regard to the MS-GFEM hypotheses (ii) interpolate a macroscopic kinematics on non-interacting patches which contain structural detail information.

In term of current limitations, two can be noted. First in the MS-GFEM strategy it is important to keep in mind the time allocated to search for the location of patches files and read the associated information during a simulation. Then the time gain due to the limited degrees of freedom is partially compensated with that additional operation. Consequently the MS-GFEM is worthwhile in the case of a rich structure modeling which would demand a very fine mesh if treated by a classic finite element method. Fortunately in the case of a CMC description many structural details must be taken into account and the MS-GFEM becomes a relevant strategy to study such materials, at the condition of an appropriate choice of the patterns and their placement. Secondly the Abaqus substructuring tool gives a representation of a linear and static behavior of a selected volume. Consequently extra researches must be led to take into account non-linear phenomena, a first step in this direction is presented in [START_REF] Friderikos | Multiscale gfem with superposition of crack enrichment functions driven by finite fracture mechanics: Theory, first computation and open problems[END_REF]. Eventually the choice in the patch positioning must be made carefully not to accumulate the perturbation functions describing a common feature. Some perspectives can be highlighted. First investigations concerning the automation of patterns selection could be carried out based on the observation of real microstructures and automatic learning algorithms. The objective would then be to obtain more complex but more representative patterns of the composite while keeping as much as possible the versatile aspect of the GFEM method.

Secondly iterative enrichments could be considered in order to update the state of each patch and to integrate degradation mechanisms such as crack propagation or creep under high temperatures. Numerical methods like POD can be considered as suggested in [START_REF] Aquino | Generalized finite element method using proper orthogonal decomposition[END_REF]. Eventually statistical and random aspects such as the variability of the fibers radius or the respect of the rates (of porosity, fibers, matrix) can also be taken into account in order to distribute the patterns while respecting the representative character of observable micro-structures.

Fig. 1 .

 1 Fig. 1. Comparison between two yarn-scaled structures obtained with a CVI or a CVI-MI process for a SiC-fiber and SiC-matrix (SiC/SiC) composite. Hybrid CVI and MI process results in a fully dense matrix. Micrography from [7].

Fig. 2 .

 2 Fig. 2. Exemple of patches constitutive of a mesh Ω h approaching a real structure Ω, subject to load boundary conditions F d on boundary ∂Ω F and displacement boundary conditions u d on boundary ∂Ωu.

Fig. 3 .

 3 Fig. 3. Crack approximation due to XFEM level-set functions.

  tentially come from explicit meshes of support ω (j) k describing structural details

Fig. 4 .

 4 Fig. 4. Construction of a SGFEM enrichment function. The piecewise linear interpolant of the enrichment function is substracted from the enrichment function ψ

  their behavior. Depending on the description strategy of the local domains ω (j) k

Fig. 5 .

 5 Fig.5. Enrichment generation by load base: patterns are chosen from a homogeneous structure. A canonic load base applied on a refined mesh including structure details is used to calculate numerical enrichment functions. Then such patterns can be superimposed to the global structure to approach the desired model.

  (j) hom (x) calculated for the global homogeneous behavior at the local scale: it is seen as a perturbation contribution of the homogeneous displacement field. By construction such enrichment functions respect the SGFEM. Eventually the created features stored in the Handbook are reinjected in the global model as shown in Fig. 6.

1 .

 1 The support ω(j) k of the enrichment function ψ (j)k is small in comparison with the support Ω k of the hat function ϕ k 2. The enrichment functions do not interact with each other.

k

  can be placed anywhere within the patch. As a convention the additional nodes a (j) k will be located at the centers C (j) k of the enrichment supports ω (j) k .

Fig. 8 .

 8 Fig. 8. Structure partitioning. The partition tool imposes the edges in the global structure so the global and local meshes are forced to fit them. Heterogeneous inclusions are highlighted in red. The characteristic size of the global mesh is about a hundred times larger than that of the local mesh.

4. 2 .

 2 Macroscopic kinematics imposed to the local meshIn the global-local GFEM and the MS-GFEM the numerical enrichment functions can be generated from the local response of the structure to the macroscopic kinematics. Then it is necessary to be able to project the global displacements onto the local mesh. In Abaqus the Submodeling tool drives the local model with the results from the global model. An interpolation is then performed to determine the solution of the submodel. For instance the point 1 of the Fig.9illustrates a hat function of a macroscopic patch. Nodes at the boundary of the patch are fixed and the central one is associated with a unitary displacement. Point 2 shows the displacements calculated for the submodel of the patch if the macroscopic kinematics is imposed on the whole surface of

Fig. 9 .

 9 Fig. 9. Projection of a global field on a local mesh through a submodeling strategy. The local mesh is assumed to be a submodel of the macroscopic structure and follows its kinematics.

Fig. 10 .

 10 Fig. 10. Numerical enrichment function seen as the perturbation between the heterogeneous response and the homogeneous response to identical boundary conditions.

k

  indicates the value of the macroscopic hat function interpolated at the coordinate x i and ψ (j,i) k denotes the value of the enrichment function at the same point x i . Nodal displacements u k and enrichment unknown nodes a (j) k must be added due to the definition of new Reference Points in order to be reused for other equation definitions. Indeed Abaqus eliminates the terms used to describe the local displacements u (local) i for a memory gain purpose. Assimilating the macroscopic nodes with Reference Points forbids Abaqus to erase them.

k

  are defined on macroscopic nodes. Nonetheless the structure behavior is contained in the local description and more particularly in the enrichment functions ψ (j) k . A static condensation will be necessary to integrate the stiffness matrix taking into account the GFEM kinematics. For that purpose the Substructuring tool enables an elimination of the internal nodes of a structure and keeps the Retained Nodal Degrees of Freedom while ensuring a linear response of the created substructure. The substructuring strategy can work in conjunction with the Equation Constraint introduced previously. The local nodes u (local) i described by the GFEM kinematics are assumed to be the internal nodes of the mesh to be substructured and the macroscopic nodes u k and a (j) k are supposed to be the retained nodal degrees of freedom. The Reference Points which are necessary in the definition of the Constraint Equations can also be associated with retained nodal degrees of freedom. If U l is the vector containing all the local displacements u (local) i and U b is the vector containing all the macroscopic displacements u k and the enrichment unknowns a (j)

  need to be solved to find the U b set. The Substructuring tool used in accordance with the Constraints (16) simply allows the construction of a macroscopic object with the desired GFEM stiffness matrix.

Fig. 11 .

 11 Fig. 11. Substructuring operation. A static condensation is performed in order to express the stiffness matrix on the macroscopic nodes of a patch. The local nodes are represented by the yellow square, the macroscopic nodes by the red circles.

Fig. 12 shows

 12 Fig. 12 shows two superimposed patches composed of four macroscopic elements each. The enrichment nodes are represented by the square points for one patch and by triangles for the other patch. Macroscopic displacement nodes are represented by the circle points. The common macroscopic element is in red and two circular heterogeneities are drawn. Superimposed nodes schematized by the same geometric shape are assumed to have the same kinematics. The Tie tool enables to impose the same displacement between two nodes with the same coordinates.

Fig. 12 .

 12 Fig.12. Two superimposed patches composed of four macroscopic element each. Superimposed nodes with the same geometric shape (circle, square or triangle) are assumed to be equal.

Fig. 13 .

 13 Fig. 13. Superimposition of an enriched patch over a homogeneous global structure.

Fig. 14 .

 14 Fig. 14. Superimposition of an enriched patch over an intermediate homogeneous patch with an opposite behavior to compensate the homgeneous one from the global structure.

Fig. 15 .

 15 Fig.[START_REF] Bobet | Study of thermal residual stresses in ceramic matrix composites[END_REF]. Reference Element Volume considered. The fibers are colored in red, the BN coating in yellow, the CVI SiC in blue and the Si MI in grey. A hypothesis of plane strain is adopted. A traction load along x is applied with u d = 1.x.

  Hexahedral elements are adopted in the reference simulation, the thickness REV Depth 5 µm REV Side length c = 60 µm Hi-Nicalon TM -S radius R = 6 µm BN interphase thickness e BN = 0.4 µm [55]SiC CVI thickness e SiC = 2 µm[START_REF] Katoh | Continuous SiC fiber, CVI SiC matrix composites for nuclear applications: Properties and irradiation effects[END_REF]

Fig. 16 .

 16 Fig. 16. First component of the strain tensor denoted E 11 written at the nodes with the Matlab software.

Fig. 17 .(E 1 , ν 1 )(E 2 , ν 2 )Fig. 18 .

 17112218 Fig.17. Examples of patterns extracted from the reference REV. Several features with two fibers must be considered because of the possible different distances between them. The fiber is in red, the interphase is yellow, the CVI matrix is blue and the MI matrix is grey. Each feature is a square of length c which must be dimensioned.

Fig. 19 .

 19 Fig. 19. Variations of the volumic strain energy along x-axis for a tensile test of amplitude σ ∞ xx . For x.x = c 2 the calculated energy reaches (1 + )ω∞.

Fig. 20 .

 20 Fig. 20. Patterns arrangement in order to recreate the reference REV. Features with one or two fibers are used, they may also partially overlap each other.

Fig. 21 .

 21 Fig. 21. Example of the conversion of the Abaqus strain field into a color map with a resolution of 512 colors. Each color scale is adapted from the strain scale.

Fig. 22 .Fig. 23 .Fig. 24 .Fig. 25 .

 22232425 Fig. 22. Strain field obtained after the post-treatment of all the contributions of the patches, for the second positioning. The color bar scale given corresponds to the strain recalculated from the color maps of the different patterns.

Fig. 26 .

 26 Fig. 26. Values of interest picked in the relative gap map (values in %).

Table 1 .

 1 material properties of interest in this study

Table 2 .

 2 Dimensions used in this model of the structure is swept by one element. The FE model contains 524 588 degrees of freedom and only one CPU is used. The first component of the strain tensor

If N P roc is the number of processors used for the computation, the number of Abaqus licenses can be deduced from the formula: 5.(Nproc)0.442

Acknowledgements

The authors wish to acknowledge the support of Safran Ceramics and the ANRT for funding this research.

Appendix 1

If u h = u F EM + u enr denotes the GFEM kinematics as the sum of the classical finite element displacement u F EM and the enriched contribution u enr and u * h = u * F EM + u * enr denotes the associated virtual field the virtual work of the internal forces can be developed as:

The Voigt convention gives (u enr + l. As a convention if the index p (respectively q) in the expression of α (respectively β) equals 1 it is assumed that α = k (respectively β = l). Since the hat functions ϕ k are not equal to zero only on the patch Ω k the domains of integration are reduced to the ones precised in [START_REF] Evans | The physics and mechanics of fibre-reinforced brittle matrix composites[END_REF].

Macroscopic shape functions ϕ

projected on the local mesh, on each local node i the value is denoted ϕ

Canonical strain basis applied on a refined mesh: