Komlan Sénam Djaka

Aurélien Villani

Vincent Taupin

Laurent Capolungo

Stéphane Berbenni

Field Dislocation Mechanics for heterogeneous elastic materials: A numerical spectral approach

Keywords: dislocation mechanics, heterogeneous media, elastic fields, spectral method, FFT, numerical algorithms

Spectral methods using Fast Fourier Transform (FFT) algorithms have recently seen a surge in interest in the mechanics of materials community. The present contribution addresses the critical question of determining accurate local mechanical fields using FFT methods without artificial fluctuations arising from materials and defects induced discontinuities. Precisely, the present work introduces a numerical approach based on intrinsic discrete Fourier transforms for the simultaneous treatment of material discontinuities arising from the presence of dislocations and from elastic stiffness heterogeneities. To this end, the elasto-static equations of the field dislocation mechanics theory for periodic heterogeneous materials are numerically solved with FFT in the case of dislocations in proximity of inclusions of varying stiffness. An optimal intrinsic discrete Fourier transform method is sought based on two distinct schemes. A centered finite difference scheme for differential rules are used for

Introduction

The elastic theory of continuously distributed dislocations initiated by Kröner [START_REF] Kröner | Kontinuumstheorie der Versetzungen und Eigenspannungen[END_REF] and others [START_REF] Bilby | Continuous distributions of dislocations: a new application of the methods of non-Riemannian geometry[END_REF][START_REF] Mura | Continuous distribution of moving dislocations[END_REF][START_REF] Willis | Second-order effects of dislocations in anisotropic crystals[END_REF] uses the Nye's dislocation density tensor [START_REF] Nye | Some geometrical relations in dislocated crystals[END_REF] and provides a rigorous basis for the description of the incompatibility of the elastic strain due to the presence of dislocations in crystalline media. This theory was recently revisited by Acharya and co-workers through the so-called Field Dislocation Mechanics (FDM) model [START_REF] Acharya | A model of crystal plasticity based on the theory of continuously distributed dislocations[END_REF]. The finite element formulation of the FDM theory and its phenomenological mesoscopic implementation (called "PMFDM") were presented in [START_REF] Roy | Finite element approximation of Field Dislocation Mechanics[END_REF] and [START_REF] Acharya | Size effects and idealized dislocation microstructure at small scales : Predictions of a phenomenological model of Mesoscopic Field Dislocation Mechanics: Part I[END_REF], respectively. One of the key features of the FDM model in contrast with classic micromechanical methods based on eigenstrain fields [START_REF] Eshelby | The determination of the elastic field of an ellipsoidal inclusion and related problems[END_REF][START_REF] Mura | Micromechanics of defects in solids[END_REF] is the Stokes-Helmholtz decomposition of the elastic distortion into incompatible and compatible parts. The incompatible part is related to the presence of a non-zero dislocation density within the body, while the compatible part is needed to ensure balance of the stress field and boundary conditions. Such a decomposition allows determining a unique solution for the elastic fields associated with a prescribed dislocation density [START_REF] Acharya | A model of crystal plasticity based on the theory of continuously distributed dislocations[END_REF][START_REF] Roy | Finite element approximation of Field Dislocation Mechanics[END_REF].

Due to computational efficiency as compared to the usual finite element method, spectral methods taking advantage of the Fast Fourier Transform algorithm (see e.g. [START_REF] Frigo | An adaptive software architecture for the FFT[END_REF]) were recently developed to numerically solve the static FDM theory [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF][START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF] and to solve the dislocation density transport equation [START_REF] Djaka | A numerical spectral approach to solve the dislocation density transport equation[END_REF]. For single dislocations, the results reported by Berbenni et al. [START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF] were successively compared to analytical solutions and finite element simulations in terms of accuracy of elastic fields near the dislocation core in linear homogeneous and isotropic elasticity. The equations for the incompatible elastic distortions and for the balance of linear momentum are solved in the Fourier space, while the resulting elastic fields are obtained in the real space by using the inverse Fourier transform. An extension of this spectral approach to field dislocation and generalized disclination mechanics was also proposed in [START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF]. To our knowledge, the numerical resolution by FFT of the elasto-static equations of FDM within a general heterogeneous linear elastic framework has not yet been explored so far, except for particular microstructures like infinite bicrystals [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF] for which the analytical piecewise uniform solution fields in each crystal can be obtained, see e.g. [START_REF] Gemperlova | Compatibility stresses in deformed bicrystals[END_REF][START_REF] Richeton | Effects of heterogeneous elasticity coupled to plasticity on stresses and lattice rotations in bicrystals: a Field Dislocation Mechanics viewpoint[END_REF]. For more complex microstructures involving dislocation/inhomogeneity interactions, the difficulty lies in the numerical resolution of the balance of linear momentum (discarding body force and inertial effects), which requires solving an implicit Lippmann-Schwinger integral equation with incompatible elastic fields arising from the presence of polar dislocation densities. Another way is to make use of FFT within a phase field micro-elasticity approach with inhomogeneities solved by an iterative-perturbation scheme [?]. One notes that alternate methods, also relying on the use of FFT, have been proposed to quantify the stress fields associated with the presence of dislocations. For example, the phase-field methods for dislocations [17,18,19,?] or Discrete Dislocation Dynamics methods [START_REF] Bertin | A FFT-based formulation for efficient mechanical fields computation in isotropic and anisotropic periodic discrete dislocation dynamics[END_REF][START_REF] Gao | Influence of misfit stresses on dislocation glide in single crystal superalloys: A three-dimensional discrete dislocation dynamics study[END_REF], where dislocations are described as eigenstrain fields [START_REF] Eshelby | The determination of the elastic field of an ellipsoidal inclusion and related problems[END_REF][START_REF] Mura | Micromechanics of defects in solids[END_REF].

The Lippmann-Schwinger integral equation for periodic heterogeneous elastic or elasto-plastic media was solved through computationally efficient schemes based on the Fast Fourier Transform (FFT) technique to determine the effective and local fields in polycrystals and composite materials [START_REF] Moulinec | A fast numerical method for computing the linear and non linear properties of composites[END_REF][START_REF] Müller | Mathematical vs. experimental stress analysis of inhomogeneities in solids[END_REF][START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF][START_REF] Eyre | A fast numerical scheme for computing the response of composite using grid refinement[END_REF][START_REF] Dreyer | An approximate analytical 2D-solution for the stresses and strains in eigenstrained cubic materials[END_REF][START_REF] Lebensohn | N-site modeling of a 3D viscoplatic polycrystal using fast Fourier transform[END_REF][START_REF] Michel | A computational scheme for linear and non-linear composites with arbitrary phase contrast[END_REF][START_REF] Vinogradov | An accelerated FFT algorithm for thermoelastic and non-linear composites[END_REF][START_REF] Lebensohn | An elasto-viscoplastic formulation based on Fast Fourier Transforms for the prediction of micromechanical fields in polycrystalline materials[END_REF][START_REF] Suquet | Multi-scale modeling of the mechanical behavior of polycrystalline ice under transient creep[END_REF][START_REF] Eisenlohr | A spectral method solution to crystal elasto-viscoplasticity at finite strains[END_REF][START_REF] Anglin | Validation of a numerical method based on fast Fourier transforms for heterogeneous thermoelastic materials by comparison with analytical solutions[END_REF]. Indeed, the FFT numerical scheme allows solving the implicit Lippmann-Schwinger integral equation of the periodic boundary-value problems, by the means of Green's function of a well-chosen reference medium, using different iterative schemes in the Fourier space, while the resulting elastic fields are obtained in the real space by the inverse Fourier transform. Well-known and pioneering iterative schemes are the "basic scheme" and the "accelerated scheme" as reported in [START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF][START_REF] Eyre | A fast numerical scheme for computing the response of composite using grid refinement[END_REF][START_REF] Lebensohn | N-site modeling of a 3D viscoplatic polycrystal using fast Fourier transform[END_REF][START_REF] Vinogradov | An accelerated FFT algorithm for thermoelastic and non-linear composites[END_REF]. The convergence of these two schemes is based on fulfillment of stress equilibrium in the Fourier space [START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF][START_REF] Lebensohn | N-site modeling of a 3D viscoplatic polycrystal using fast Fourier transform[END_REF] and strongly depends on the elasticity moduli contrast and on the choice of the homogeneous reference medium. It was observed that the "basic scheme" and the "accelerated scheme" generally fail to converge when the material contains voids or too rigid inclusions. Therefore, in order to overcome this difficulty, Michel et al. [START_REF] Michel | A computational scheme for linear and non-linear composites with arbitrary phase contrast[END_REF] introduced an iterative method based on "augmented Lagrangians" that is more efficient when the composite exhibits a strong mechanical contrast between phases. More recently, new schemes like the "conjugate gradient" method [START_REF] Zeman | Accelerating a FFT-based solver for numerical homogenization of periodic media by conjugate gradients[END_REF][START_REF] Brisard | FFT-based methods for the mechanics of composites: A general variational framework[END_REF] or the "polarization-based scheme" [START_REF] Monchiet | A polarization-based FFT iterative scheme for computing the effective properties of elastic composites with arbitrary contrast[END_REF] were proposed to solve the linear elastic heterogeneous problem. Some comparisons between the "accelerated", "augmented Lagrangians" and "polarization-based" schemes were provided in [START_REF] Moulinec | Comparison of three accelerated FFT-based schemes for computing the mechanical responses of composite materials[END_REF] with a discussion on the different convergence criteria to use for each scheme. It it was shown that the "accelerated" and the "augmented Lagrangians" schemes are particular cases of the "polarization-based" scheme. The conjugate gradient method was later adapted to non-linear elastic behavior using the Newton-Raphson algorithm [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]. Efficient fixed-point and Newton-Krylov solvers for FFT-based homogenization of elasticity at finite strains (hyperelasticity) were recently developed by Kabel and co-workers [START_REF] Kabel | Efficient fixed point and Newton-Krylov solvers for FFT-based homogenization of elasticity at large deformations[END_REF][START_REF] Kabel | Mixed boundary conditions for FFT-based homogenization at finite strains[END_REF]. Furthermore, the Lippmann-Schwinger formulation for arbitrary mixed boundary conditions solved by FFT-based homogenization at finite strains were also reported [START_REF] Kabel | Mixed boundary conditions for FFT-based homogenization at finite strains[END_REF].

Another important issue inherent to Fourier-based numerical methods is to predict accurate local fields near materials discontinuities like phase/grain boundaries, inclusion corners, or materials defects by preventing from the emergence of spurious oscillations. To control and remove this numerical artifact inherent to spectral methods, different non-classical methods for the computation of partial spatial derivatives of first and second orders have been proposed. For example, "intrinsic" discrete Fourier transforms (DFT) approximation, based on 9-pixel centered finite difference formulas, were used to compute partial derivatives in the Fourier space to solve the elasto-static field equations of linear heterogeneous solids [START_REF] Müller | Mathematical vs. experimental stress analysis of inhomogeneities in solids[END_REF][START_REF] Dreyer | An approximate analytical 2D-solution for the stresses and strains in eigenstrained cubic materials[END_REF]. Later, this method was successfully used to solve the elasto-static fields of FDM in a linear homogeneous isotropic elastic medium with incompatibilities [START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF], in a couple stress medium with generalized disclinations [START_REF] Berbenni | A fast Fourier transformbased approach for generalized disclination mechanics within a couple stress theory, Generalized Continua as Models for Classical and Advanced Materials[END_REF], and, in the elasto-viscoplastic FFT algorithm for polycrystals with non local plasticity [START_REF] Lebensohn | Numerical implementation of non-local polycrystal plasticity using fast Fourier transforms[END_REF] or with a couple stress elasto-viscoplastic formulation [START_REF] Upadhyay | A higher order elasto-viscoplastic model using fast fourier transforms: Effects of lattice curvatures on mechanical response of nanocrystalline metals[END_REF]. It was observed that this "intrinsic" DFT technique allows avoiding spurious oscillations that occur with the use of the classic FFT approximation initially proposed in [START_REF] Moulinec | A fast numerical method for computing the linear and non linear properties of composites[END_REF][START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF], especially when narrow defect cores, possibly assigned on a single pixel, are considered [START_REF] Berbenni | A fast Fourier transformbased approach for generalized disclination mechanics within a couple stress theory, Generalized Continua as Models for Classical and Advanced Materials[END_REF]. Other "intrinsic" discrete Fourier transforms (DFT) schemes based on forward and backward finite differences were first introduced by Willot and Pellegrini [START_REF] Willot | Fast Fourier transform computations and build-up of plastic deformation in 2D, elastic-perfectly plastic, pixelwise disordered porous media, Continuum Models and Discrete Systems[END_REF] for elastic-perfectly plastic porous media, and, modified discrete Green operators for heterogeneous conductivity problems were reported in [START_REF] Willot | Fourier-based schemes with modified Green operator for computing the electrical response of heterogeneous media with accurate local fields[END_REF]. A recent scheme for the computation of partial derivatives, based on centered finite differences on a rotated grid, was proposed by Willot [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF] to compute the modified discrete Green operator. Such a scheme, referred to as "rotated scheme", was shown to be efficient to give accurate local fields, devoid of spurious oscillations when computing the local stress/strain fields in very crude situations, such as cubic inclusions or voids embedded in a matrix phase. Another very recent finite difference discretization scheme performed on a staggered grid was developed by Schneider and co-workers [START_REF] Schneider | Computational homogenization of elasticity on a staggered grid[END_REF]. This refined scheme combined with DFT was also efficiently used for three dimensional porous materials to give solutions devoid of oscillations compared with the solutions obtained from the classic Moulinec-Suquet's discretization. Furthermore, a new discretization method based on linear hexaedral finite elements on Cartesian grids adapted for the basic scheme and the conjugate gradient method was reported [START_REF] Schneider | FFT-based homogenization for microstructures discretized by linear hexahedral elements[END_REF]. Interestingly, an equivalence has been demonstrated between the "rotated scheme" and the use of hexahedral finite elements with reduced integration. The present contribution addresses the critical question of determining accurate local mechanical fields using FFT methods without artificial fluctuations arising from materials and defects induced discontinuities like dislocations interacting with elastic heterogeneities such as inclusions or pores. Precisely, the present work introduces a numerical approach based on intrinsic discrete Fourier transforms for the simultaneous treatment of material discontinuities arising from the presence of dislocations and from elastic stiffness heterogeneities. Motivated by the accuracy and efficiency of the existing FFT-based iterative schemes for heterogeneous elastic materials, we propose in this paper to study the interactions between dislocations and elastic heterogeneities. A special attention will be paid to the control of spurious oscillations inherent to numerical spectral methods regarding local field solutions by comparing different discrete Fourier transforms (DFT) schemes for the approximation of spatial derivatives in the field equations. In particular, we will adopt a DFT scheme, based on 9-pixel centered finite differences for the resolution of the incompatible elastic distortion, and based on the "rotated scheme" for the resolution of the Lippmann-Schwinger equation with polar dislocation densities to determine the compatible elastic strain.

The paper is organized as follows. Notations are introduced in section 2. The elastostatic field equations of continuum dislocation mechanics (static FDM equations) in small strains are introduced in Section 3 as well as the Lippmann-Schwinger integral equation with incompatibilities. For convenience, the complete set of field equations of the FDM theory including the time evolutions of dislocation density and plastic distortion tensors is recalled in Appendix A. In Section 4, the static field equations are solved using the FFT algorithm with prescribed macroscopic uniform stress over the periodic unit cell. More specifically, general three-dimensional solutions for incompatible and compatible elastic fields in the presence of dislocations and elastic heterogeneities are derived to compute strain/stress fields. The FFTbased numerical algorithm is also introduced for the resolution of the Poisson-type equation and of the Lippmann-Schwinger integral equation with incompatibilities.

In Section 4, the different numerical approximations for the treatment of partial spatial derivatives in the Fourier space are introduced using DFT in the case of general three-dimensional problems. Section 5 is dedicated to the comparison of the different approaches to determine the ideal combination of DFT schemes leading to the best accuracy near material discontinuities. Some comparisons will be performed with analytical solutions to simulate the interaction between an edge dislocation and a circular inhomogeneity [START_REF] Dundurs | Interaction between an edge dislocation and circular inclusion[END_REF][START_REF] Dundurs | Edge dislocation inside a circular inclusion[END_REF]. In section 5, the chosen numerical approach is also applied to simulate different physical cases such as dislocation dipoles located near the matrix/inclusion interface or dislocation loop distributions surrounding cubic inclusions in a 3D two-phase composite material.

A bold symbol denotes a tensor or a vector. The symmetric part of tensor A is denoted A sym . Its skew-symmetric part is A skew and its transpose is denoted by A t .

The tensor A•B, with rectangular Cartesian components A ik B kj , results from the dot product of tensors A and B, and A ⊗ B is their tensorial product, with components A ij B kl . The vector A•V, with rectangular Cartesian components A ij V j , results from the dot product of tensor A and vector V. A symbol ":" represents the trace inner product of the two second order tensors A : B = A ij B ji , in rectangular Cartesian components, or the product of a higher order tensor with a second order tensor, e.g., A : B = A ijkl B lk . The div and curl operations for second order tensors are defined row by row, in analogy with the vectorial case. For any base vector e i of the reference frame:

(div A) t • e i = div(A t • e i ), (curl A) t • e i = curl(A t • e i ).
In rectangular Cartesian components, the cross product "×" of a second-order tensor A and a vector V, the div and curl of A are given by:

(A × V) ij = e jkl A ik V l , (divA) i = A ij,j , (curl A) ij = e jkl A il,k ,
where e jkl is a component of the third-order alternating Levi-Civita tensor X and the spatial derivative with respect to a Cartesian coordinate is indicated by a comma followed by the component index. 8

Dislocation density tensor and incompatibility

The theory is developed in the small strain framework (i.e. linear kinematics). Let us assume that the elastic distortion tensor U e is a periodic field within a unit cell V . Let us define a closed circuit C delimiting a surface S of unit normal n so that dislocation lines cross the surface S. In the case where the resulting Burgers vector of all threading dislocation lines is non-zero, there is lattice incompatibility, in the sense that the displacement field becomes multi-valued. A measure of this incompatibility is precisely the net Burgers vector b, which can be measured by integrating the elastic distortion along the circuit C:

b = C U e • dr. (1)
The above equation is a point-wise measure of lattice incompatibility. The continuous density that is associated with incompatible elastic distortion is the Nye's dislocation tensor α, also referred to as "geometrically necessary", "polar" or "excess" dislocation density tensor in the literature. In a Cartesian reference frame (e 1 , e 2 , e 3 ), its component α ij = b i t j provides the net Burgers vector in direction e i per unit surface of S, b i , and the line direction along e j , t j . With this definition, the net Burgers vector can be expressed as:

b = S α • ndS , ( 2
)
where n is the normal to the surface S delimited by the circuit C. Applying now the Stokes's theorem to Eq. ( 1) and comparing with Eq. ( 2) yields the well-known incompatibility equation [START_REF] Kröner | Kontinuumstheorie der Versetzungen und Eigenspannungen[END_REF]:

curlU e = α. (3)
The above equation means that the presence a dislocation density is associated with an incompatible (non-gradient) part of the elastic distortion, whose curl is non-zero.

From Eq. ( 3), there follows immediately a conservation law:

div α = 0. (4)
This equation means that dislocation lines cannot terminate within the body. Here, let us now assume that the dislocation density α is a periodic field over the unit cell V . In this case, the field equations were introduced in [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF][START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF]. The aim of this paper is to determine the elastic fields arising from the presence of dislocation density α coupled with elastic heterogeneities in periodic media.

Poisson-type equation for incompatible elastic distortion

The elastic distortion tensor can be decomposed according to the Stokes-Helmholtz orthogonal decomposition. Indeed, there exists a unique periodic tensor field χ (up to a constant second order tensor) and a unique periodic vector field w (up to a constant vector) such that U e could be written as the sum of the rotational of χ and the gradient of w:

U e = U e,⊥ + U e, = curl χ + grad w, (5)
where, U e,⊥ = curl χ and U e, = grad w represent respectively the incompatible and compatible parts of U e . By applying this decomposition, since curl grad w = 0, Eq. ( 3) yields:

α = curl U e,⊥ , (6)
and since div curl χ = 0:

div U e,⊥ = 0.

Invoking now the identity: curl curl U e,⊥ = grad div U e,⊥ -div grad U e,⊥ , and using Eq. ( 7), the incompatible elastic distortion is solution of the following Poissontype equation:

div grad U e,⊥ = U e,⊥ = -curl α. (8)
In component form, Eq. ( 8) reads

U e,⊥ ij,kk = -e jkl α il,k . (9)
Eq. ( 8) will be transformed in the Fourier space in section 4 and numerically solved using discrete Fourier transforms together with the FFT algorithm in section 5. Note that if the incompatible elastic distortion U e,⊥ is strictly related to a given distribution of the dislocation density α within the periodic unit cell V , the compatible elastic distortion U e, is not involved in this incompatibility problem. As shown below, the compatible part serves to satisfy the balance of the stress field as well as periodic boundary conditions, in the presence of dislocation densities.

Lippmann-Schwinger equation for compatible elastic strain

In the absence of body force and inertia effects, an overall uniform stress T is applied as traction boundary conditions. This corresponds to the spatial average of stress T over the unit cell V using spatial average theorem. For any periodic stress field T, the stress equilibrium equation reads:

div T = 0 in V . (10)
In the following, we adopt a linear elastic constitutive law (Hooke's law) of the form:

T = C : U e,sym = C : ε e , ( 11
)
where C is the fourth order tensor of spatially heterogeneous elastic moduli (with classic symmetries for components:

C ijkl = C jikl = C ijlk = C klij )
and ε e is the elastic strain tensor, i.e. the symmetric part of the elastic distortion tensor U e .

Using Eqs. ( 5) and ( 11), the balance of linear momentum (Eq. ( 10)) can be rewritten in the form of a heterogeneous Navier-type equation in V :

div C : ε e, + f ⊥ = 0, (12)
where U e, = grad w, and the fictitious body force density

f ⊥ = div C : ε e,⊥
reflects the incompatibility arising from the presence of dislocations. Assuming a homogeneous reference medium with linear elastic moduli C 0 , such that C(x) = C 0 + δC(x) (x being a position vector in the unit cell), and using U e, = grad w, Eq. ( 12) yields, in component form:

C 0 ijkl w l,kj + τ ij,j = 0, (13)
where the stress polarization tensor τ is defined as follows:

τ = C : ε e,⊥ + δC : grad w, (14)
where the first term is T ⊥ = C : ε e,⊥ and represents the stress due to incompatible elastic strains obtained from the resolution of the incompatibility problem only, see Eq. ( 8) or Eq. ( 9). Since τ contains the unknown compatible elastic strain ε e, through the term δC : grad w, Eq. ( 13) is to be solved through an integral equation as described now. Introducing the Green's function technique, Eq. ( 13) can be solved in the form of an integral Lippmann-Schwinger equation for the unknown compatible elastic strain ε e, , with the additional presence of an incompatible term due to polar dislocation densities:

ε e, (x) = ε e, -Γ 0 τ (x), (15)
where denotes spatial convolution product and Γ 0 is the modified Green tensor associated with the homogeneous reference elastic medium C 0 . ε e, represents the spatial average of ε e, over the unit cell V . As shown in Appendix B, ε e, can also be written as ε e, = C 0 -1 : Tτ where τ is the averaged stress polarization over the unit cell V . Let us note that Eq. ( 15) is consistent with the Lippmann-Schwinger formulation recently derived by Kabel et al. [START_REF] Kabel | Mixed boundary conditions for FFT-based homogenization at finite strains[END_REF] with mixed boundary conditions in the case where only traction boundary conditions are prescribed. Then, the field solution of Eq. ( 15) is given by series expansion [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF] as:

ε e, (x) = +∞ n=0 -Γ 0 δC (x) n : ε e, -Γ 0 T ⊥ (x) . (16)
The integral equation (Eq. ( 15)) and the traction boundary conditions set an elasticity problem for the unknown field ε e, = (grad w) sym . This equation is similar to the classic integral Lippmann-Schwinger equation for elastic composites. The latter was generally solved through a computationally efficient numerical method based on the Fast Fourier Transform (FFT) using different possible iterative schemes. In our work, we will use in the next section a fixed-point algorithm (the so-called "basic scheme") as in the pioneering work by Moulinec and Suquet [START_REF] Moulinec | A fast numerical method for computing the linear and non linear properties of composites[END_REF] to solve Eq. ( 16)

except for the boundary conditions expressed in terms of tractions. Two novelties are highlighted: first, the Poisson-type equation (Eq. ( 8)) for U e,⊥ is solved by FFT from a prescribed dislocation density tensor α. This makes the numerical scheme different from thermoelastic-like and eigenstrain-like problems as solved by Vinogradov and Milton [START_REF] Vinogradov | An accelerated FFT algorithm for thermoelastic and non-linear composites[END_REF] and Anglin et al. [START_REF] Anglin | Validation of a numerical method based on fast Fourier transforms for heterogeneous thermoelastic materials by comparison with analytical solutions[END_REF], respectively. Second, the series expansion (Eq. ( 16)) contains an extra term which is a convolution product between the modified Green tensor Γ 0 and T ⊥ . Therefore, the numerical FFT algorithm proposed in the next section needs two different procedures which are the calculation of U e,⊥ and the iterative resolution of ε e, , respectively.

4 Spectral approach and numerical algorithm

Fourier transform-based method

The Poisson-type equation (Eqs. ( 8), ( 9)) and the series expansion (Eq. ( 16)) will be solved in the following by using the Fourier Transform method. Indeed, as shown below, the incompatible field U e,⊥ and the compatible one ε e, can be solved in the Fourier space and the obtained elastic fields are finally computed in the real space by using the inverse Fourier Transform.

In the Fourier space, let ξ be the Fourier vector of magnitude ξ = √ ξ • ξ and with components ξ i in Cartesian coordinates. The complex imaginary number is denoted by i and defined as i = √ -1. Let α(ξ), U e,⊥ (ξ), be the Fourier transforms of α(x)

and U e,⊥ (x). Then, the Poisson-type equation is solved using the differentiation theorem in Fourier space. Using component notations, Eq. ( 9) writes in the Fourier space [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF][START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF]:

U e,⊥ ij (ξ) = i ξ k ξ 2 e jkl α il (ξ) ∀ξ = 0, U e,⊥ ij (0) = 0. (17)
Let ε e, (ξ) and Γ 0 (ξ) be, respectively, the continuous Fourier transform of ε e, (x) and Γ 0 (x). The Fourier transform of the integral Lippmann-Schwinger-type equation (Eq. ( 15)) yields:

ε e, (ξ) = -Γ 0 (ξ) : τ (ξ) ∀ξ = 0, ε e, (0) = ε e, . (18)
The Fourier Transform of the modified Green operator Γ 0 can be calculated in Fourier space for anisotropic materials [START_REF] Mura | Micromechanics of defects in solids[END_REF]. In the particular case of isotropic elasticity, the components of Γ 0 (ξ) are given in [START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF]:

Γ 0 ijkl (ξ) = 1 4µ 0 ξ 2 (δ ik ξ j ξ l + δ il ξ j ξ k + δ jk ξ i ξ l + δ jl ξ i ξ k ) - (λ 0 + µ 0 ) µ 0 (λ 0 + 2µ 0 ) ξ i ξ j ξ k ξ l ξ 4 . (19)

FFT algorithm

The direct and the inverse Fourier transforms are computed by using the Fast Fourier Transform (FFT) algorithm. The unit cell V is assumed to have spatial dimensions T 1 , T 2 and T 3 in the x 1 , x 2 and x 3 directions, respectively, and is discretized by a regular rectangular grid with

N 1 × N 2 × N 3 voxels with position vector x = (i 1 δ 1 , i 2 δ 2 , i 3 δ 3 ), where i 1 = 0 → N 1 -1, i 2 = 0 → N 2 -1, i 3 = 0 → N 3 -1 and δ 1 , δ 2 ,
δ 3 are the voxel sizes in the x 1 , x 2 and x 3 directions. The total number of FFT grid points is

N tot = N 1 × N 2 × N 3 and δ = δ 1 = δ 2 = δ 3 . The discrete Fourier transform (DFT) of a given spatial function f is f = FFT (f ). Its inverse Fourier transform is f = FFT -1 f .
The numerical algorithm used to solve the balanced elastic fields associated with any dislocation density distribution within periodic media is reported in the present section. It is constituted of two major procedures: (i) the initialization of the iterative scheme corresponds to the initialization of ε e, for a macroscopic imposed stress T and the computation of the incompatible elastic distortion through the resolution of the Poisson-type equation using Eq. ( 17) in the Fourier space, and, (ii) the global iterative procedure based on the basic scheme (fixed-point algorithm) to solve Eq. ( 16), where ε e, is calculated after convergence is reached.

Hence, the algorithm (see Algorithm 1) works as follows: Once, the dislocation den-sity α(x) is prescribed in the real space, the initialization procedure begins with the computation of α(ξ) in the Fourier space by using direct FFT (step 1). Then, Eq. ( 17) is used to obtain U e,⊥ (ξ) in the Fourier space (step 2), and, U e,⊥ (x) in real space by the use of the inverse FFT (step 3). The initial compatible elastic strain ε e, 0 is taken as the homogeneous elastic solution C 0 -1 : T (step 4). The initialization procedure is concluded by the computation of the initial stress field T 0 (x) set to

C(x) : ε e, 0 + ε e,⊥ (x) (step 5)
, where ε e,⊥ is the symmetric part of U e,⊥ .

In the global iterative loop at iteration (n+1), the stress field known from iteration (n) denoted T n (x), is cast in the Fourier space by direct FFT to obtain T n (ξ) (step 6), and is used to test the convergence criterion based on stress equilibrium in the Fourier space (step 7). If the convergence is reached then the iterative procedure is stopped, else the polarization stress tensor τ n (ξ) is computed in the Fourier space by calculating the FFT of C 0 : ε e, n , and subtracting it from T n (ξ) (step 8). The compatible elastic strain is computed in the Fourier space at step 9 using Eq. ( 18) to obtain ε e, n+1 (ξ) (∀ξ = 0) and ε e, n+1 (0) is taken as ε e, from iteration (n) considering traction boundary conditions (see also Appendix B). Then, the compatible elastic strain is given in the real space from inverse FFT (step 10). The latter is used to update the stress field in the real space (step 11). The iterative loop starts again with the updated compatible elastic strain and stress fields until convergence is reached (see Eq. ( 20)).

Algorithm 1. FFT-based algorithm for elasto-static FDM equations

Initialization: (α(x) being known) 1: α(ξ) ←FFT(α) 2: Solve Poisson equation (Eq.( 17)) to obtain U e,⊥ (ξ) :

U e,⊥ ij (ξ) = i ξ k ξ 2 e jkl α il (ξ) ∀ξ = 0 and U e,⊥ ij (0) = 0 3: U e,⊥ (x) ←FFT -1 ( U e,⊥ ) 4: ε e, 0 ← ε e, 0 = C 0 -1 : T 5: T 0 (x) ← C(x) : ε e, 0 + ε e,⊥ (x) Iterate : n + 1 (ε e, n (x
) and T n (x) being known) 6: T n (ξ) =FFT(T n ) 7: Convergence test based on Eq. ( 20)

8: τ n (ξ) ← T n (ξ) -C 0 : ε e, n (ξ) 9: ε e, n+1 (ξ) ← -Γ 0 (ξ) : τ n (ξ) ∀ξ = 0 and ε e, n+1 (0) ← ε e, n 10: ε e, n+1 (x) ←FFT -1 ( ε e, n+1 (ξ)) 11: T n+1 (x) ← C(x) : ε e, n+1 (x) + ε e,⊥ (x)
The stress equilibrium convergence criterion used at step 7 reads:

e n = div(T n ) 2 | T n | = ξ. T n (ξ) 2 | T n (0)| ≤ , (20)
where . 2 denotes the L 2 norm, |.| is the Euclidian norm of a second order tensor and e n is the error at iteration n. The convergence is reached when the error is smaller than a given precision (typically 10 -6 -10 -8 in our simulations). Then, the elastic strain (incompatible and compatible parts) and stress fields are obtained in the real space.

Different approximations for differentiation rules

Classic approximation

According to [START_REF] Moulinec | A fast numerical method for computing the linear and non linear properties of composites[END_REF][START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF], the classic approximation of the partial derivative of first order is defined as:

iξ j = i 2πm j T j , (21)
where the m j with j = 1 → 3 are arranged in Fourier space as follows:

m j = - N j 2 + 1 , - N j 2 + 2 , ..., -1, 0, 1, ..., N j 2 -1 , N j 2 , (22)
if N j is even, and

m j = - N j -1 2 , ..., -1, 0, 1, ..., N j -1 2 , (23)
if N j is odd.

In this paper, the classic approximation (Eq. ( 21)) will be used to compute partial derivatives of first and second orders. It will be applied to study the elastic fields near material discontinuities. This method will be denoted by P and LS for solving the Poisson-type equation and Lippmann-Schwinger equation with incompatibilities, respectively.

Centered finite difference approximation

The differentiation rules for first and second order partial derivatives onto a discrete grid are based on a 9-pixel approximation of partial derivatives using centered differences. Using this approximation, the corresponding multipliers in the Fourier space between continuous and discrete Fourier transforms for partial derivatives are the following (see [START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF]):

iξ j C = i δ sin 2πm j N j , (24)
-(ξ j ) 2 C = 2 δ 2 cos 2πm j N j -1 , (25)
-ξ j ξ l C = 1 2δ 2 cos 2π m j N j + m l N l -cos 2π m j N j - m l N l , (26)
where, m j and m l are defined with Eq. ( 22) or Eq. ( 23), j = 1 → 3 and l = 1 → 3.

In the following, the accuracy of such centered finite difference approximation will be tested in comparison with the classic approximation (Eq. ( 21)) regarding the numerical resolution of the Poisson-type equation. This centered finite difference method applied to the Poisson-type equation will be denoted by P C . In this abbreviation, the subscript C means "centered finite difference".

4.3.3

Centered finite difference approximation on a rotated grid: "rotated scheme"

The first order partial derivative operator in the Fourier space of the "rotated scheme" is obtained as follows. The details of the method are given in the original paper by Willot [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF]. Here, let us briefly summarize this procedure. First, the displacement field is expressed in the Fourier space at the four corners of a voxel in a 45 • -rotated basis with respect to the original Cartesian basis. The corresponding stress/strain fields are obtained in the center of the voxel in the 45 • -rotated basis using the centered finite difference approximation. Such stress and strain fields are expressed back in the original Cartesian basis, and then cast in the real space through the inverse DFT. Hence, the corresponding multiplier in the Fourier space between continuous and discrete Fourier transform for partial derivatives of first order is given by: iξ

R j = i 4δ tan πm j N j 1 + exp i 2πm 1 N 1 1 + exp i 2πm 2 N 2 1 + exp i 2πm 3 N 3 , (27)
where m j is defined with Eq. ( 22) or Eq. ( 23).

In the following, the accuracy of such centered finite difference approximation on a rotated grid will be tested in comparison with the classic approximation (Eq. ( 21)) to solve the Lippmann-Schwinger equation with incompatibilities. This "rotated scheme" applied to the Lippmann-Schwinger equation will be denoted by LS R .

In this abbreviation, the subscript R means "rotated scheme".

In section 5, the three different approximations will be studied and discussed to solve the Poisson-type equation with P or P C approximations and the integral Lippmann-Schwinger equation with incompatibilities through the basic iterative FFT algorithm decsribed in section 4.2 together with the LS or LS R approximations. First, the comparison of the three combinations, P -LS, P C -LS and P C -LS R will be conducted in the case of a pure edge dislocation density restricted to a single pixel, and interacting with a circular inclusion embedded in a 2D infinite matrix. Second, the P C -LS R method will be applied to physical problems such as dislocation dipoles located near the matrix/inclusion interface, or 3D dislocation loop distributions surrounding cubic inclusions in a two-phase composite material.

5 Numerical results

Dislocation densities and materials parameters

In the forthcoming numerical simulations with 2D and 3D unit cells, different spatial distributions of polar dislocation densities will be assigned within two-phase com-posite materials constituted of inclusions (denoted I) and a matrix phase (denoted M ). For numerical simulations, the Burgers vector magnitude b of the dislocation density is set to 4.05 × 10 -10 m. The dislocations will be located either in the inclusion or outside the inclusion (i.e. within the matrix). Dislocation dipoles and loops constrained at the interfaces between inclusion and matrix will be also studied.

The linear elastic behavior of the matrix phase is isotropic with Young's modulus With such elastic constants, Michel et al. [START_REF] Michel | A computational scheme for linear and non-linear composites with arbitrary phase contrast[END_REF] suggested that an optimal value for the Young's modulus of the reference medium for the iterative basic scheme should be E 0 = β(E M + E I ) with β = 0.5. However, our numerical simulations (see subsection 5.2.2) have demonstrated that choosing a slightly higher value for β such as β = 0.505 improves the convergence of the FFT algorithm reported in section 4.2 with the P C -LS R approximation for k ranging from k = 10 -10 to k = 10 3 . This improvement is also consistent with the convergence study by Willot [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF] who found numerically optimized values of β in the range 0.5003 ≤ β ≤ 0.509 for simulations with the iterative basic scheme together with a pure LS R procedure.

E M =

Internal stresses of edge dislocation interacting with circular inclusion (2D cases)

For these simulations, a large unit cell with dimension 900b × 900b has been discretized using a 2D FFT grid constituted of 2048 × 2048 pixels, with a pixel size δ 0.44b. A circular inclusion I is located in the center of the unit cell with a radius R = 80δ. First, it is assumed that the composite material contains a single edge dislocation represented by a Burgers vector b = be 1 and unit line vector t = e 3 .

The corresponding dislocation density is restricted to a single pixel with magnitude α 13 = b/δ 2 . The dislocation density is defined such that Eq. ( 2) is satisfied. These studied cases correspond to α = 0, which means that the non periodic part of U e can not be obtained [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF]. Different cases such as dislocation density located either in inclusion center or in the matrix phase will be considered. The internal stresses due to mechanical interactions between the circular inclusion and the single edge dislocation density will be computed to compare the accuracy of stress fields generated by the P -LS, P C -LS and P C -LS R procedures with analytical solutions reported in [START_REF] Dundurs | Interaction between an edge dislocation and circular inclusion[END_REF] and [START_REF] Dundurs | Edge dislocation inside a circular inclusion[END_REF]. In order to avoid a null denominator in the stress-based convergence criterion, see Eq. ( 20), a slight non zero macroscopic stress T 11 =1MPa is imposed as traction boundary conditions without any consequence on the reported values for internal stress profiles. An admissible error = 10 -6 was seen to be sufficient to get accurate results in comparison with the analytical solutions. Therefore, a strainbased convergence criterion [START_REF] Kabel | Efficient fixed point and Newton-Krylov solvers for FFT-based homogenization of elasticity at large deformations[END_REF] was not used for the present numerical applications.

5.2.1

Edge dislocation in the center of the inclusion: comparisons between the P -LS, P C -LS and P C -LS R procedures

Here, the dislocation density α 13 is restricted to a single pixel in the center of the circular inclusion I. For this case, the stress profiles T 11 and T 12 normalized by µ M obtained with the numerical iterative FFT scheme are compared to the analytical solutions given by Dundurs and Sendeckyj [START_REF] Dundurs | Edge dislocation inside a circular inclusion[END_REF]. The results for a mechanical contrast k = 0.01 are first reported in Fig. 1. Fig. 1(a) reports the results obtained from the present FFT scheme using the classic P -LS procedure. Spurious oscillations are observed with the P -LS procedure especially for T 11 which demonstrates that this procedure is not able to capture the strong variations of stress field in an accurate fashion at discontinuity interfaces and near the defect core. In contrast, the results obtained from the P C -LS R and the P C -LS procedures reported in Fig. 1(b) do not exhibit oscillations at the same resolution scale and match well the analytical solution.

The results for a mechanical contrast k = 100 are reported in Fig. 2. First, Fig. 2(a) shows the results obtained with the P -LS procedure. Like in Fig. 1(a), spurious oscillations are also observed with the P -LS procedure for T 11 due to strong material discontinuity inside the inclusion and at the matrix/inclusion interfaces.

In contrast, the results obtained from the P C -LS R and the P C -LS procedures reported in Fig. 2(b) dot not exhibit oscillations at the same resolution scale and match the analytical solution well.

From these first results, it appears that the P -LS procedure leads to non accurate stress fields for this problem due to large oscillations near discontinuities reported for the T 11 component. Therefore, this procedure is not efficient and will not be considered in the rest of the paper. In contrast, the P C -LS and P C -LS R procedures appear to give more accurate results even inside the inclusion and near the dislocation core region, as compared to the analytical solutions. However, careful comparisons using a zoom on the T 11 profile in Fig. 3 reveals the occurrence of little oscillations for both contrasts (Fig. 3(a) and Fig. 3(c)) in the case of the P C -LS approximation. The P C -LS R solution is the only one able to reproduce the analytical solution with accuracy and without any oscillation (Fig. 3(b) and Fig. 3(d))

. These observations reinforce and extend the conclusions stressed in [START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF], where the DFT based on 9-pixel centered difference approximation was efficiently used to solve incompatible elastic distortions (P C procedure) in a linear elastic homogeneous isotropic Navier-type equation without any spurious oscillation, especially when the dislocation density is restricted to a single pixel (see figure 18 in [START_REF] Berbenni | A numerical spectral approach for solving elasto-static Field Dislocation and G-Disclination Mechanics[END_REF]). On the other hand, the use of the "rotated scheme" to compute the modified Green operator in the Lippmann-Schwinger equation through the LS R procedure allows to give solutions devoid of any numerical oscillations, which is consistent with the results already reported in [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF] for elastic heterogeneous materials.

As an example of accurate 2D maps of stress fields arising from the interaction between the edge dislocation (located in the center of the inclusion) and the inclusion itself, Fig. 4 shows the components T 11 and T 12 (normalized by µ M ) as a function of and k = 10 3 with a precision = 10 -6 and a slight non zero macroscopic stress T 11 =1MPa. In Fig. 5, it is seen that the P C -LS R is able to converge even for an inclusion akin to a void (i.e. k = 10 -10 ) for a number of iterations (Niter) of 404, which is not the case for the P C -LS procedure for which the convergence study was stopped at k = 10 -3 with Niter=788. Therefore, in addition to a very good accuracy on stress fields in the presence of strong material discontinuities, the convergence of the P C -LS R procedure is faster than that of the P C -LS for k < 1. Interestingly, the P C -LS R can also treat the interactions between dislocations and holes. Now, for an inclusion with k = 10 3 , both P C -LS R and P C -LS procedures give similar number of iterations (in the order of 5000 iterations, i.e. Niter=4729 for P C -LS R and Niter=5710 for P C -LS). This trend for low and high values of k was also observed by [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF] for heterogeneous elastic composites. The advantage of the P C -LS R is the possible treatment with a reasonable number of iterations of the interaction between a hole (then possibly pore or crack distributions) and a dislocation without making use of a more refined iterative scheme than the basic fixed-point algorithm.

As an example, the stresses T 11 and T 12 (normalized by µ M ) obtained with the P C -LS R procedure are plotted in Fig. 6 and Fig. 7 for the edge dislocation interacting with a circular hole (k = 10 -10 ). Both figures highlight the mechanical interaction of the polarized dislocation density with the circular hole, with large stress gradients between the hole (at zero stress) and the dislocation. Furthermore, the T 11 and T 12 stress profiles (normalized by µ M ) obtained near the center of the circular hole are successfully compared to the analytical solutions provided by Dundurs and Mura [START_REF] Dundurs | Interaction between an edge dislocation and circular inclusion[END_REF].

Dislocation dipoles located in the inclusion near the matrix/inclusion interface (2D case)

The present FFT method with the P C -LS R procedure is now applied to a physical case often encountered in crystalline materials, such as a dipole of dislocations located in the inclusion and constrained at the matrix/inclusion interface of a composite material. A 2D dislocation dipole is a plane representation of a dislocation loop constrained at the inclusion/matrix interface in a real three-dimensional configuration. The edge dislocation densities α 13 that compose the dipole are restricted to a single pixel and are assigned with opposite signs. This situation means that α = 0 then a periodic elastic distortion solution can be found [START_REF] Brenner | Numerical implementation of static Field Dislocation Mechanics theory for periodic media[END_REF]. For this physical case, we set k to 0.1. For these calculations, the precision in the convergence criterion is set to = 10 -8 . Fig. 8 shows the resulting stress profiles (blue solid lines).

These profiles are compared to those obtained in the same case without dislocation dipole (pure heterogeneous elasticity problem), as given by the LS R procedure only (red solid lines). reveal that for both loading cases, the presence of the dislocation dipole remarkably modifies the stress field at the inclusion/matrix interface and inside the inclusion, with large stress variations from the center of the inclusion to the interface. Indeed, the solution given by the LS R leads to a uniform stress inside the circular inclusion which is consistent with the classic Eshelby equivalent inclusion method [START_REF] Eshelby | The determination of the elastic field of an ellipsoidal inclusion and related problems[END_REF][START_REF] Mura | Micromechanics of defects in solids[END_REF].

In contrast, the stresses are similar at some distance from the inclusion for both solutions with or without dipoles, which shows a screening effect. The stresses of opposite dislocations cancel each other at long-range distances of single dislocations so that the stresses due to heterogeneous elasticity only are recovered.

Single dislocation loop around a cubic inclusion (3D case)

As a first three-dimensional configuration to apply the P C -LS R procedure, a rectangular shear dislocation loop around a cubic inclusion embedded in the matrix phase is considered. The cubic inclusion is ten times stiffer than the matrix phase

Dislocation loop distributions around cubic inclusions in a two-phase composite material: study of different mechanical contrasts (3D cases)

In these last three-dimensional simulations, the numerical scheme is applied to a more realistic two-phase microstructure with dislocation loop distributions located around cubic-shaped precipitates embedded in a matrix phase. As shown in Fig. 11, the unit cell is constituted of a 3D FFT grid of 128×128×128 = 2097152 voxels with a voxel size of 10 -9 m. The unit cell contains 8 × 8 × 8 = 512 cuboidal inclusions.

The elastic constants of the matrix phase are µ M = 124.2GPa and ν M =0.3. The elastic constants of the precipitates (inclusions) are parametrized as before in this section, i.e. with k the mechanical contrast between both phases and with ν I = ν M .

The size of the cuboidal inclusions is 8 ± 3 × 10 -9 m and the volume fraction of precipitates in the unit cell V is f I = V I /V = 0.15. Each inclusion is then surrounded by 3 ± 1 prismatic dislocation loops with Burgers vector magnitude as set before (b=4.05×10 -10 m) but now directed along the e 3 direction. The dislocation loops are introduced as in the previous 3D simulation (dislocation density on single voxels) to study the efficiency of the P C -LS R procedure in order to accurately compute local fields near dislocations and cuboidal inclusions with sharp corners. A slight non zero macroscopic stress T 33 = 1MPa is imposed such that the obtained elastic fields will be equivalent to the internal fields due to this complex dislocation loop distribution. Comparisons with analytical solutions in the case of an edge dislocation located on a single pixel in the center of a 2D circular inclusion embedded in a matrix phase revealed that the P C -LS R procedure is more accurate and stable than the P C -LS procedure. In addition, it was shown that the P -LS procedure enhances strong oscillations and is not suited for this problem. The efficiency and the accuracy of the P C -LS R procedure was confirmed in a more extreme case such as an edge dislocation density near a 2D circular hole by successful comparisons with analytical solutions. Furthermore, the comparison of the number of iterations between the P C -LS R and the P C -LS procedures, for various mechanical contrasts, has revealed a faster convergence of the P C -LS R procedure for k < 1 (see Fig. 5).

The P C -LS R procedure allows capturing accurately stress fields for physical cases in crystalline materials such as an edge dislocation dipole located in the inclusion and constrained at the matrix/inclusion interface of a 2D composite material. Furthermore, the present approach can be easily extended to 3D simulations, as shown in the case of a cubic inclusion stiffer than the matrix phase and surrounded by a rectangular dislocation loop with edge/screw segments.

Further studies will be focused on the extension of the present numerical methodology for dislocation mechanics to other more refined FFT schemes than the basic fixed-point algorithm like the conjugate gradient scheme. Optimal choices for rigid inclusions will be considered in a future work using the refined discretization procedure using staggered grid [START_REF] Schneider | Computational homogenization of elasticity on a staggered grid[END_REF] together with a dual scheme already developed for heterogeneous media [START_REF] Bhattacharya | A model problem concerning recoverable strains of shape-memory polycrystals[END_REF]. Furthermore, a FFT crystal plasticity spectral framework for the complete elastoviscoplastic Phenomenological Mesoscopic Field Dislocation Mechanics (PMFDM) theory [START_REF] Acharya | Size effects and idealized dislocation microstructure at small scales : Predictions of a phenomenological model of Mesoscopic Field Dislocation Mechanics: Part I[END_REF] is intended to large polycrystalline aggregates by coupling the present numerical approach to the recent FFT resolution of the dislocation density transport equation by Djaka et al. [START_REF] Djaka | A numerical spectral approach to solve the dislocation density transport equation[END_REF]. a spectral approach and starting from the Lippmann-Schwinger equation described in section 3.3. The Cauchy stress tensor T and U e, are obtained when solving the balance of linear momentum without body forces and inertia effects, i.e. div T = 0, together with the Hooke's law T = C : ε e . Similarly, the plastic distortion U p is not a gradient either, and its incompatible part U p,⊥ is defined as the opposite of U e,⊥ . The total plastic distortion rate Up results from the mobility of polarized dislocations, with velocity vector V, constitutively depending on T, α. Completed by boundary conditions on traction and displacement vectors, Eq. (A.1) forms a complete set of partial differential equations with a propagative (i.e. hyperbolic) character for the evolution of dislocation densities α = -curl (α × V). The last equation of Eq. (A.1) was recently numerically solved with a filtered spectral approach by Djaka et al. [START_REF] Djaka | A numerical spectral approach to solve the dislocation density transport equation[END_REF]. The definitions of the div, curl operators and the cross-product in Eq. (A.1) for a rectangular Cartesian basis were introduced in section 2.

B Averaged compatible elastic strain over the unit cell

The stress tensor T is constitutively defined from the Hooke's law (Eq. ( 11)) and using Eq. ( 5) as follows:

62780

  MPa and Poisson ratio ν M = 0.3647. The Young's modulus and Poisson ratio of the inclusion phase are E I and ν I , respectively. These ones will be set up from the definition of the mechanical contrast k = E I /E M = µ I /µ M assuming ν I = ν M with µ I and µ M denoting the shear moduli of inclusion and matrix phases.

Figure 1 .Figure 2 .

 12 Figure 1. Stress profiles T 11 and T 12 normalized by µ M obtained in the case of an edge dislocation density α 13 located on a single pixel in the center of a circular inclusion of radius R = 80δ. For the simulation, a 2D FFT grid of 2048 × 2048 pixels is adopted and µ I = 2.3001 × 10 2 MPa, ν I = 0.3647 with µ I /µ M = 0.01, ν M = ν I . The simulations are performed with a precision = 10 -6 . (a): P -LS solution (green solid lines) vs. analytical solutions AS (dashed lines). (b): P C -LS R solution (blue solid lines) and P C -LS solution (red solid lines) vs. analytical solutions denoted AS (dashed lines). The locations of the inclusion-matrix interfaces are represented by pink dotted lines.

x 1 and x 2 (

 2 normalized by b) obtained with the P C -LS R procedure for three different contrasts k = 0.1 (Figs. 4(a)(b)), k = 10 (Figs. 4(c)(d)) and k = 1 (Figs. 4(e)(f)) (i.e. homogeneous matrix phase, no inclusion). These figures show the influence of the inclusion stiffness (for k = 0.1 and k = 10 compared to the homogeneous case k = 1) on the shape of the stress field due to the dislocation, which reflects different image stresses in the composite due to different elastic moduli contrasts. Furthermore, traction continuity or discontinuity at inclusion/matrix interfaces is well reproduced.

Figure 3 .Figure 4 .

 34 Figure 3. Zoom of stress T 11 normalized by µ M in the inclusion reporting little oscillations with the P C -LS solution (red solid lines) in contrast with the P C -LS R solution (blue solid lines). (a): P C -LS solution with µ I /µ M 0.01, (b): P C -LS R solution with µ I /µ M = 0.01, (c): P C -LS solution with µ I /µ M = 100, and (d): P C -LS R solution with µ I /µ M = 100. The locations of the inclusion-matrix interfaces are represented by pink dotted lines.

Figure 5 .Figure 6 .

 56 Figure 5. Comparisons of the number of iterations needed by the P C -LS and the P C -LS R procedures to converge, according to different mechanical contrasts k in the case of a circular inclusion interacting with an edge dislocation located in the matrix phase of a composite material. The simulations are performed with a precision = 10 -6 and a slight non zero macroscopic stress T 11 =1MPa.

Figure 7 .

 7 Figure 7. Stress profiles T 11 and T 12 (normalized by µ M ) obtained in the case of an edge dislocation density α 13 located on a single pixel at a distance of 120δ to the right from the center of circular hole of radius R = 80δ, embedded in matrix (µ M = 2.3001 × 10 4 MPa, ν M = 0.3647) discretized with a 2D FFT grid of 2048 × 2048 pixels. P C -LS R solution (blue solid lines) vs. analytical solutions denoted AS (dashed lines). The simulations are performed with a precision = 10 -6 and a slight non zero macroscopic stress T 11 =1MPa. The locations of the internal surfaces of the hole are represented by pink dotted lines.

Fig. 8 (

 8 a) shows the T 11 stress component obtained with a non zero macroscopic stress T 11 = 0.01E M and Fig. 8(b) shows the T 12 stress component obtained with non zero macroscopic stresses T 12 = T 21 = 0.01E M . These comparisons

Figure 8 .

 8 Figure 8. Stress profiles T 11 and T 12 (normalized by µ M ) obtained in the case of edge dislocation density dipole located in the inclusion, near the matrix/inclusion interface of a composite material. A circular inclusion of radius R = 80δ with µ I = 2.3001 × 10 3 MPa, ν I = 0.3647 is embedded in a matrix phase with µ I /µ M = 0.1, ν M = ν I . A 2D FFT grid of 2048 × 2048 pixels is considered. For these calculations, the precision is set to = 10 -8 . (a): Stress component T 11 obtained with a non zero macroscopic stress T 11 = 0.01E M . (b): Stress component T 12 obtained with non zero macroscopic stresses T 12 = T 21 = 0.01E M . P C -LS R solution (blue solid lines) vs. pure heterogeneous solution obtained with a LS R procedure only without dislocation dipoles (red solid lines).

(k = 10 ) 2 2 b/δ 2 ,

 1022 with size C = 20δ. The dimension of the unit cell is 900b × 900b × 900b and is discretized using a 3D FFT grid of 128 × 128 × 128 voxels. The dislocation loop is defined by dislocation line and Burgers vectors given by t and b, respectively. The Burgers vector of the loop is in the e 1 direction. The loop is successively constituted of a positive screw segment with line direction along e 1 , a positive edge segment with line direction along e 2 , a negative screw segment with line direction along e 1 and a negative edge segment with line direction along e 2 as shown in Fig.9. Hence, the screw and edge dislocation density components α 11 and α 12 take the value of ±b/δ 2 on each corresponding segment. The corners have both edge and screw densities with a value of ± √ in order to have a uniform norm of the dislocation density tensor everywhere on the loop.The 3D composite material is subjected to non zero macroscopic stresses T 12 = T 21 = 0.01E M . Convergence is reached after 60 iterations for a precision of = 10 -8 .The results obtained with the P C -LS R procedure are shown to be accurate, devoid of any numerical oscillation and the stress components T 11 , T 12 and T 13 (normalized by µ M ) are reported in Fig.10. The stress field generated by the edge segments can be observed through the T 11 component in Fig.10(a) in the (e 1 ,e 3 ) plane. The T 12 and T 13 components are shown in the (e 2 ,e 3 ) plane in Fig. 10(b) and Fig. 10(c), respectively. These ones show the stress fields enhanced by the screw segments.

Figure 9 .Figure 10 .

 910 Figure 9. Schematic illustration of a rectangular dislocation loop surrounding a cubic inclusion with size C = 20δ embedded in a matrix phase of a 3D composite material. The dislocation loop is defined by dislocation line and Burgers vectors given by t and b, respectively. The loop is successively composed of a positive screw segment with line direction along e 1 , a positive edge segment with line direction along e 2 , a negative screw segment with line direction along e 1 and a negative edge segment with line direction along e 2 .

  The loops are made of densities α 31 (dislocation line along the e 1 direction) and α 32 (dislocation line along the e 2 direction). Five contrasts k are performed, namely k = 0.9, k = 0.1, k = 10, k = 0.01 and k = 100. The stress-based criterion for numerical convergence is set to a precision of = 10 -6 and the iteration numbers at convergence are 5, 30, 35, 156 and 283 for the five contrasts, respectively. Figs.11(b),(c), (d), (e), (f) show the predicted internal stress field T 33 normalized by µ M for the five contrasts, respectively. No numerical oscillations are observed although strong stress gradients are observed for the highest contrast k = 100. Due to the polarity of prismaticdislocation loops, the cubic precipitates are in a compression stress state. As it can be observed in these figures, the numerical scheme allows capturing the strong dependence of the internal stress profile on the mechanical contrast for more complex distributions of dislocations with low computational costs in terms of iterations.

Figure 11 .

 11 Figure 11. Calculation of the internal stress field in a matrix containing cuboidal inclusions surrounded by prismatic dislocation loops. (a) the unit cell is constituted of a 3D FFT grid of 128 × 128 × 128 = 2097152 voxels with a voxel size of 10 -9 m containing 8 × 8 × 8 = 512 cuboidal inclusions (green colored). The inclusion volume fraction is 0.15. The inclusions are surrounded by prismatic dislocation loops (gray colored) with Burgers vector along the e 3 direction. The size of the cuboidal inclusions is 8 ± 3 × 10 -9 m and each inclusion is surrounded by 3 ± 1 prismatic dislocation loops with Burgers vector magnitude b=4.05 × 10 -10 m along the e 3 direction. Internal stress field T 33 normalized by µ M for different mechanical contasts k: 0.9 (b), 0.1 (c), 10 (d), 0.01 (e) and 100 (f), respectively.

T

  = C : ε e, + ε e,⊥ . (B.1)Using the definition of the stress polarization tensor τ (see Eq. (14)):T = C 0 : ε e, + τ . (B.2)Therefore, using Eq. (B.2) and traction boundary conditions with macroscopic stress T, the average stress over the unit cell V denoted T yields:T = T = C 0 : ε e, + τ , (B.3)and finally, ε e, can be identified as:ε e, = C 0 -1 : Tτ . (B.4)Then, using Eq. (B.4), the integral Lippmann-Schwinger equation for the unknown compatible elastic strain ε e, (see Eq. (15)) with traction boundary conditions (see section 3.3) also writes:ε e, (x) = C 0 -1 : Tτ -Γ 0 τ (x). (B.5)

Corresponding author.Email address: Stephane.Berbenni@univ-lorraine.fr (StéphaneBerbenni)

Acknowledgements

This work is supported by the French State through the National Research Agency under the program Investment in the future (Labex DAMAS referenced as ANR-11-LABX-0008-01) and the Region Lorraine.

A Field Dislocation Mechanics (FDM) equations

The main Field Dislocation Mechanics equations are presented in the infinitesimal deformation assumption:

The FDM theory uses the continuum description of dislocation densities based on Nye's dislocation tensor α (see section 3.1). The total distortion U is compatible and can be decomposed in a small deformation theory as U e + U p , where U e , U p are the elastic and the plastic distortion tensors, respectively. The elastic distortion U e is not a gradient because of the displacement discontinuity due to the presence of dislocations. Therefore, the elastic distortion can be decomposed into incompatible and compatible parts denoted respectively U e,⊥ and U e, , invoking the Stokes-Helmholtz decomposition (see Eq. ( 5)). The incompatible part U e,⊥ is the solution of the incompatibility equation curl U e,⊥ = α, which can be transformed to a Poisson-type equation (see Eq. ( 8)). The compatible part U e, first introduced in section 3.2 is the gradient of an unknown vector field w, which is solved in the present paper using


Information Systems Security Applications

Protection of the scientific potential and technology of the Nation: Utopia or Reality Finding the right balance : foster innovation, protect knowledge

 The research potential of an laboratory gives a strategic character to the protection of the most sensitive scientific and technical assets.  The security breaches can easily affect scientific or technological data as well as their scientific, technical or human potential.  The laboratory often lives in a complex environment by the diversity of its guardianship and the diversification of its resources, while being confronted with an increasing scientific competition.

Rules for the Protection of the Scientific and Technical Heritage of the Nation, edicted by the Government should be followed, but without hampering research, competitiveness, international exchanges and cooperation: an appropriate balance has to be found.

Laboratories:

  Their equipments, their scientific and technical production, their know-how and their information system (IS) are targets.  Information systems present vulnerabilities connected to diverse threats which can be environmental, intrinsic, human, etc.  Consequences are varied, sometimes dramatic, preventing the pursuit of activities Purpose of protection of the scientific and technical heritage  Classification in protected Sectors and Restrictive Access Zones  Contradiction with the ethics of the public research  Rules of confidentiality, thesis defense, etc.  Impossibility of respecting submission deadlines if authorization is required  Hardening of access permissions (PhD and postdoctoral students, visitors, foreign researchers, etc.)  Additional work due to new management procedures  Loss of PhD students due to Ministerial delays  Members are implied in the valuation of their work, but precautions are necessary  The raising awareness and training of all the researchers for a good level of safety and confidence  Equipment, software, etc.  Databases, strategic data  Industrial and government contracts  Know-how of researchers, etc  Patents pending, current publications  Current and future projects  Research and experience  Reputation and image Scientific Heritage External or internal  Passive or Active  Physical, electronic  Organizational, human  Industrial espionage  Cyber spying and hacking  Interference with persons or property  Fraudulent use of computer resources  Legal risks Issues A failure of the security of the information system can affect all activities and assets of the laboratory Loss of  Credibility  Degradation of the image  Trust of partners  Market share  Confidentiality, Availability, Integrity of IS  Disappearance of scientific advantages,  Recognition of peers  Reputation of the laboratory Security criteria to minimize the impact  IS Confidentiality, IS availability, IS integrity  Authenticity of identities (human, machine, etc.)  Protection of sensitive data  Legal protection  Ownership  Responsibility  Use Threats But the everyday world is quite different, Constraints lived as an infringement on the freedom of research Ideally Researchers must have a strategy of protection of both their knowledge and their know-how Striking a balance between the spreading of knowledge, exchanges and international cooperations on one side, and the security constraints on the other  The adopted measures should not hamper the research and the competitiveness  Simple and common sense solutions for the whole laboratory  Communication, training and demonstrations to find the best way to limit constraints  Take into account professional needs in order to find necessary measures and not the reverse  But don't forget simple security and pragmatic measures URSI-France Workshop Radio Science for Humanity 2017, Feb 2017, Sophia Antipolis, France Jean-Pierre DAMIANO, UCA CNRS LEAT