K Du

L Cheng

J F Barthélémy

I Sevostianov

A Giraud

A Adessina

Effective elastic properties of transversely isotropic materials with concave pores

Keywords: concave shape, supersphere, superspheroid, compliance contribution tensor, effective elastic properties, transversely isotropic matrix

The aim of this paper is to extend recent works devoted to the study of the effect of 3D pores of concave shape embedded in isotropic matrix to the case of transversely isotropic (TI) matrix. In the first part of the paper, approximate relations for the compliance contribution tensor of pores of two reference shapes, supersphere and axisymmetrical superspheroid, are developed on the basis of 3D Finite Element Modelling, recently presented in a companion paper, and known exact solutions for the limiting cases of spherical pores and circular crack. In the second part application to effective elastic coefficients of transversely isotropic materials such as clay rocks, in the frame of homogenization theory is presented to illustrate the impact of concavity parameter on overall properties.

Introduction

In the present paper, we analyze the effect of the concavity of pores on the overall elastic properties of a porous material with transversely-isotropic solid phase. For this goal, we use two homogenization techniques: Mori-Tanaka-Benveniste scheme and Maxwell scheme. Both of them are based on the solution for a single inhomogeneity problem and can be easily formulated for ellipsoidal inhomogeneities using Eshelby results [START_REF] Eshelby | Elastic inclusions and inhomogeneities[END_REF]. Non-ellipsoidal shapes of the inhomogeneities are not so well studied and most of the results are obtained in 2 -D by conformal mapping [START_REF] Kachanov | Effective Moduli of Solids With Cavities of Various Shapes[END_REF]. For three-dimensional case, the problem of irregular (non-ellipsoidal) inhomogeneities reduces to integral equations and generally requires computational approaches although, in some cases, solution can be obtained in the form of infinite series, see, for example [START_REF] Krasnitskii | Effect of a rigid toroidal inhomogeneity on the elastic properties of a composite[END_REF] . They can be generally subdivided onto two groups: (i) direct computation of stress and strain fields for a given (deterministic) microstructure by discretizing the domain and using the FEM, and then post-processing the averages of the stress and strain fields (see, for example, [START_REF] Garboczi | Elastic moduli of composites containing a low concentration of complex-shaped particles having a general property contrast with the matrix[END_REF]) and (ii) computation of the contribution of one isolated inhomogeneity into the effective elastic properties as a function of its shape. The latter results constitute basic building blocks for theoretical models that cover diverse orientation distributions and concentrations of inhomogeneities. [START_REF] Böhm | Effects of particle shape on the thermoelastoplastic behavior of particle reinforced composites[END_REF] analyzed shape effects on the effective elastic and thermal properties of the composites containing randomly oriented and distributed spherical, octahedral, cubical and tetrahedral particles. [START_REF] Drach | Prediction of the effective elastic moduli of materials with irregularly-shaped pores based on the pore projected areas[END_REF] proposed to evaluate effect of pores of irregular shape on the overall elastic moduli using pore projected areas. This approach works well for prediction of the overall Young's moduli in different directions. [START_REF] Drach | Comparison of full field and single pore approaches to homogenization of linearly elastic materials with pores of regular and irregular shapes[END_REF] performed comprehensive numerical analysis of the pore shape on the overall properties of solids with porosity levels up to 25%. [START_REF] Trofimov | Effective elastic properties of composites with particles of polyhedral shapes[END_REF] compared predictions of overall elastic properties of composites reinforced with particles of a different polyhedral shapes by FEM and micromechanical schemes. The results of the two approaches are in good agreement for volume fractions up to 30% for all studied material combinations. The inverse problem -design of material microstructure has been has been done by Zohdi [START_REF] Zohdi | Genetic design of solids possessing a random-particulate microstructure[END_REF], who determine optimal geometrical and mechanical properties of inhomogeneities for prescribed overall elastic moduli. Effect of the concavity factor of superspheres and axisymmetric concave pores was analyzed in the works of [START_REF] Chen | Evaluation of the effective elastic and conductive properties of a material containing concave pores[END_REF][START_REF] Chen | Combined effect of pores concavity and aspect ratio on the elastic properties of a porous material[END_REF][START_REF] Sevostianov | Compliance and resistivity contribution tensors of axisymmetric concave pores[END_REF]. The authors supplemented finite element modelling with analytical approximations for compliance contribution tensors of pores of such shapes. These results were used to calculate overall elastic properties of materials with multiple concave pores: oolitic rock [START_REF] Kalo | Effective elastic moduli of a heterogeneous oolitic rock containing 3-d irregularly shaped pores[END_REF] and 3 -D printed S i 3 N 4 ceramics [START_REF] Lurie | Mechanical behavior of porous S i 3 N 4 ceramics manufactured with 3-D printing technology[END_REF]. All the mentioned results have been obtained for materials with isotropic matrix. The number of explicit results on elastic properties of heterogeneous materials with anisotropic matrix is substantially smaller. Piezoelectric properties of transversely isotropic materials containing circular fibers aligned with the axes of symmetry of the matrix have been calculated using various homogenization techniques by [START_REF] Sevostianov | On the modeling and design of piezocomposites with prescribed properties[END_REF]. [START_REF] Sevostianov | Effective elastic properties of the particulate composite with transversely isotropic phases[END_REF] calculated compliance contribution tensor for a spheroidal inhomogeneity of arbitrary aspect ratio embedded in a transversely-isotropic material. [START_REF] Levin | Elastic properties of inhomogeneous transversely isotropic rocks[END_REF] calculated effective elastodynamic properties of transversely isotropic rocks containing aligned spherical and strongly oblate spheroidal pores. Effective porothermoelastic properties of transversely isotropic rocks such as mudstones, argillites, shales have been studied in the frame of Effective Media Theory (EMT), by [START_REF] Giraud | Effective poroelastic properties of transversely isotropic rocks-like composites with arbitrarily oriented ellipsoidal inclusions[END_REF]. [START_REF] Cosenza | Effect of the local clay distribution on the effective elastic properties of shales[END_REF] calculated overall properties of transversely-isotropic clay containing spherical inhomogeneities. [START_REF] Vasylevskyi | On micromechanical modeling of orthotropic solids with parallel cracks[END_REF] calculated overall properties of a transversely-isotropic material containing parallel circular cracks. [START_REF] Seyedkavoosi | Randomly oriented cracks in a transversely isotropic material[END_REF] used approach developed by [START_REF] Guerrero | On a possible approximation of changes in elastic properties of a transversely isotropic material due to an arbitrarily oriented crack[END_REF] to calculate overall properties of a transversely-isotropic material containing arbitrarily oriented cracks. This result was used by [START_REF] Seyedkavoosi | Modeling of the overall elastic behavior of a transversely isotropic material reinforced with arbitrarily oriented transversely isotropic platelets[END_REF] to estimate properties of a transverselyisotropic material with multiple arbitrarily oriented oblate inhomogeneities and applied to calculation of the overall properties of dentin.

In the text to follow, we use the recent numerical results of [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF] where compliance contribution tensors of concave pores in a transversely-isotropic material are obtained. We suggest an analytical approximation of the components of these tensor using approach of [START_REF] Trofimov | Overall elastic properties of a material containing inhomogeneities of concave shape[END_REF] and evaluate overall elastic properties of transversely-isotropic materials containing such pores. The results are illustrated by example of shale rock containing concave pores.

Compliance contribution tensor of a concave pore

Refer to appendices 7 and 8 for background on tensors and property contribution tensors. Recently [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF] calculated components of the compliance contribution tensors of superspheroidal and axisymmetric superspheroidal pores embedded in a transversely-isotropic material numerically. • axisymmetrical superspheroidal pore In what follows we will only consider ς = 1, the first shape is then a supersphere, and the second shape obtained by a rotation about symmetry axis x 3 . Supersphere and axisymmetrical superspheroid coincide with sphere in the case p = 1 but strongly differ in the limiting case p → 0: supersphere tends to three orthogonal needles along coordinates axes and superspheroid tends to a circular crack of unit radius crossed by a perpendicular needle along symmetry axis x 3 .

• superspheroidal pore | x 1 a | 2p + | x 2 a | 2p + | x 3 ςa | 2p = 1 (1)
x 2 1 + x 2 2 a 2 p + | x 3 ςa | 2p = 1 (2)
Compliance contribution tensor H E 0 of a superspherical pore aligned with the direction of a TI matrix, with symmetry axis x 3 , respects tetragonal symmetry. With three orthogonal planes of symmetry (with normal e i ) and equivalence between x 1 and x 2 axes, its [START_REF] Trofimov | Effective elastic properties of composites with particles of polyhedral shapes[END_REF] as

H E 0 = 6 i=1 h i E i (3)
Average compliance contribution tensors H E 0 of superspherical and axisymmetrical superspheroidal pores have been numerically calculated by using 3D Finite Element Method (FEM), see [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF] considering TI matrix with elastic properties of a shale (see Tables 12). Numerical results of [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF] are summarized in Tables (7-8) (see appendix 9). In the particular case of an ellipsoidal pore E embedded in an infinite matrix 0 of stiffness C 0 and compliance S 0 tensors, compliance H E 0 and stiffness N E 0 contribution tensors are analytical and write (see [START_REF] Kachanov | Micromechanics of materials, with applications[END_REF] for details) :

H E 0 = (S E -S 0 ) -1 + Q E 0 -1 , N E 0 = (C E -C 0 ) -1 + P E 0 -1 (4)
where P E 0 and Q E 0 denote the fourth order Hill's tensors [START_REF] Hill | A self-consistent mechanics of composite materials[END_REF] of the inhomogeneity, related by relation

Q E 0 = C 0 : I -P E 0 : C 0 ( 5
)
Strain Hill tensor P E 0 of a spheroidal inhomogeneity aligned in a TI matrix may be found in [START_REF] Sevostianov | Effective elastic properties of the particulate composite with transversely isotropic phases[END_REF][START_REF] Barthélémy | Simplified approach to the derivation of the relationship between Hill polarization tensors of transformed problems and applications[END_REF] and it is recalled for convenience in appendix 8. In the next Section we approximate these results analytically and then use them to calculate overall elastic properties of transversely isotropic matrix containing multiple concave pores.

3. Approximation formula for compliance contribution tensor of a superspherical or axisymmetrical superspheroidal pore embedded in a transversely isotropic host matrix

We investigate in this section the extension to transverse isotropy of approximation formula for the compliance contribution tensor of 3D pores of particular shapes previously presented. We restrict the study to the following assumptions

• same directions of symmetry between matrix and pore inclusion (aligned case)

• study is focused on the concavity parameter p • in the case of superspherical pore, compliance contribution tensor respects tetragonal symmetry (6 independent components) but we will consider for applications random orientation distributions in the isotropic plane x 1x 2 (x 3 denotes the symmetry axis) on the one hand, in 3D space on the other hand. Related compliance contribution tensors respects transversely isotropic symmetry.

Volume and surface area of superspherical and axisymmetrical superspheroidal pores

Approximation formula may be obtained by using basic geometric information related to the considered reference shapes, supersphere and axisymmetrical superspheroid, defined in relations (1-2), with aspect ratio ς = 1. These informations are volume, total surface area and projected areas onto planes 0x i x 3 (with i = 1, 2).. and corresponding volumes write (Γ denotes Euler Gamma function, see [START_REF] Chen | Evaluation of the effective elastic and conductive properties of a material containing concave pores[END_REF][START_REF] Sevostianov | Compliance and resistivity contribution tensors of axisymmetric concave pores[END_REF][START_REF] Trofimov | Overall elastic properties of a material containing inhomogeneities of concave shape[END_REF] for details)

V se (p) = 2 3 Γ 1 2p 3 p 2 Γ 3 2p , V so (p) = 4 π 3 Γ 1+2p 2p Γ 1 p Γ 3 2p ( 6
)
where superscripts se and so respectively refer to supersphere and axisymmetric superspheroid. In the range 0 < p ≤ 1, supersphere and axisymmetrical superspheroid with unit semi-length a are superscribed by unit sphere of volume V 0 = 4π/3. Ratios V se (p)/V 0 , V so (p)/V 0 and V se (p)/V so (p) are presented in figure [START_REF] Krasnitskii | Effect of a rigid toroidal inhomogeneity on the elastic properties of a composite[END_REF].. Supersphere and axisymmetrical superspheroid coincide with sphere in the case p = 1 but strongly differ in the limiting case p → 0: supersphere tends to three orthogonal needles along coordinates axes and superspheroid tends to a circular crack of unit radius crossed by a perpendicular needle along symmetry axis x 3 . . Except for some particular values of concavity parameter (p = 1 4 , 1 2 , 1), the total surface area needs to be calculated by numerical integration. As in [START_REF] Trofimov | Overall elastic properties of a material containing inhomogeneities of concave shape[END_REF], we use the surface area of the supersphere A se (p) given by [START_REF] Trott | The Mathematica Guidebook for Numerics[END_REF], and the surface area A so (p) of the axisymmetrical superspheroid is given by the single integral accounting for symmetry of revolution

A so (p) = 1 0 1 -x 2p 1 2p 1 + x -2(1+2p) 1 -x 2p 1-2p p 1 2 dx (7)
The projection area S proj (p) of both 3D shapes onto planes x i x 3 (i = 1, 2, x 3 denotes symmetry axis of the axisymmetrical superspheroid) writes (Beta denotes Euler Beta function) /V so V se (p)/V so (p)

S proj (p) = 2 p Beta 1 + 1 2p , 1 2p (8
Figure 3: Left: ratios V se (p)/V 0 and V so (p)/V 0 functions of concavity parameter p, right: ratio V se (p)/V so (p)

Approximation formula for superspherical pore

We restrict this study to the range 0.2 ≤ p ≤ 1, on the basis of [START_REF] Trofimov | Overall elastic properties of a material containing inhomogeneities of concave shape[END_REF] we propose approximation formula

H i jkl (p) =      S proj (p)/(V se (p)) 2/3 S octa proj /(V se octa ) 2/3 f se-a i jkl (p) H octa i jkl 0.2 ≤ p < 0.5 , no sum over i and j A se (p)/(V se (p)) 2/3 A octa se /(V se octa ) 2/3 f se-b i jkl (p) H octa i jkl 0.5 ≤ p ≤ 1 , no sum over i and j (9)
where octa denotes octahedron (particular case of supersphere at p = 0.5 represents an octahedron). Functions f se-a i jkl (p) and f se-b i jkl (p) are given in appendix. A quadratic fit has been considered for components H ii j j (no sum over i and j) whereas a fourth degree polynomial has been necessary to fit shear components H 1212 , H 2323 . Comparisons between approximate relations [START_REF] Zohdi | Genetic design of solids possessing a random-particulate microstructure[END_REF]

FEM approx

Figure 4: The 6 independent components H i jkl of the tetragonal compliance contribution tensor of a superspherical pore embedded in TI matrix, as a function of concavity parameter p. Comparison between FEM results (dashed lines) and approximate relations (plain lines).

Approximation formula for axisymmetrical superspheroidal pore

Similar approximation formula are proposed for components H 1111 , H 1122 , H 1133 but total surface area has been used instead of projection area, in the concave range p < 0.5.

H 11ii (p) =      A so (p)/(V so (p)) 2/3 A octa so /(V so octa ) 2/3 f so-a 11ii (p) H octa 11ii 0.2 ≤ p < 0.5 , no sum over i, i ∈ [1, 2, 3] A so (p)/(V so (p)) 2/3 A octa so /(V so octa ) 2/3 f so-b 11ii (p) H octa 11ii 0.5 ≤ p ≤ 1 , no sum over i, i ∈ [1, 2, 3] (10)
where octa corresponds to the case (p = 0.5) which is not an octahedron but a double-conical shape. Functions f so-a 11ii (p) and f so-b 11ii (p) are given in appendix. Polynomials of degree 4 have been considered for both concave and convex domains in the range 0.2 < p < 1. It may be noticed that semi-analytical approximations using the limiting cases of aligned circular crack p → 0 and sphere p → 1 may be used for components H 3333 and H 2323 (with x 3 symmetry axis of axisymmetrical superspheroid and TI matrix ). Approximate solutions writes

H 3333 (p) = V sphere V so (p) 1 -p 1 -α o H c 3333 + p -α o 1 -α o H sphere 3333 , α o = 0.19, 0.2 ≤ p ≤ 1 ( 11
)
H 2323 (p) =    f 2323 (p) 0.2 ≤ p < 0.5 H sphere 2323 0.5 ≤ p ≤ 1 ( 12
)
f 2323 (p) = A + B exp (-ω(p -p 0 )) , p 0 = 0.2 (13)
Constants A and B are determined by using analytical solutions for the limiting cases of circular crack and sphere, which are imposed at p = 0.2 and p = 0.5 (difference between solutions of axisymmetrical octahedron p = 0.5 and sphere are not significant for H 2323 component).

f 2323 (p = 0.2) = H c

2323

V sphere V so (0.2) , f 2323 (p = 0.5) = H sphere 2323 [START_REF] Lurie | Mechanical behavior of porous S i 3 N 4 ceramics manufactured with 3-D printing technology[END_REF] and then

A + B = H c 2323 V sphere V so (0.2) , A + B exp (-0.3 ω) = H sphere 2323 ( 15
)
Constant ω is determined by fit of finite element results (see appendix 9 and [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF] for details on finite element modellings). One obtains

A = 0.221795, B = 0.591786, ω = 22 (16)
Analytical solution for the compliance contribution tensor of a spheroidal pore aligned with the directions of a TI matrix is recalled in appendix. It may be noticed that approximation [START_REF] Chen | Combined effect of pores concavity and aspect ratio on the elastic properties of a porous material[END_REF] for H 3333 component numerically coincides with the corresponding component of a spheroidal pore with same volume than superspheroidal pore :

H 3333 (p) ≈ H spheroid 3333 (γ(p)), 0.2 ≤ p ≤ 1 ( 17
)
with

γ(p) = V so (p) V sphere = Γ 1+2p 2p Γ 1 p Γ 3 2p ( 18
)
The analytical solution for the aligned spheroidal pore

H spheroid i jkl (it

Evaluation of the effective elastic properties of materials with transversely isotropic matrices

In this section, we calculate effective elastic properties using three homogenization techniques: Non Interaction Approximation, Mori Tanaka-Benveniste and Maxwell schemes (respectively referred with superscripts NIA, MTB and MX), see [START_REF] Mori | Average stress in matrix and average elastic energy of materials with misfitting inclusions[END_REF][START_REF] Benveniste | A new approach to the application of Mori-Tanaka's theory in composite materials[END_REF][START_REF] Kachanov | Micromechanics of materials, with applications[END_REF])

S NIA = S 0 + ϕH E 0 , S MTB = S 0 + ϕ 1 -ϕ H E 0 , S MX = S 0 + 1 ϕ H E 0 -1 -Q Ω 0 -1 (19)
where ϕ denotes the porosity. Q Ω 0 denotes the second Hill tensor of the effective inclusion of the Maxwell scheme, which is supposed of spheroidal shape (with aspect ratio γ Ω ) and aligned with the directions of the TI host matrix. Q Ω 0 is related to the strain Hill tensor P Ω 0 by the relation (see appendix 8 for details):

Q Ω 0 = C 0 : I -P Ω 0 : C 0 ( 20
)
For numerical examples, we use elastic constants of shale and mudstone given in Tables 1 and2 ( [START_REF] Giraud | Effective porothermoelastic properties of transversely isotropic rock-like composites[END_REF][START_REF] Cosenza | Effect of the local clay distribution on the effective elastic properties of shales[END_REF]. In what follows, effective elastic properties of porous clay matrix at mesoscopic scale are estimated by homogenizing micropores. We do not consider solid mineral inclusions of calcite and quartz which would need to be added for the transition from mesoscopic to the macroscopic scale (the centimeter scale of standard geomechanical laboratory tests).

The porosity of clay matrix ϕ (denoted f I p in [START_REF] Giraud | Effective porothermoelastic properties of transversely isotropic rock-like composites[END_REF]) is comprised in the range ϕ ≤ 0.30 which has been considered for the sensitivity study.. It must be emphasized that most of the existing homogenization results do not account for anisotropy of the host matrix attributing the overall anisotropy to the microstructure of the pore space. This assumption is invalid for shale rock, in particular. For geomaterials, the superspherical shape of pores is more realistic than the axisymmetric one since it approximately represents intergranular pores. We consider only a random orientation distribution of superspherical pores which does not violate the orientation of the symmetry axes of the transversely isotropic matrix. Transverse isotropic projection of compliance contribution tensor Π T I H E 0 will be used instead of the compliance contribution of the superspherical pore relative distance between H tensor and its TI projection as a function of concavity parameter. This relative distance is equal to zero in the case of the spherical inclusion p = 1.

H E 0 .

Aligned axisymmetric superspheroidal pores

Effective elastic coefficients

E MTB 1 , E MTB 3 , G MTB 31
obtained with Mori-Tanaka-Benveniste (MTB) approximation are presented in figure [START_REF] Trofimov | Effective elastic properties of composites with particles of polyhedral shapes[END_REF]. Approximation formula (10-11-12) deduced from FEM are compared to approximation of compliance contribution tensor of an oblate spheroidal pore of same volume (semi axis length of axisymmetrical superspheroid is equal to the greater semi axis length of the oblate spheroid) ):

H i jkl (p) ≈ H spheroid i jkl (γ(p)), 0.2 ≤ p ≤ 1 (23)
with aspect ratio of oblate spheroid γ(p) defined by relation ( 18 It may be observed that:

• normal Young's modulus E MTB 3 may be estimated by a very simple approximation using compliance contribution of an oblate spheroidal pore of same volume than the axisymmetrical superspheroidal pore, in the concavity range p < 0.5. Numerically it is mainly related to normal component of compliance contribution tensor H 3333 . Comparison between axisymmetrical superspheroidal and oblate spheroidal pores shows that both approximations lead to the same effective coefficient E MTB 3 . In other words, concavity or convexity has no significant effect on this coefficient, when comparing same pores of same volume. It must be emphasized that this result is very specific and cannot be generalised. It only holds for this particular shape and the normal Young's modulus..

• oppositely, transverse Young's modulus E MTB 1 and shear coefficient G MTB 31 are strongly related to the concav-ity parameter p. Comparison of estimates based on concave and convex pores of same volume (respectively axisymmetrical superspheroid and oblate spheroid), in the range p < 0.5, see figure [START_REF] Trofimov | Effective elastic properties of composites with particles of polyhedral shapes[END_REF], shows significant differences. This result is expected as an approximation based on an oblate spheroid of same volume is not precise for all components H 1111 , H 1122 , H 1133 , H 2323 , particularly in the concave range 0.2 < p < 0.5. It confirms, for a 3D shape embedded in an anisotropic matrix that the concavity parameter is of major importance when estimating effective elastic properties. 4.2. Aligned axisymmetric superspheroidal pores compared with a random orientation distribution of superspherical pores in the transverse plane Comparisons of effective elastic moduli obtained with (MTB) approximation for aligned axisymmetrical superspheroidal pores and randomly oriented superspherical pores in the transverse plane are presented in figure [START_REF] Zohdi | Genetic design of solids possessing a random-particulate microstructure[END_REF]. Obtained effective porous material is transversely isotropic with same symmetry axis than matrix. It may be observed that • effects of these two shapes on elastic effective properties are strongly different in the concave range 0.2 < p < 0.5. It is expected as the supersphere tends to three orthogonal needles (with zero volume and zero surface) when p tends to zero, whereas the axisymmetrical superspheroid tends to a circular crack with one central orthogonal needle (the latter having zero volume but non zero surface).. As previously indicated, the most relevant shape compared to microstructure of porous materials is certainly supersphere.

• a significant anisotropic degree in the case of aligned axisymmetrical superspheroidal pores, in the limit p → 0. as a fonction of porosity ϕ, MTB approximation, superspherical (plain lines) and aligned axisymmetrical superspheroidal (dashed lines) pores.

Comparisons between NIA, MTB and Maxwell homogenization schemes

Effective elastic properties predicted by Maxwell, MTB, NIA are presented in figures [START_REF] Chen | Evaluation of the effective elastic and conductive properties of a material containing concave pores[END_REF][START_REF] Chen | Combined effect of pores concavity and aspect ratio on the elastic properties of a porous material[END_REF] for respectively for axisymmetrical superspheroidal and superspherical pores randomly oriented in transverse plane. The shape of the effective inclusion of the Maxwell scheme is still an open issue when host matrix is anisotropic (see [START_REF] Sevostianov | On the shape of effective inclusion in the Maxwell homogenization scheme for anisotropic elastic composites[END_REF][START_REF] Giraud | Effective electrical conductivity of transversely isotropic rocks with arbitrarily oriented ellipsoidal inclusions[END_REF]). The sensitivity study on the shape of the effective inclusion, and oblate spheroid of aspect ratio γ Ω = 0.5 -1 confirms that it is a parameter of major importance.

Concluding remarks

In the present work, effective properties of a transversely-isotropic material containing concave pores are discussed and illustrated on the example of porous clay matrix. For this goal we used non-interaction approximation, Mori-Tanaka-Benveniste and Maxwell homogenization schemes. All techniques require the explicit analytical representation of the compliance contribution tensor for a single pore. These tensors were calculated for the set of superspherical and axisymmetrical superspheroidal pores with concavity parameter p in the range 0.2 ≤ p ≤ 1 which covers both concave (0.2 ≤ p < 0.5) and convex shapes (0.5 ≤ p ≤ 1) using FEM. Based on the numerical solution for octahedron (p = 0.5) and analytical solution for sphere we built analytical approximations of the compliance contribution tensor in terms of the pore concavity parameter. The accuracy of this approximation is better than 5% for all the tensor components. We show that the concavity parameter p is a parameter of major importance on the overall elastic behavior of transversely isotropic materials containing such pore shapes. It is impossible to match effect of concave pores by oblate spheroidal pores except in the specific case of normal Young's modulus and axisymmetrical shape. The main novelty of this study is the account of concavity effects related to 3D shapes embedded in an anisotropic matrix while previous studies where done in the case of isotropic matrix. On the basis of presented results, it appears that it is not possible to separate effect of anisotropy from the effect of concavity to extend previous results obtained in the isotropic case. Accordingly, it is not possible to take effect of concavity from isotropic matrix and put it into TI matrix.

Acknowledgements

JFB thanks Cerema and Université Gustave Eiffel for allowing the creation of the joint research team entitled " Équipe de Recherche Commune sur les Matériaux pour une Construction Durable (ERC MCD)", within which these research works were carried out. IS acknowledges financial support from National Aeronautics and Space Administration (NASA) Cooperative Agreement NNX15AL51H.

Appendix. Background on tensors

Notations : Barred letters A, C, D, Q refer to fourth order tensors, bold letters ε, σ, i refer to second order tensors, underlined letters x, x refer to first order tensors. Einstein's summation convention over repeated indices is used unless otherwise indicated. ⊗, : and :: respectively represent tensor product, (dot product), contracted products on two and four indices. i, I, J and K = I -J respectively represent the second-rank identity tensor, the fourth-rank symmetric identity tensor, and fourth-rank spherical and deviatoric isotropic projectors (δ i j denotes Kronecker delta symbol,

δ i j = 1 if i = j, δ i j = 0 otherwise). a ⊗ b = a i b j e i ⊗ e j , a s ⊗ b = 1 2 a i b j + a j b i e i ⊗ e j , a ⊗ b = a i j b kl e i ⊗ e j ⊗ e k ⊗ e l (24)
a : b = a i j b ji , A : B = A i jop B pokl e i ⊗ e j ⊗ e k ⊗ e l , A ::

B = A i jkl B lk ji (25) a⊗b = 1 2 a ik b jl + a il b jk e i ⊗ e j ⊗ e k ⊗ e l (26)
J = 1 3 i ⊗ i, I = i⊗i, i = δ i j e i ⊗ e j , J i jkl = 1 3 δ i j δ kl , I i jkl = 1 2 δ ik δ jl + δ il δ jk (27)
It may be interesting to introduce standard notation and the corresponding simplified algebra for fourth-order transversely isotropic tensor (see [START_REF] Walpole | Fourth-rank tensors in the thirty-two crystal classes: multiplication tables[END_REF]). See also post of Sébastien Brisard on github, http://sbrisard.github.io/, intitled Decomposition of transverse isotropic, fourth-rank tensors. By denoting n the unit vector of symmetry axis of the material, let us introduce the second-order tensors

i N = n ⊗ n = n i n j e i ⊗ e j , i T = i -i N (28)
In the particular case of n = e 3 , (28) writes

i N = e 3 ⊗ e 3 , i T = e 1 ⊗ e 1 + e 2 ⊗ e 2 (29)
One introduces fourth-order tensors

E 1 = i N ⊗ i N , E 2 = 1 2 i T ⊗ i T , E 3 = 1 √ 2 i N ⊗ i T , E 4 = 1 √ 2 i T ⊗ i N (30)
E 5 = i T ⊗i T - 1 2 i T ⊗ i T , E 6 = i T ⊗i N + i N ⊗i T (31)
It may be shown that any transversely isotropic fourth-order tensor can be decomposed as

L = 6 i=1 l i E i = l i E i (32)
Considering symmetry axis equal to n = e 3 , Walpole matrix representation of tensor H = h i E i writes

H =                h 2 +h 5 2 h 2 -h 5 2 h 4 √ 2 0 0 0 h 2 -h 5 2 h 2 +h 5 2 h 4 √ 2 0 0 0 h 3 √ 2 h 3 √ 2 h 1 0 0 0 0 0 0 h 6 0 0 0 0 0 0 h 6 0 0 0 0 0 0 h 5                (33)
Relations with usual H i jkl components write [START_REF] Giraud | Effective electrical conductivity of transversely isotropic rocks with arbitrarily oriented ellipsoidal inclusions[END_REF] We detail herefafter calculation of compliance contribution H E 0 of a spheroidal pore E in TI basis, in terms of components c i and p i of matrix stiffness tensor C 0 and strain Hill polarisation tensor P E 0

h 1 = H 3333 , h 2 = H 1111 + H 1122 , h 3 = √ 2H 3311 , h 4 = √ 2H 1133 , h 5 = H 1111 -H 1122 , h 6 = 2H 2323
H E 0 = Q E 0 -1 , Q E 0 = 6 i=1 q i E i , P E 0 = 6 i=1 p i E i , C 0 = 6 i=1 c i E i , q 4 = q 3 , p 4 = p 3 , c 4 = c 3 ( 35
)
q 1 = c 1 -c 2 1 p 1 -c 2 3 p 2 -2 c 1 c 3 p 3 , q 2 = c 2 -c 2 3 p 1 -c 2 2 p 2 -2 c 2 c 3 p 3 ( 36
)
q 4 = q 3 = c 3 -c 1 c 3 p 1 -c 2 c 3 p 2 -c 1 c 2 + c 2 3 p 3 ( 37
)
q 5 = c 5 (1 -c 5 p 5 ) , q 6 = c 6 (1 -c 6 p 6 ) (38)
and components h i of compliance contribution tensor are deduced from calculation rules in transversely isotropic basis (with

q 4 = q 3 )   h 1 h 3 h 4 h 2   =   q 1 q 3 q 3 q 2   -1 = 1 q 1 q 2 -q 2 3   q 2 -q 3 -q 3 q 1   , h i = 1 q i = 1 c i (1 -c i p i ) (i = 5, 6) (39)
See similar derivations in [START_REF] Sevostianov | Effective elastic properties of the particulate composite with transversely isotropic phases[END_REF] , with different TI tensor basis (see formula 2.56).

A fourth order tensor T may be transversely isotropised by projection onto a transversely isotropic basis to obtain Π T I (T)

Π T I (T) = 2 i=1 (E i :: T) E i + (E 4 :: T) E 3 + (E 3 :: T) E 4 + 1 2 6 i=5 (E i :: T) E i (40)
In case of a symmetry axis is x 3 , we apply the transformation Q α representing a rotation of angle ϕ about x 3 axis to fourth order tensor T and non zero components of matrix Q α write

T α i jkl = Q α ip Q α jq Q α kr Q α ls T pqrs (
Q α 11 = Q α 22 = cos(ϕ), Q α 12 = -Q α 21 = sin(ϕ), Q α 33 = 1 (43)
Transverse isotropic projection Π T I (T) (40) corresponds to the following average over orientations, in the transverse plane

x 1 -x 2 Π T I (T) i jkl = 1 2 π 2π 0 T α i jkl (ϕ) dϕ = 1 2 π 2π 0 Q α ip (ϕ) Q α jq (ϕ) Q α kr (ϕ) Q α ls (ϕ) T pqrs dϕ (44)
Compliance contribution tensor of a superspherical pore aligned with the directions of symmetry of a transversely isotropic matrix (with x 3 symmetry axis of matrix) is tetragonal, with 6 independent components H 1111 , H 1122 , H 1133 , H 3333 , H 1212 , and H 2323 . Average compliance contribution tensor related to an isotropic orientation distribution of superspherical pores in the transverse isotropic plane x 1x 2 is transversely isotropic and obtained by transverse isotropic projection

h T I 1 = H 3333 , h T I 2 = H 1111 + H 1122 , h T I 3 = h T I 4 = √ 2 H 1133 h T I 5 = H 1111 -H 1122 + 2H 1212 2 , h T I 6 = 2H 2323 (45)

Appendix. Background on property contribution tensors

Property contribution tensors are used in micromechanics to describe the contribution of a single inhomogeneity to the property of interest ( [START_REF] Kachanov | Micromechanics of materials, with applications[END_REF]). Compliance contribution tensors have been first introduced in the context of pores and cracks by [START_REF] Hori | Overall moduli of solids with microcracks: Load induced anisotropy[END_REF] as the extra average strain produced by a pore. The average strain, over representative volume | Ω |, can be represented as a sum

ε = S 0 : Σ + ∆ε (48)
where S 0 where is the compliance tensor of the matrix and Σ is uniform remotely applied stress. The material is assumed to be a linear elastic; hence, the extra strain ∆ε due to presence of an inhomogeneity E is a linear function of the applied stress:

∆ε = f H E 0 : Σ, with f = | E | | Ω | (49) γ 1 = C 1111 C 3333 -C 2 1133 -2 C 1133 C 2323 + √ ∆ 2 C 2323 C 3333 ; γ 2 = C 1111 C 3333 -C 2 1133 -2 C 1133 C 2323 - √ ∆ 2 C 2323 C 3333 (57) with ∆ = (C 1111 C 3333 -C 2 1133 ) (C 1111 C 3333 -C 2 1133 -4 C 1133 C 2323 -4 C 2 2323 ) (58)
The square root of a complex argument is defined with a positive real part and it is consistent with the cmath library of Python. The inverse hyperbolic cosine of a complex argument, denoted by arcosh , has one branch cut, extending left from 1 along the real axis to -∞, continuous from above. Note that arcosh Z = ln (Z + √ Z 2 -1) where the principal value of the logarithm is chosen such that the imaginary part has the smallest value and belongs to ] -π 2 , π 2 ] with the same branch cut as the square root.

I k (η) = 1 z=-1 z k z 2 + η 2 (1 -z 2 ) dz and J k (η 1 , η 2 ) = 1 z=-1 z k z 2 + η 2 1 (1 -z 2 ) z 2 + η 2 2 (1 -z 2 ) dz ( 59
)
I 0 (η) I 2 (η) I 4 (η) if η = 1 2 2 3 2 5 if η 1 2 arcosh η η √ η 2 -1 2 η arcosh η- √ η 2 -1 (η 2 -1) 3 2 2 3 3 η 3 arcosh η+(1-4 η 2 ) √ η 2 -1 (η 2 -1) 5 2 Table 5: Calculation of integrals I k (η) J 0 (η 1 , η 2 ) J 2 (η 1 , η 2 ) J 4 (η 1 , η 2 ) if η 1 = η 2 = 1 2 2 3 2 5 if η 1 = η 2 1 arcosh η 1 +η 1 √ η 2 1 -1 η 3 1 √ η 2 1 -1 η 1 √ η 2 1 -1-arcosh η 1 η 1 (η 2 1 -1) 3
2 (2+η 2 1 ) √ η 2 1 -1-3 η 1 arcosh η 1 (η 2 1 -1) 5 2 if η 1 η 2 (η 2 1 -1) I 0 (η 1 )-(η 2 2 -1) I 0 (η 2 ) η 2 2 -η 2 1 (η 2 1 -1) I 2 (η 1 )-(η 2 2 -1) I 2 (η 2 ) η 2 2 -η 2 1 (η 2 1 -1) I 4 (η 1 )-(η 2 2 -1) I 4 (η 2 ) η 2 2 -η 2 1 Table 6: Calculation of integrals J k (η)
Related compliance contribution of the spheroidal pore is then deduced from relations [START_REF] Walpole | Fourth-rank tensors in the thirty-two crystal classes: multiplication tables[END_REF][START_REF] Hori | Overall moduli of solids with microcracks: Load induced anisotropy[END_REF][START_REF] Sevostianov | Explicit cross-property correlations for anisotropic two-phase composite materials[END_REF][START_REF] Sevostianov | Effective viscoelastic properties of short-fiber reinforced composites[END_REF][START_REF] Markov | On the applicability of replacement relations to tetrahedron-like inhomogeneities[END_REF]. The exact solution of compliance contribution tensor of the penny shaped crack, H = lim ω→0 ω (Q 0 ) -1 , is detailed in [START_REF] Barthélémy | Simplified approach to the derivation of the relationship between Hill polarization tensors of transformed problems and applications[END_REF]. We only recall final results. One uses

                                     p 1 1 = π 2 (γ 1 + γ 2 ) C 1111 -(γ 2 1 + γ 1 γ 2 + γ 2 2 ) C 2323 C 2323 C 3333 p 1 2 = π 4 (γ 1 + γ 2 ) 1 γ 1 γ 2 C 3333 + 1 C 2323 p 1 3 = - √ 2 π 4 (γ 1 + γ 2 ) C 1133 + C 2323 C 2323 C 3333 p 1 5 = π 8 1 γ 1 γ 2 (γ 1 + γ 2 ) C 3333 + 1 γ 1 + γ 2 + 1 γ 3 1 C 2323 p 1 6 = π 8 C 1111 -2 γ 1 γ 2 C 1133 -γ 1 γ 2 (γ 2 1 + γ 1 γ 2 + γ 2 2 ) C 3333 γ 1 γ 2 (γ 1 + γ 2 ) C 2323 C 3333 - γ 3 C 2323 (60a) (60b) (60c) (60d) (60e)
γ 1 , γ 2 , γ 3 have been previously defined. Compliance contribution tensor of a penny shaped crack aligned in a transversely isotropic host matrix writes

H = lim ω→0 ω (Q 0 ) -1 = - 1 C 2 3333 p 1 1 + 2 C 2 1133 p 1 2 + 2 √ 2 C 3333 C 1133 p 1 3 E 1 - 1 4 C 2 2323 p 1 6 E 6 (61)
where p 1 1 , p 1 2 , p 1 3 and p 1 6 are given in (60a), (60b), (60c) and (60e). One deduces the components denoted h c 1 , h c 6 and used in approximate relations of compliance contribution tensor of an aligned axisymmetric superspheroidal pore: Finite element results of a superspherical pore embedded in a TI matrix, are given in table 4 of paper [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF], and recalled in table [START_REF] Drach | Comparison of full field and single pore approaches to homogenization of linearly elastic materials with pores of regular and irregular shapes[END_REF]. Finite element results of an axisymmetrical superspheroidal pore embedded in a TI matrix are given in table 3 of paper [START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF] and recalled in table [START_REF] Trofimov | Effective elastic properties of composites with particles of polyhedral shapes[END_REF]. Approximate relations for the concavity parameter range 0.2 < p < 1 write f so-a 1111 (p) = 54.285p 4 -85.8182p 3 + 44.7893p 2 -6.31294p + 0.294191 f so-a 1122 (p) = 50.5731p 4 -84.3085p 3 + 46.7634p 2 -7.38197p + 0.378348 f so-a 1133 (p) = 39.0659p 4 -42.9843p 3 + 7.

h c 1 = - 1 C 2 3333 p 1 1 + 2 C 2 1133 p 1 2 + 2 √ 2 C 3333 C 1133 p 1 3 , h c 6 = -
p
p

3 Figure 1 :

 31 Figure 1: Superspherical pore (relation (1) with ς = 1)

Figure 2 :

 2 Figure 2: 2D representation in diametral plane of a 3D axisymmetrical superspheroidal pore, with ς = 1 and symmetry axis x 3

1 Figure 8 : 1 ( 3 (

 1813 Figure 8: Effective transverse E MTB 1 (top, left) and normal E MTB 3 (top, right) Young's elastic moduli, effective axial shear modulus G MTB 31 as a function of porosity ϕ, MTB approximation, aligned axisymmetrical superspheroidal pores. App 1 (plain lines): approximation formula (relations 10-11-12) for axisymmetrical spheroidal pores, App 2 (dashed lines) : approximation oblate spheroid with same volume (relations 18-23). Note that the two approximations coincide only for normal Young's modulus E MTB 3 .

1 Figure 9 : 1 ( 3 (

 1913 Figure 9: Effective transverse E MTB 1 (top, left) and normal E MTB 3 (top, right) Young's elastic moduli , effective shear coefficient G MTB 31

e = 1 Figure 10 :e = 1 Figure 11 :

 110111 Figure 10: Effective transverse E ef 1 (top, left) and normal E ef 3 (top, right) Young's moduli, effective shear coefficient G ef 31 (bottom) as a fonction of porosity ϕ, for aligned axisymmetric superspheroidal pores randomly oriented in the isotropic transverse plane.

41) with e α 1 =

 1 cos(ϕ)e 1 + sin(ϕ) e 2 , e α 2 =sin(ϕ)e 1 + cos(ϕ) e 2 , e α 3 = e 3 (42)

9 .

 9 Appendix: numerical results for approximation formula of compliance contribution tensors

  ) = -0.197702p 2 + 0.220097p + 0.939377 f se-b 1122 (p) = -0.24264p 2 + 0.0415192p + 1.0399 f se-b 1133 (p) = -0.357454p 2 + 0.52452p + 0.827103 f se-b 3333 (p) = -0.326089p 2 + 0.439731p + 0.861657 f se-b 1212 (p) = -2.04372p 4 + 7.78585p 3 -11.3964p 2 + 7.68465p -0.839061 f se-b 2323 (p) = -2.65601p 4 + 9.6441p 3 -13.4224p 2 + 8.59804p -0.982841 (64)

  43095p 2 + 5.44332p -0.646506 (65) f so-b 1111 (p) = -3.74036p 4 + 13.3244p 3 -18.0592p 2 + 11.0818p -1.45743 f so-b 1122 (p) = -4.07865p 4 + 14.7289p 3 -20.227p 2 + 12.5889p -1.82352 f sp-b 1133 (p) = -1.92702p 4 + 6.39843p 3 -7.8095p 2 + 3.95979p + 0.293327 (66)

Table 1 :

 1 Reference transversely isotropic elastic parameters

	E 0 1 (GPa) E 0 3 (GPa)	ν 0 12	ν 0 31	G 0 31 (GPa)
		20.44	11.306		0.1027 0.1798	1.5851
	Table 2: Reference transversely isotropic elastic parameters: C 0 i jkl components of C 0 tensor and related c i components in transversely isotropic
	tensor basis E i
	C 0 1111 (GPa) C 0 3333 (GPa)	C 0 1122 (GPa)	C 0 1133 (GPa) C 0 2323 (GPa)
	22.3639	12.9994		3.8275	4.7092	1.5851
	c 1 (GPa)	c 3 (GPa)	c 3 = c 4 (GPa)	c 5 (GPa)	c 6 (GPa)
	12.9994	26.1914		6.65983	18.5363	3.1702

Table 3 :

 3 and finite element results are presented in figure 4. Maximal relative errors of approximate relations (9) compared to FEM results are given in table 3, they are lower than 4.%. Maximal relative errors of approximate relations compared to FEM results, (H FEM i jkl -H FEM i jkl )/H FEM

	H E 1111	H E 1122	H E 1133	H E 3333	H E 1212	H E 1313
	0.0349 0.03677 0.0350 0.0356 0.0061 0.0060
	p ∈ [0.2, 1]					i jkl	∞	for the superspherical pore with

Table 4 :

 4 -11-12-17) and finite element results are presented in figure5. Compliance contribution tensor of aligned axisymmetrical superspheroidal and spheroidal pores of same volume have quasi the same normal component H 3333 . In other words, only the volume characterizes this component, not affected by concavity (superspheroid) of convexity (spheroid). It must be emphasized that it is not the case for all the other components H i jkl , including the shear component H 2323 , for which the concavity parameter p is of major importance (volume is not sufficient to characterize compliance contribution tensor). Shear component H 1212 in the plane of transverse isotropy has been used to check accuracy of the symmetry of revolution by comparing to H 1111 -H 1122 )/2. (see figure6), it may be noticed that both coincide as expected.Maximal relative errors of approximate relations (10-12-17) compared to FEM results are given in table4, they are lower than 5.%. Maximal relative errors of approximate relations compared to FEM results, (H FEM Figure 6: H 1212 component of the T I compliance contribution tensor of an axisymmetric superspheroidal pore embedded in TI matrix, as a function of concavity parameter p. H 1212 is used to check accuracy of transverse isotropy by comparing to (H 1111 -H 1122 )/2..

	H E 1111	H E 1122	H E 1133	H E 3333	H E 1313
	0.00133 0.00174 0.00297 0.04890 0.03425
					i jkl -H FEM i jkl )/H FEM i jkl	∞	for the axisymetrical super-
	spheroidal pore with p ∈ [0.2, 1]

includes the particular case of the sphere, H sphere i jkl with γ(p = 1) = 1) is deduced from the exact Hill tensor recalled in appendix 8. Comparisons between approximate relations (

10

Figure

5

: The 5 independent components H i jkl of the T I compliance contribution tensor of an axisymmetric superspheroidal pore embedded in TI matrix, as a function of concavity parameter p. Comparison between FEM results (dashed lines) and approximate relations (plain lines). Note that H 1212 is used to check accuracy of transverse isotropy by comparing to (H 1111 -H 1122 )/2..

  See appendix for detail, only components H 1111 , H 1122 and H 1212 are modified, other components are equal,

				H T I 1212 =	H 1111 -H 1122 + 2 H 1212 4		(22)
	Figure (7) illustrates the numerical difference between the tetragonal tensor H E 0 and its TI projection Π T I H E 0 .
			supersphere: H tetragonal						H -Π T I (H) / H
		0.35						0.1
				H 1212
		0.3		(H 1111 -H 1122 )/2	0.08
		0.25						0.06
	H 1212	0.2
								0.04
		0.15
								0.02
		0.1
								0
		0.2	0.4	0.6 p	0.8	1		0.2	0.4	0.6 p	0.8	1
	H T I 1133 = H 1133 , H T I 3333 = H 3333 , H T I 2323 = H 2323 (H T I i jkl denotes Π T I H E 0 i jkl )
			H T I 1111 =	3 H 1111 + H 1122 + 2 H 1212 4	, H T I 1122 =	H 1111 + 3 H 1122 -2 H 1212 4	(21)

Figure 7: Left figure: comparison between H 1212 component and (H 1111 -H 1122 )/2 to illustrate tetragonal symmetry of H tensor. Right figure:

Table 7 :

 7 Numerical estimation of H E i jkl for the superspheroidal pore embedded in a transversely isotropic corrected model with different values of concavity p ∈ [0.2, 1], from[START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF].Approximate relations for the concavity parameter range 0.2 < p < 1 writef se-a 1111 (p) = 2.72394p 2 -2.6248p + 1.63141 f se-a 1122 (p) = 1.23332p 2 -1.20253p + 1.29293 f se-a 1133 (p) = 3.75077p 2 -3.50795p + 1.81628 f se-a 3333 (p) = 2.77572p 2 -2.46747p + 1.5398 f se-a 1212 (p) = 17.1508p 4 -43.7713p 3 + 38.3878p 2 -13.1282p + 2.36576 f se-a 2323 (p) = 50.3296p 4 -90.0984p 3 + 59.3072p 2 -15.8086p + 2.1938

		H E 1111	H E 3333	H E 1122	H E 1133	H E 1212	H E 1313
	0.2 0.5401 0.9682 -0.1386 -0.1093 0.2043 0.4972
	0.25 0.3391 0.6106 -0.0895 -0.0676 0.1288 0.3327
	0.3 0.2521 0.4567 -0.0682 -0.0496 0.0984 0.2715
	0.35 0.2052 0.3749 -0.0567 -0.0400 0.0844 0.2454
	0.4 0.1770 0.3269 -0.0497 -0.0343 0.0779 0.2330
	0.45 0.1589 0.2967 -0.0450 -0.0308 0.0749 0.2265
	0.5 0.1460 0.2747 -0.0413 -0.0284 0.0735 0.2229
	0.6 0.1337 0.2542 -0.0371 -0.0263 0.0727 0.2199
	0.7 0.1269 0.2419 -0.0342 -0.0252 0.0728 0.2194
	0.8 0.1229 0.2345 -0.0321 -0.0246 0.0731 0.2200
	0.9 0.1204 0.2297 -0.0305 -0.0242 0.0735 0.2212
	1	0.1188 0.2266 -0.0293 -0.0239 0.0740 0.2226

Table 8 :

 8 Numerical estimation of H E i jkl for the superspheroidal pore embedded in a transversely isotropic corrected model with different values of concavity p ∈ [0.2, 1], from[START_REF] Du | Numerical computation of compliance contribution tensor of a concave pore embedded in a transversely isotropic matrix[END_REF].

		H E 1111	H E 1122	H E 1133	H E 3333	H E 1212	H E 1313	H E 1111 -H E 1122 2	Error 1212	1
	0.2 0.1004 -0.0192 -0.0520 3.6670 0.0597 0.8424	0.0598	0.16%
	0.25 0.1026 -0.0212 -0.0442 1.5948 0.0619 0.4132	0.0619	0.04%
	0.3 0.1068 -0.0232 -0.0395 0.9354 0.0650 0.2928	0.0650	0.01%
	0.35 0.1102 -0.0248 -0.0359 0.6445 0.0675 0.2495	0.0675	0.00%
	0.4 0.1128 -0.0260 -0.0331 0.4929 0.0694 0.2324	0.0694	0.00%
	0.45 0.1147 -0.0270 -0.0310 0.4047 0.0708 0.2248	0.0708	0.00%
	0.5 0.1159 -0.0276 -0.0293 0.3492 0.0718 0.2211	0.0718	0.00%
	0.6 0.1175 -0.0284 -0.0271 0.2900 0.0730 0.2187	0.0730	0.01%
	0.7 0.1182 -0.0288 -0.0258 0.2604 0.0735 0.2188	0.0735	0.01%
	0.8 0.1186 -0.0291 -0.0249 0.2436 0.0738 0.2197	0.0738	0.01%
	0.9 0.1187 -0.0292 -0.0243 0.2333 0.0739 0.2211	0.0739	0.01%
	1	0.1187 -0.0293 -0.0239 0.2265 0.0740 0.2226	0.0740	0.01%
	1 Relative error of	H E 1111 -H E 1122 2	with respect to H E 1212

References References

where | E | is the pore volume and H E 0 is fourth-rank compliance contribution tensor of the pore. The H E 0 tensor is determined by the shape and size of the inhomogeneity, as well as properties of the matrix and of the inhomogeneity material. This tensor is also affected by elastic interactions. In the non-interaction approximation, it is taken by treating the inhomogeneities as isolated ones. In the case of multiple inhomogeneities, the extra strain produced by m-th inhomogeneity is ∆ε (m) = f (m) H E(m) 0 : Σ so that the extra compliance due to all the inhomogeneities is given by

Formula (50) highlights the fundamental importance of the compliance contribution tensors: these tensors have to be summed up and averaged in the context of the effective elastic properties. The sum

properly reflects compliance contributions of individual inhomogeneities and constitutes the general microstructural parameters in whose terms the effective compliance should be expressed. Components of this tensor were calculated for 2 -D pores of various shape and 3 -D ellipsoidal pores in isotropic material by [START_REF] Kachanov | Effective Moduli of Solids With Cavities of Various Shapes[END_REF]. For the general case of elastic inhomogeneities, these tensors were introduced and calculated (for ellipsoidal shapes) by Sevostianov and Kachanov ([37]). Components of the compliance contribution tensor for various concave pores in isotropic matrix have been calculated by [START_REF] Chen | Evaluation of the effective elastic and conductive properties of a material containing concave pores[END_REF] (supersphere), [START_REF] Sevostianov | Effective viscoelastic properties of short-fiber reinforced composites[END_REF] (axisymmetric pore obtained by rotation of a supersphere around one of its principal diagonals), [START_REF] Chen | Combined effect of pores concavity and aspect ratio on the elastic properties of a porous material[END_REF] (combined effect of concavity and aspect ratio), and [START_REF] Markov | On the applicability of replacement relations to tetrahedron-like inhomogeneities[END_REF] (tetrahedron-like pores). [START_REF] Sevostianov | Effective elastic properties of the particulate composite with transversely isotropic phases[END_REF][START_REF] Barthélémy | Simplified approach to the derivation of the relationship between Hill polarization tensors of transformed problems and applications[END_REF] calculated components of this tensor for a spheroidal inhomogeneity embedded in a transversely-isotropic material. We recall the compact solution of strain Hill polarization tensor of a spheroidal inclusion aligned in a transversely isotropic host matrix recently presented in [START_REF] Barthélémy | Simplified approach to the derivation of the relationship between Hill polarization tensors of transformed problems and applications[END_REF] (reader may refer to this reference for the detailed derivation and python script of the complete solution). The corresponding compliance contribution tensor may be deduced from relations presented in section [START_REF] Drach | Comparison of full field and single pore approaches to homogenization of linearly elastic materials with pores of regular and irregular shapes[END_REF]. In what follows, the aspect ratio of the spheroidal inclusion is denoted ω, and the symmetry axis n of the host matrix and the spheroidal inclusion is taken as n = e 3 . The particular case of the spherical pore is then deduced by setting ω = 1. The solution writes in Walpole TI tensor basis

with


Marine Schott

Laurent Ranno

Hélène Béa

Claire Baraduc

Stéphane Auffret

Anne Bernand-Mantel

Electric field control of interfacial Dzyaloshinskii-Moriya interaction in Pt/Co/AlOx thin films

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

I. INTRODUCTION

Electric-field (EF) control of magnetism in metals is currently an actively expanding field of study, motivated by potential applications in low power spintronics devices. The magnetization of metallic magnetic materials in spintronics devices can be manipulated using magnetic field or spin-polarized currents flow. The ability to control the magnetic properties of a system using a gate voltage is a promising technique which could lead to the development of EF-assisted magnetization switching devices. Since the first experimental demonstration of EF control of coercive field in a FePt thin film by Weisheit et al. [START_REF] Weisheit | Electric field-induced modification of magnetism in thin-film ferromagnets[END_REF], lots of studies have been conducted with the final intention of finding the best conditions for functional devices. Significant variations in the magnetic anisotropy energy under EF application have been observed by many groups [START_REF] Maruyama | Large voltage-induced magnetic anisotropy change in a few atomic layers of iron[END_REF][START_REF] Wang | Electric-field-assisted switching in magnetic tunnel junctions[END_REF][START_REF] Niranjan | Electric field effect on magnetization at the Fe/MgO(001) interface[END_REF][START_REF] Nakamura | Giant Modification of the Magnetocrystalline Anisotropy in Transition-Metal Monolayers by an External Electric Field[END_REF][START_REF] Endo | Electric-field effects on thickness dependent magnetic anisotropy of sputtered MgO/CoFeB/Ta structures[END_REF][START_REF] Yamada | Electric Field Modulation of Magnetic Anisotropy in MgO/Co/Pt Structure[END_REF], as well as EF variation of T c in ultrathin films [START_REF] Ohno | Electricfield control of ferromagnetism[END_REF][START_REF] Chiba | Electrical control of the ferromagnetic phase transition in cobalt at room temperature[END_REF][START_REF] Shimamura | Electrical control of Curie temperature in cobalt using an ionic liquid film[END_REF]. Recently, the EF impact on the exchange stiffness parameter A ex has been addressed experimentally [START_REF] Ando | Modulation of the magnetic domain size induced by an electric field[END_REF][START_REF] Dohi | Effect of electric-field modulation of magnetic parameters on domain structure in mgo/cofeb[END_REF] and theoretically [START_REF] Oba | Electric-Field-Induced Modification of the Magnon Energy, Exchange Interaction, and Curie Temperature of Transition-Metal Thin Films[END_REF]. In these systems, a dielectric oxide layer is present at the top surface of the ferromagnet to allow the application of an EF in a capacitor geometry. This induces a broken spatial inversion symmetry for the ferromagnetic thin film. The presence of these non-identical interfaces is the source of an antisymmetric type of exchange, the Dzyaloshinskii-Moriya interaction (DMI), which has been shown to be a source of chiral spin textures in Pt/Co/AlO x for instance [START_REF] Belmeguenai | Interfacial dzyaloshinskii-moriya interaction in perpendicularly magnetized ultrathin films measured by brillouin light spectroscopy[END_REF]. DMI has been recently shown to be tuneable using strain [START_REF] Gusev | Manipulation of dzyaloshinskii-moriya interaction in co/pt multilayers with strain[END_REF], or using the current flowing in the heavy metal [START_REF] Karnad | Modification of dzyaloshinskii-moriyainteraction-stabilized domain wall chirality by driving currents[END_REF][START_REF] Kato | Currentinduced modulation of the interfacial dzyaloshinskii-moriya interaction[END_REF] . The interfacial origin of this DMI and its consequences in spin textures in ultrathin films make also DMI an interesting EF tunable parameter. However, distinguishing EF modulation on DMI from other contributions is not straightforward as direct access to DMI is difficult. In addition to theoretical predictions for ultra-thin samples [START_REF] Yang | Controlling Dzyaloshinskii-Moriya Interaction via Chirality Dependent Atomic-Layer Stacking, Insulator Capping and Electric Field[END_REF], several experimental studies have reported a change in the DMI factor under EF application [START_REF] Yang | Controlling Dzyaloshinskii-Moriya Interaction via Chirality Dependent Atomic-Layer Stacking, Insulator Capping and Electric Field[END_REF][START_REF] Nawaoka | Voltage induction of interfacial DzyaloshinskiiMoriya interaction in Au/Fe/MgO artificial multilayer[END_REF][START_REF] Srivastava | Large-Voltage Tuning of Dzyaloshinskii-Moriya Interactions: A Route toward Dynamic Control of Skyrmion Chirality[END_REF][START_REF] Koyama | Electric field control of magnetic domain wall motion via modulation of the Dzyaloshinskii-Moriya interaction[END_REF][START_REF] Zhang | Electrical field enhanced interfacial Dzyaloshinskii-Moriya interaction in MgO/Fe/Pt system[END_REF][START_REF] Suwardy | Voltage-controlled magnetic anisotropy and Dzyaloshinskii-Moriya interactions in CoNi/MgO and CoNi/pd/MgO[END_REF]. This control was shown for relatively thick layers of Fe (20 nm) [START_REF] Nawaoka | Voltage induction of interfacial DzyaloshinskiiMoriya interaction in Au/Fe/MgO artificial multilayer[END_REF][START_REF] Suwardy | Voltage-controlled magnetic anisotropy and Dzyaloshinskii-Moriya interactions in CoNi/MgO and CoNi/pd/MgO[END_REF],

or a material presenting a weak DMI, Ta/FeCoB/TaO x , as demonstrated by Srivastava et al. [START_REF] Srivastava | Large-Voltage Tuning of Dzyaloshinskii-Moriya Interactions: A Route toward Dynamic Control of Skyrmion Chirality[END_REF]. More recently, Zhang et al. [START_REF] Zhang | Electrical field enhanced interfacial Dzyaloshinskii-Moriya interaction in MgO/Fe/Pt system[END_REF] and Koyama et al. [START_REF] Koyama | Electric field control of magnetic domain wall motion via modulation of the Dzyaloshinskii-Moriya interaction[END_REF] showed larger EF induced variation of DMI in Pt/Fe/MgO and Pt/Co/Pd/MgO systems presenting intermediate

DMI values.

In this work, we analyzed labyrinthine domains (also referred to as stripe domains) in samples of Pt/Co/AlO x presenting large DMI values [START_REF] Belmeguenai | Interfacial dzyaloshinskii-moriya interaction in perpendicularly magnetized ultrathin films measured by brillouin light spectroscopy[END_REF][START_REF] Vanatka | Velocity asymmetry of dzyaloshinskii domain walls in the creep and flow regimes[END_REF]. A reversible evolution of the labyrinthine magnetic domain configurations under EF is observed. We also report an EF variation of the saturation magnetization and anisotropy field. The analysis of the variation of equilibrium stripe width allowed us to estimate the EF variation of the DMI term D. We stress here that this estimation strongly depends on the assumptions made on the exchange constant value.

II. SYSTEM DESCRIPTION AND CHARACTERIZATION

Our system is a Pt(3 nm)/Co(0.6 nm)/AlO x sputtered trilayer, presenting a gradient of oxidation at the Co/AlO x interface. This was induced by the post-oxidation of a wedgedshaped Al top layer (Fig. 1(a)). Using this technique, we get access to a low thicknesses range for Co due to the Co partial oxidation in the region where the deposited Al is thin.

The sample has been covered by about 50 nm of a high-k insulator, HfO 2 , deposited using atomic layer deposition. Then an Indium Tin Oxide (ITO) layer was DC sputtered and patterned in electrodes. These transparent electrodes allow to obtain magneto-optical Kerr effect (MOKE) images of the magnetic domains in the Co layer under a gate voltage (Fig. 1(d)). Voltages of different magnitudes were applied between this top electrode and the Co layer. Here a positive voltage application induces electrons accumulation at the surface of the Co layer (Fig. 1(a)). We used MOKE in polar geometry and recorded images and magnetic hysteresis loops on a single position on the sample through the ITO electrode under different applied voltages (Fig. 1(d)). To determine the EF variation of the saturation magnetization M s , we measured the variation of the p-MOKE signal between the two opposite saturated states, for out-of-plane applied magnetic field. The average value of magnetization of the wedged sample has been obtained by measuring M s with a Vibrating Sample Magnetometer-Superconducting Quantum Interference Device (VSM-SQUID). A T c of 366 K was estimated from VSM-SQUID measurements on the sample under study which possesses an averaged Co thickness of 0.49 nm due to partial oxidation of the deposited Co (see Fig. 1

(a)).

A strong EF dependence of M s is shown in Fig. 1(e). Such sensitivity of the saturation magnetization is due to the proximity of the T c to the measurement temperature (room temperature). Similar EF control of T c associated with strong M s variations were observed in recent studies [START_REF] Ando | Modulation of the magnetic domain size induced by an electric field[END_REF][START_REF] Chiba | Control of magnetism in Co by an electric field[END_REF]. The EF variation of the anisotropy field µ 0 H a is presented in Fig. 1(f). The anisotropy field µ 0 H a is extracted using Stoner-Wohlfarth fit of the magnetization rotation toward the hard axis direction, on p-MOKE hysteresis loops, for in-plane applied magnetic field (see Fig. 1(b)). From the measurements of µ 0 H a (V ), we deduced the variations of the effective magnetic anisotropy energy from the formula µ 0 H a = 2K eff /M s ,

where

K eff = K s /t -1 2 µ 0 M 2
s is the effective magnetic anisotropy and K s is the surface magneto cristalline anisotropy. A voltage controlled magnetic anisotropy (VCMA) coefficient [START_REF] Dieny | Perpendicular magnetic anisotropy at transition metal/oxide interfaces and applications[END_REF] β PMA = ∆K s /∆E = 562 ± 102 fJ/(V•m) is obtained. This is larger compared to the usual range obtained for a charge accumulation effect (10-290 fJ/(V•m)) [START_REF] Wang | Electric-field-assisted switching in magnetic tunnel junctions[END_REF][START_REF] Endo | Electric-field effects on thickness dependent magnetic anisotropy of sputtered MgO/CoFeB/Ta structures[END_REF][START_REF] Kita | Electric-field-control of magnetic anisotropy of Co0.6Fe0.2B0.2/oxide stacks using reduced voltage[END_REF][START_REF] Nozaki | Large voltage-induced changes in the perpendicular magnetic anisotropy of an mgo-based tunnel junction with an ultrathin fe layer[END_REF]. The enhanced effect is again likely due to the proximity of T c to the measurement temperature.

III. EF VARIATION OF EQUILIBRIUM LABYRINTHINE DOMAINS

As discussed earlier, the analysis of the equilibrium domain configuration has been shown to be one of the ways to study the influence of the EF on different magnetic properties [START_REF] Ando | Modulation of the magnetic domain size induced by an electric field[END_REF][START_REF] Dohi | Effect of electric-field modulation of magnetic parameters on domain structure in mgo/cofeb[END_REF]. At zero and negative gate voltages, the domains are blurry due to a very strong thermally induced domain wall motion, as observed in previous works [START_REF] Bergeard | Dynamic fluctuations and two-dimensional melting at the spin reorientation transition[END_REF]. Thus, in the following we focus our study on positive gate voltages, for which the labyrinthine domains are stable, between 6 V and 14 V. To record the images of labyrinthine states and their variation with the EF, we proceed as follow: under a constant EF, we first saturate the sample with an applied magnetic field, then we turn off the magnetic field and record an image after a waiting time of a few seconds. The images presented in Fig. 2(a)-(d) correspond to the labyrinthine domains which spontaneously appear due to thermally activated nucleation and domain wall propagation. The labyrinthine domains are clearly and reversibly influenced by the application of an EF, similarly to what we observed in our previous experiments [START_REF] Schott | The skyrmion switch: Turning magnetic skyrmion bubbles on and off with an electric field[END_REF]. The characteristic labyrinthine domain width L eq is estimated using a 2D Fourier transform of the images in Fig. 2(a)-(d). In Fig. 2(e) we present L eq as function of the applied voltage, which shows an increase as a positive voltage is applied to the top ITO electrode. To understand in which way the EF is influencing the domain periodicity, we used an analytical model to describe the domain width L eq in ultrathin films:

L eq = C • t • exp πσ ω 2K d t , (1)
0 5 1 0 1 5 0 2 4 ( f ) ( e ) 1 4 V 1 2 V 6 V ( d ) ( c ) ( b ) L e q ( µm )
V o l t a g e ( V ) where t is the cobalt thickness (which is fixed here to 0.49 nm), K d the dipolar constant [START_REF] Kaplan | The domain structure in ultrathin magnetic films[END_REF][START_REF] Gehanno | Magnetic susceptibility and magnetic domain configuration as a function of the layer thickness in epitaxial fepd(0 0 1) thin films ordered in the l1o structure[END_REF] of order unity. The experimental domain wall energy σ ω in Fig. 2(f) is deduced using the measured L eq and M s values. We find domain wall energies σ ω around 1.3 -1.6 mJ/m 2 for the range of voltages under study, which is consistent with earlier studies [START_REF] Schott | The skyrmion switch: Turning magnetic skyrmion bubbles on and off with an electric field[END_REF]. In the presence of sufficiently strong DMI, the domain wall is of Néel type with a sense of rotation which lowers the domain wall energy [START_REF] Heide | Dzyaloshinskii-moriya interaction accounting for the orientation of magnetic domains in ultrathin films: Fe/w(110)[END_REF]:

K d = 1 2 µ 0 M 2 s and C is a model-dependent constant
σ w = σ Bloch w + σ DM I w = 4 A ex K eff -πD,
where D = D s /t stands for the bulk DMI constant and is expressed in J/m 2 , while the surface DMI value D s is expressed in J/m. We see that in order to extract D using the measured σ w and K ef f values it is necessary to know the exchange constant A ex . We have used two values for A ex . The first value A ex =7.5 pJ/m (for V=0 V) was deduced from a fit using Kuzmin formula [START_REF] Kuz | min, Shape of Temperature Dependence of Spontaneous Magnetization of Ferromagnets: Quantitative Analysis[END_REF] of the temperature dependence of M s in our sample measured by VSM-SQUID. This reduced value is coherent with recent studies [START_REF] Yastremsky | Thermodynamics and Exchange Stiffness of Asymmetrically Sandwiched Ultrathin Ferromagnetic Films with Perpendicular Anisotropy[END_REF]. The second value

A ex =16 pJ/m (for V=0 V) is a value corresponding to bulk thin Co films [START_REF] Metaxas | Creep and flow regimes of magnetic domain-wall motion in ultrathin Pt/Co/Pt films with perpendicular anisotropy[END_REF]. The fact that the value of A ex extracted from our VSM-SQUID measurements is lower than the bulk Co value is linked to the low T c of our sample (366 K) and coherent with the low M s value.

For the same reason we can expect A ex to be influenced by EF as it was already suggested by previous studies [11-13, 37, 38]. We present the different deduced energy terms and their EF variation in Fig. 3. The total domain wall energy σ w , also shown in Fig. 2 A ex discussed ealier. A ex is considered to vary with EF as ∆A ex /A ex = 2∆M s /M s , as the general tendency described by mean field approximation is that A ex ∝ M 2 s near T c [START_REF] Atxitia | Multiscale modeling of magnetic materials: Temperature dependence of the exchange stiffness[END_REF]. Other theoretical studies have found a scaling law of A ex ∝ M 2- s , with being close to 0.2 for low temperatures [START_REF] Moreno | Temperature-dependent exchange stiffness and domain wall width in co[END_REF][START_REF] Kvashnin | Electronic topological transition and noncollinear magnetism in compressed hcp co[END_REF][START_REF] Turek | Ab initio theory of exchange interactions in itinerant magnets[END_REF]. Fixing ∆A ex /A ex = 2∆M s /M s is then a good approximation and will give us an upper limit for the variation we can expect for D(V ). From these assumptions we deduce the DMI value and EF variations which are presented in Fig. 3.

The D value at 6 V is ranging between 1.3-2.4 mJ/m 2 depending on the chosen A ex value.

We see that EF variation of σ Bloch w (Fig. 3(a)) is stronger than the variation of the total wall energy σ w (Fig. 3(b)) and that an EF variation of D is necessary to compensate this variation (Fig. 3(c)).

To discuss the DMI variations under EF, we use the voltage-control DMI (VCDMI) coefficient [START_REF] Yang | Controlling Dzyaloshinskii-Moriya Interaction via Chirality Dependent Atomic-Layer Stacking, Insulator Capping and Electric Field[END_REF][START_REF] Srivastava | Large-Voltage Tuning of Dzyaloshinskii-Moriya Interactions: A Route toward Dynamic Control of Skyrmion Chirality[END_REF]), defined as β DM I = ∆D/∆E (in J/(V•m), with ∆E = ∆V /t Ox , where ∆V is the voltage variation and t Ox the dielectric tickness. The coefficient β DMI is plotted versus A ex in Fig. 4. It is in the range of β DM I = 1100 -2000 ± 700 fJ/(V•m). In order to compare this result with previous works we give a summary of the different values in Table 1. origin, and consequently it is the variation of surface DMI constant D s which should be compared. In addition, as high-k oxides have been used in order to allow smaller EF for a similar effect on DMI we also have to take into account the effective permittivity of the (possibly multilayered) oxide and compare the effect for a given induced surface charge density and not a given electric field (assuming that the EF effect is induced by electron displacement). After performing this renormalization of the EF changes on DMI, we can now compare the η DM I values corresponding to different studies from several teams. The obtained values lie within the range 0.6 to 14.5 ×10 -12 J•m/C. We first see that the two studies with long time scale measurements (made by Brillouin light Spectroscopy (BLS) [START_REF] Srivastava | Large-Voltage Tuning of Dzyaloshinskii-Moriya Interactions: A Route toward Dynamic Control of Skyrmion Chirality[END_REF][START_REF] Zhang | Electrical field enhanced interfacial Dzyaloshinskii-Moriya interaction in MgO/Fe/Pt system[END_REF]) provide much larger η DM I . This is potentially due to the ion migration contribution from these long measurements. Long (days or hours) and short time (seconds or minutes) measurements present a factor of 5 ratio in η DM I [START_REF] Srivastava | Large-Voltage Tuning of Dzyaloshinskii-Moriya Interactions: A Route toward Dynamic Control of Skyrmion Chirality[END_REF].

FIG. 1 :

 1 FIG. 1: (a) Schematics of the wedge sample and of the set-up used for EF application: the Pt/Co/AlO x trilayer is covered by a ∼ 50 nm of dielectric layer (HfO 2 ) and a transparent top electrode (ITO). (b) p-MOKE hard axis (in-plane field) hysteresis loop (black squares) and Stoner-Wohlfarth fit (red line). (c) p-MOKE easy axis (out-of-plane field) hysteresis loop (black squares). (d) p-MOKE image of the demagnetized state taken at 0 T, after saturation with an out-of-plane field of 0.4 mT for +12 V applied gate voltage. (e) Evolution of saturation magnetization with gate voltage. A systematic error of 7%, was estimated from the error on the VSM-SQUID measurement. (f) Evolution of anisotropy field with voltage measured from p-MOKE hysteresis loops with inplane field. A systematic error of 10%, was estimated from fluctuations in the field calibration of the p-MOKE measurement with in-plane applied magnetic field. The blue solid lines in c and d show the regions where the EF-efficiencies have been extracted.

FIG. 2 :

 2 FIG. 2: (a)-(d) p-MOKE images of the evolution of equilibrium stripe width L eq with gate voltage. (e) Equilibrium stripe width L eq as a function of gate voltage. We estimated from repeated measurements a ∼ 10% error on L eq due to imperfect demagnetization. (f) Deduced domain wall energy as function of gate voltage. The error bars correspond to a 16% error calculated from the propagation of the errors in M s , t and L eq .

2 )FIG. 3 :

 23 FIG. 3: Different components of the wall energy σ w as a function of applied EF : (a) Bloch wall energy σ Bloch w

= 4 A

 4 (f) is presented in Fig. 3(b). In Fig. 3(a) we present two values and variations of the Bloch wall energy σ Bloch ω ex K eff deduced from the measured H a and M s and the two values of

FIG. 4 :

 4 FIG. 4: β DMI coefficients deduced from the experimental M s , t, L eq and K eff values as a function of the chosen A ex value. The error bars are calculated from the standard error on the slope on a linear fit of D(V ).

  the EF calculated corresponds to an applied EF or takes into account the effect of r in the oxide, which is still debated in the community, in particular due to the underestimation of the dielectric permittivity of MgO. These two values are thus giving ranges for theoretical η DM I in the Pt/Co/MgO system. The relatively small dispersion in the extracted η DM I values is noticeable, as all these studies have been performed in different teams with samples with different ferromagnetic materials, seed layers and oxides, grown with different techniques (magnetron sputtering or molecular beam epitaxy), measured with different methods (BLS, DW asymmetric motion under in-plane magnetic field, analysis of the labyrinthine domain structures, excitation of magneto-static surface spin-waves) with EF ranging from 18 to 800 MV/m. This indicates that the underlying physics is similar for all these samples despite the very different β DM I proposed. In order to carefully compare interface effect and optimize the EF effect on DMI, one should thus perform calculation of this renormalized, thickness and oxide-independent η DM I .ConclusionIn conclusion, we were able to deduce a strong EF variation of the magnetic anisotropy energy (β PMA = 562 ± 102 fJ/(V•m)) and of the DM I term (β DMI ∼ 1100 -2000 ± 700 fJ/(V•m)). The observed variation of domain size L eq is a result of the combination of variations of every magnetic parameter. The strong β coefficients we were able to deduce for PMA and DMI energies are the result of significantly increased EF effects for measurement temperatures close to T c .

  Our present result obtained within minute-timescale in ultrathin Co gives intermediate η DM I . For the case of the theoretical paper by Yang et al., we have reported two values for the η DM I , considering or not that

	seed	FM	MOx	∆D	∆E	β DM I	η DM I	ref
	nm	nm	nm	mJ/m 2 MV/m fJ/(V•m) 10 -12 J•m/C
	Pt	Co	AlO x /HfO x	0.14-0.26	133 1100-2000	3.9-7.2	here
	3	0.49	6-53	±0.2		±700
	Ta	FeCoB TaO x /AlO x /HfO x	0.1	167	600	3.2	[20]*
	3	0.65	1-10-50
	Ta	FeCoB TaO x /AlO x /HfO x	0.105	667	158	0.61	[20]
	3	0.65	1-10-50
	Au	Fe	MgO/SiO 2	4 × 10 -5	12.5	3.2	1.72	[19]
	50	20	10-270
	Pt	Fe	MgO	0.06	800	75	2.65	[22]*
	4	2	367
	Ta/Pt Co/Pd	MgO/HfO x	0.038	577	67	0.38	[21]
	2.6-2.4 0.78-0.4	2-50
	V/Fe	Co	MgO/SiO 2	1.8 × 10 -3	18	100	5.5	[23]
	20-20	0.14	5-50
	V/Fe	Co	MgO/SiO 2	1.2 × 10 -3	18	65	14.5	[43]
	20-20	0.26	5-50
	Pt	Co	MgO			26	0.6-1.76	[18]
		0.6

TABLE I :

 I

Summary of DMI variations under EF. The long time scale studies (hours or days) have a star * in the last column.

Acknowledgements

The authors thank M. Chshiev and A. Thiaville for fruitful discussions and C. Mahony for his contribution to the editing of the manuscript. This work was supported by the French ANR via Contract No. ELECSPIN ANR-16-CE24-0018), the DARPA TEE program through Grant MIPR No. HR0011831554 and the EUR Grant NanoX No. ANR-17-EURE-0009 in the framework of the Programme des Investissements d'Avenir.