Effects of nano?voids and nano?cracks on the elastic properties of a host medium: XFEM modeling with the level?set function and free surface energy

M. R. Kired · B. E. Hachi · D. Hachi · M. Haboussi

Abstract This work deals with the influences of nano-heterogeneities in the form of voids/cavities or cracks on the elastic properties of a host medium. with a relatively large ratio of apparent surface to volume and particularly strong physical interactions with the surrounding medium at nano-scale, nano-heterogeneities can potentially affect the elastic properties of the parent medium (matrix) containing the m in a significant manner. This has been reported by various theoretical and experimental studies, some of which are discussed in the present paper. To describe the positive (reinforcement) or negative (degradation)effect of the nano-heterogeneities from the modeling perspective, it is necessary to take into account the energy of interfaces/surfaces between nano-heterogeneities and the matrix, which, because of the relatively large extent of their apparent surface and their strong physical interaction with their neighborhood, can no longer be neglected compared to those of the volume energy. Thus, to account for the effects of interfaces/surfaces in a nanostructured heterogeneous medium, the coherent interface model is considered in the present investigation within a periodic homogenization procedure. In this interface/surface model, the displacement vector is assumed to be continuous across the interface while the stress vector is considered to be discontinuous and satisfying the Laplace–Young equations. To solve these equations coupled to the classical mechanical equilibrium problem, a numerical simulation tool is developed in a two-dimensional (2D) context using the extended finite element method and the level-set functions. The developed numerical tool is the n used to carry out a detailed analysis about the effect of nano-heterogeneities on the overall mechanical properties of a medium. The nano-heterogeneities are present in the medium initially as cylindrical cavities (circular in 2D) before being reduced to plane cracks (line in 2D) by successive flattenings.

Keywords Interface/surface energy · XFEM · Level-Set function · Periodic homogenization · Nano-voids/nano-cavities ·

1 Introduction

The fabrication of nanostructured materials such as nanocrystalline metals, thin films, nano-wires, nano-beams,nano-porous materials, and nano-composites wherein nanoparticles are used as reinforcements or doping agents, is now a reality in nanotechnology. As the number of atoms near the surface/interface in these nanostructured materials is relatively large (compared to the total number of atoms),the surface/interface starts to have a specific behavior and affects significantly the properties of the nano-particle or the nanostructured material. This has been confirmed by several the oretical, numerical, and experimental works such as those cited in Refs. [1–17], and the many others cited in the review paper [1]. At the nano-scale, the surface/interface may have its own physics, materialized by a proper constitutive law like the elastic law determined by Miller and Shenoy [2], at the free surface of nano-pores in isotropic aluminum matrix (with shear modulus = 34.7 GPa and Poisson’s ratio = 0.3). Depending on the crystallographic orientation of the polycrystal at the aluminum free surface,these authors identified two sets of surface Lamé’s constants(λs= 3.48912 N·m?1, μs= ? 6.2178 N·m?1for the orientation A l[100]; and λs= 6.842 N·m?1, μs= ? 0.3755 N·m?1for the orientation A l[111]), leading to positive and negative bulk modulii, respectively. The calculations of Miller and Shenoy[2] suggest that the modification of the crystallographic orientation leads the oretically to different free surface behavior.In fact, it is possible to tailor the surface structure in order to get particular behavior by chemical functionalization [3]or controlled in situ crystallization [4]. In this case, it is possible to obtain novel and unusual bulk properties such as nano-porous or nano-cellular materials whose stiffness matches or even exceeds that of the parent materials. This may enable considerable reduction in size and weight of structural elements without sacrificing their strength and other important physical properties. This statement, which has numerically been confirmed in Ref. [5], will also be examined in the present work.

The main objective of this work is the numerical study of the effect of nano-heterogeneities on the effective stiffness of a parent material (aluminum) hosting such nanoheterogeneities. The latter can be either a single (or multiple) nano-void(s)/nano-cavity(ies) or nano-crack(s). In fact,nano-cracks are the ultimate state of flattened nano-voids/nano-cavities in this work. Hence, the initially cylindrical nano-voids/nano-cavities are transformed into increasingly flattened elliptical-cylinder like voids/cavities, before shrinking to the state of plane nano-cracks. As the developments are addressed in two-dimensional (2D) context, the cylindrical void/cavity is represented by a circular domain and the plane crack is represented by a line. The nanoscopic nature of the examined heterogeneities is accounted for by considering the energy contribution of their surfaces, through the use of the Laplace–Young equations of surface equilibrium[18, 19], as it has been done in various publications, some of which are cited in Refs. [5–9].

Regarding the particular problem of cracking effect on effective properties, in the past, several studies have been dedicated to the analysis of the effects of micro- and macrocracks on the effective stiffness of a cracked medium as mentioned in the comprehensive review by Kachanov et al.[20] and in other works [21–23]. In particular, the cracks’distribution and interaction were the main concerns. More recently, as an example of the very few contributions that are devoted to the study of the effective elastic properties of nano-cracked medium, the re is a work [24] that investigates the surface tension effect of nano-cracks on the size-dependent effective electroelastic properties of a piezoelectric, within the framework of non-interaction approximation, under continuously damaged interface formed by nano-cracks.

In accordance with the announced goal, a numerical simulation tool is developed in this work in order to evaluate the elastic properties of a medium containing single (or multiple) nano-void(s), which become crack(s) by successive flattenings. This simulation tool uses an approach combining extended finite element method (XFEM) and the technique of level-set functions, [6, 8, 25–29]. Once validated by comparison with existing numerical [6] and analytical [12]results, the simulation tool is deployed to solve the boundary value problems corresponding to the adopted periodic homogenization procedure [30, 31]. Hence, the effect of nano-voids and /or nano-cracks on the effective stiffness of a medium containing such heterogeneities is investigated.Regarding the surface model considered here, salient results of this study are: (1) identification of a theoretical value of the void radius below which, the porous medium becomes stiffer than the parent matrix; (2) highlighting of the evanescent effect of the surface energy once the void becomes a crack-like.

This paper is organized as follows: in Sect. 2, the theoretical foundations of the present work are given. The n the XFEM discretization of the mathematical problem is presented combined to a level-set in Sect. 3. In Sect. 4, some numerical applications are proposed. The paper is ended with some conclusions and remarks.

Fig. 1 Illustration of the studied heterogeneous nanostructured domain

2 Theoretical formulation of the problem

Consider a linear elastic 2D material domainin static equilibrium, composed of N isotropic sub-domains(Fig. 1). The sub-domain (N) corresponds to the solid phase whereas the first(N?1)sub-domains correspond to inclusions or voids/cavities Ω(i). For the case of voids or cavities,the stiffness tends towards zero. Each sub-domain Ω(i)obeys the following equilibrium field equation:

where ing(i)and b(i)correspond to the bulk Cauchy stress tensor and the body force, respectively. Each inclusion Ω(i)(i=1,2,…,N?1)possesses an interface (with zero thickness) Γ(i)obeying the coherent Laplace–Young constitutive field equation:

w here ing(i)is the surface divergence of the surfaces Cauchy stress at each point x of the interface Γ(i). Both surface quantities, namely the surface divergence and the surface Cauchy stress, are defined in terms of the bulk Cauchy stress = ?ings∶P and ings=P ing P by using the following surface projection operator P:

where n(i)(x)is the outer unit normal vector at each point x of the interface Γ(i).

The field equations (2) and (3) are completed by the boundary conditions on the external boundary ing of the domain Ω:

with ingF∪ ingu= ing and ingF∩ ingF= .

Additionally, the displacement continuity at the different interfaces is required:

For the finite element calculation purpose, the local mathematical problem [Eqs. (1–6)] is reformulated by using the divergence the or em under the following variational form:

In Eq. (7), the unknown is the displacement vector that is assumed to be sufficiently regular and kinematically admissible, i.e. verifying ingu = 0 oningu.

After introducing the volume and interface behavior laws,

and and using again the divergence transformation theorem,Eq. (7) becomes:

In Eqs. (8) and (9), Ω denotes the union of all Ω(i)P?u=usdenote the surface strain tensor ingsand the surface displacement vector us, respectively.

In Eq. (10),F?is the applied force on αingI. The corresponding integral on αingIvanishes if the interface ingIis closed.? and ?sare space-dependent isotropic elastic fourth order tensors of the volume Ω and surface ingIdomains. Precisely,?is null inside a void phase whereas ?schanges from an interface to another according to the relation (see Ref. [6]):

where λsand μsare the surface Lamé’s constants and P is the surface projection operator expressed in Eq. (3).

Also in Eqs. (9) and (10), the stress σ0is such that σ0= τ0P, where τ0represents a residual surface tension at the interface. τ0is taken equal to 0 in the coming applications.

3 Numerical solution procedure: XFEM with the level?set technique

To solve the boundary value problem presented above, an approach using the XFEM combined with the level-set technique is deployed. The domain (Ω) is the n discretized into a set of 3-node triangular elements (P1), whereas no node is specifically assigned to the interface ing(i). This latter is geometrically described by a level-set function φ(i), and is treated from the kinematic point of view by making use of the neighboring nodes, enriched for this purpose.

Each point of the interface ing(i)is parameterized by a proper level-set function (see Ref. [28])

which is equal to zero when the point is on the interface ing(i).When the point is inside Ω(i), the function φ(i)is negative,outside,is positive. Depending on the constants xic, ai, and pi, it is possible to consider different shapes of the interface (circle, ellipse, square, rhombus, etc.).

Each node in the neighborhood of the interface ing(i)is assigned a value of the function. It is the n possible to define a polynomial interpolated function as,where Nj(x)is a 1-degree polynomial shape function. The function(x)will be used in the calculations, instead of the original one φ(i)(x). Accordingly,for each point x= of the interface ing(i), the unit outer normal vector is defined as follows:

For the 3-noded P1 triangular element that is adopted here, the derivativesare uniform functions.

The polynomial approximation of the elementary displacement vector is defined by [25]:

where ψ(x) is the enrichment function defined on the domain(Ω(i)). When Ω(i)corresponds to an inclusion (not a void/cavity), we adopt for ψ(x) the following form,

Initially proposed by Mo?s et al. [27], the enrichment function in Eq. (16) is used here to meet the continuity conditions at the interfaces (continuity of displacements and discontinuities of the traction vector or deformations) by using additional degrees of freedom, without disturbing/falsifying the displacements of the enriched neighboring nodes.

In the presence of a void/cavity instead of an inclusion,another enrichment function, namely V(x), is used, defined as:

For the needs of the different numerical integrations,the interface is geometrically discretized. The position of each discretization point depends on the position of the neighboring nodes and their level function φ(i)values,according to the formula:

Fig. 2 Geometrical discretization of the interface

This is a linear interpolation as shown in Fig. 2. The interface subdivision thus defined enables us to calculate the contour integrals in Eq. (9). It should also be noted that the intersections between the elements and the interface generate triangular sub-elements as shown in Fig. 3. These subelements are also used in numerical integrations by adopting the appropriate number of Gauss integration points.

As announced formerly, the present study that concerns the effect of nano-voids/nano-cavities on the effective behavior of a medium, also covers the case of nano-voids/nanocavities shrinking to nano-cracks. For the comparison purpose, the case of cracked medium is also specifically solved by XFEM, by considering special enrichment for the elements cut by crack. This enrichment is of Heaviside type,H(x). It is used in the interpolated displacement that is of similar form as Eq. (15):

where

Besides the Heaviside function H(x)that is used instead of Ψ(x)in Eq. (20), another level-set function φ(x), not given here, is used in Eq. (21) instead of φ(x)in Eq. (12).

Fig. 3 Splitting of an element crossed/cut by an interface

Taking into account the polynomial approximation[Eq. (15)] of the displacement in the weak formulation[Eq. (10)], leads to the following system of linear algebraic equations:

with in which α = 0 and α = 1 for plane strain and axisymmetric configuration, respectively.

The behavior of domain being linearly elastic and isotropic, the spatially variable stiffness tensor C in Eq. (23) takes the form:

where λ and υ are the Lamé constants. These constants tend towards zero for nano-cavities and nano-cracks.

The surface elasticity tensor in Eq. (24) has an expression that depends of the surface projection operator P defined in Eq. (3), as follows:

Relations (23)–(25) define the bulk stiffness matrix,the surface stiffness matrix and the nodal force vector. In Eq. (22),U is the nodal displacement vector U=u v…u v a b…a bTwhere N=6×N and11NN11MMeM≤6×Newith Nebeing the number of elements.

Equations (22)–(25) are valid for plane strain [γ = 1 in Eqs. (23)–(25)] and axisymmetric problems. For both problems, the adopted expression of strain and stress are:

The matrix Mpused in Eqs. (24) and (25) is given by:

The strain interpolation sub-matrix used in Eqs. (23)–(25)st and s as:

for a non-cut element, and as,

Remarks It is worth noting that the matrix B in Eq. (30) is defined without considering any enrichment. It is suitable for any element of the domain (Ω) that is not cut by the contour of an inclusion and for an element, which is cut by the contour of a void/cavity, as well. In fact, since the inclusion is a void/cavity, the element nodal displacement is not enriched.

When an element is cut by an inclusion contour, the matrix B is enriched by as many columns as additional used degrees of freedom [Eq. (15)]. In Eq. (31), it was assumed that all three nodes of the triangular element were enriched.

Equations (23)–(25) contain volume integrals (on Ω) that are analytically computed for uncut elements over which the matrix B is uniform (shape functions P1: polynomial function of degree 1). In case of cut elements, the volume integrals are numerically evaluated by subdividing cut elements into triangular sub-elements (see Fig. 3). Three Gauss integration points are considered for the numerical integration over each sub-element. For the surface integrals (on Γ(i)),only cut elements are concerned in this calculation. Indeed,the interface (Γ(i)) is discretized by segments ΔΓ(i)defined by φ(i)(x)=0, whose extremities are determined according to the procedure described in Eqs. (18) and (19) and Fig. 2.

Equations above have been implemented and solved numerically for the particular problem of homogenization of a material containing nano-heterogeneities. The concept of representative elementary volume (REV) was used the rein, with different possible arrangements (square, hexagonal, r and om) for the inclusions or voids. Also, different boundary conditions have been considered, such as imposed average strain or periodic displacement. In this paper, only the periodic displacements conditions have been considered.This gave rise to a simulation tool that enables to study the effect of nano-heterogeneities on the effective stiffness of a medium. It was possible to consider a single or multiple nano-heterogeneities that are stiffer or s ofter than the parent material (matrix), or even without stiffness as in the case of voids/cavities. For the latter, we chose to flatten the m until the y shrink to nano-cracks. A ll these calculations have been performed by accounting of interface/free surface effects through a specific energy contribution.

The next section is dedicated to the simulations performed to validate the developed numerical tool and to show its predictive capabilities. Different applications have been addressed, in relation with the effect of cylindrical (circular or elliptical in 2D) nano-voids/nano-cavities and plane (line in 2D) nano-cracks on the effective stiffness of a material.

4 Applications and numerical studies

4.1 Case of hexagonally distributed cylindrical voids in an aluminium matrix

Hexagonal arrangement of voids (see Fig. 4) is considered in an aluminum matrix with surface effects. The properties of the linearly elastic and isotropic aluminum are EM= 70 GPa and νΜ= 0.32.

For the surface properties of the interface, three cases are consid(ered rega)rding the surface bulk modulusdefined as=λs+2υs, as in Refs. [5, 6, 12]. The three cases are:

Fig. 4 Hexagonal arrangement of voids. a Considered REV.b Meshed part

The Lamé’s constants λsand υshave been obtained by Miller and Shenoy thanks to ab initio calculations [2].

The obtained numerical results concerning the normalized in-plane bulk modulus,=(λ +2υ)∕3),are reported in Fig. 5. The y are compared with the results obtained by Y vonnet et al. [6] and Quang and He [12] for the validation purpose of the developed code.

The curves depicted in Fig. 5a, b of the normalized inplane bulk modulus versus void radius and volume fraction,show a good agreement between the present calculations and the ones of Yvonnet et al. [6] and Quang and He [12]obtained numerically and analytically, respectively. In particular, the results show the role of interface/surface energy(≠0), which differs from the case without interface/surface energy (Ks=0), when the void radius is relatively small. Depending on the interface/surface properties, one observes either an increase of the medium effective modulus (Ks>0), or a decrease (<0)compared to the case without interface/surface energy=0). An increase of the radius yields to a decrease of the effect of the surface energy on the effective modulus. Moreover, for a constant radius (R0= 1 nm), an increase of the volume fraction of the nano-void yields to a decrease of the effective modulus in all cases. This decrease, however, is less for>0, and more for<0 than the case=0.

Fig. 5 Normalized in-plane bulk modulus versus. a Void radius. b Void volume fraction

Table 1 Normalized in-plane bulk modulus for different mesh densities (void radius R0 = 1 nm;void volume fraction f0 = 0.2)

Different mesh densities have been used to examine the mesh-dependency of the numerical results in both cases(negative and positive surface bulk moduli). The obtained results reported in Table 1 show that the re is no meshdependency as the computed effective modulus stabilizes after a sufficiently refined meshing.

4.1.1 Effect of the contrast EM/EF

Fig. 6 Normalized in-plane bulk modulus (k? =k ∕max(k M,k F)versus the contrast lg(E M/E F))

Actually, the role of the surface energy defined bystrongly depends on the ratio between the rigidity of the matrix and that of the inclusion or void/cavity. To show this effect, we performed a numerical study concerning the influence of EM/EF, on the medium effective properties, by considering the different surface behaviors Eqs. (32)–(34).

Fig. 7 Normalized in-plane bulk modulus versus the void radius under 0.5 nm for different void volume fractions

The obtained numerical results are summarized in Fig. 6.The y show that once the matrix is more flexible than the inclusion, no influence of the interface/surface behavior,materialized by, is noticed on the calculated in-plane effective modulus. In addition, for a relatively very flexible matrix, the calculated effective modulus tends to zero,which is physically a realistic trend. For a stiffer matrix, as in the presence of void/cavity, we notice an effect of the surface energy through K′, which stabilizes beyond a certain contrast.

After this initial validation, the numerical tool is used in the subsequent sections to simulate different physical configurations. The first one concerns the case when the void radius becomes very small (less than 0.5 nm), in which a very interesting emerging tendency, regarding the effective behavior, has been captured by the simulation tool. The second treated situation concerns multiple randomly distributed voids. For this situation we also consider the case when the nano-voids/nano-cavities flatten until to shrink to nano-crack-like. All these situations are treated by considering the surface energy contribution, for accounting of the nano-metric size and the physical interactions matrix/nanoheterogeneities, at this scale.

4.1.2 An interesting case of very small void radius

We keep the same evaluation procedure as before concerning EM, νMand K′. The effective bulk modulus is calculated for different assumed values of the void radius (R0) in the interval [0, 0.5 nm] and for a void volume fraction (f) between 0.1 and 0.6.

From Fig. 7, we observe that, in the presence of nanovoid/nano-cavity, the effective stiffness of the medium is lower than that of the matrix as long as the radius is greater than a critical ( and purely theoretical) value=0.04 nm.Below that, the overall effective rigidity exceeds that of the matrix. In addition, the increase of the void volume fraction becomes uncommonly favorable to the increase of the stiffness contrary to what happens for R>=0.04 nm.is equal to that of the matrix whatever the volume fraction

4.2 Case of randomly distributed multiple nano?voids with constant volume fraction

We keep the same aluminum matrix, which contains now randomly distributed nano-voids, as shown in Fig. 8.

In order to study the effects of the size of nano-voids on the effective properties of a medium, we consider 30 circular voids randomly distributed. The size of the square medium domain is fixed such that the volume fraction is kept equal to f = 0.3, while varying the radius of the voids. To get statistically representative values of the computed properties,several realizations are performed, up to 40 in order to stabilize the average results as shown in Fig. 9a. The graph shows an increase in the effective stiffness with the decrease of the nano-void radius. At 40 realizations, the results obtained are grouped and compared with those found by Yvonnet et al.[6] in Fig. 9b. We observe a very good agreement between the two results.

4.3 Case of a nano?void flattened towards a nano?crack

We keep the same aluminum matrix, and we consider the case of nano-void, which flattens and gradually tends towards a nano-crack. The questions we answer in the following are: how does the effective stiffness of the medium evolve when the void volume fraction becomes null for the crack configuration? Does the surface energy as it was modeled by Young–Laplace equations, continue to have an influence for nano-cracks?

We define two geometry characteristics for this application, a surface fraction index f=a∕w, and the flattening coefficient c = a/b, (see Fig. 10). The volume fraction of the void is the n equal to π f∕(2c). In the following, the surface fraction index is maintained equal to f′= 0.4, the flattening coefficient is variable, between 0 and 40 (for c = 40 and f′= 0.4, f ~ 1.6%).

Fig. 8 Randomly distributed 30 voids with volume fraction f = 0.3

Fig. 9 a Statistical convergence of the normalized effective in-plane bulk modulus versus the number of realizations for different radius of the nano-void. b Normalized effective in-plane bulk modulus after 40 realizations versus nano-void radius

4.3.1 Without surface energy

We first calculate the normalized effective elastic modulus by considering=0, i.e. without surface energy.The results of the calculations reported in Fig. 11 show an effective elastic modulus that increases very rapidly with the flattening coefficient c, i.e. when the nano-void closes.The elastic property corresponding to a flattened void tends to wards the case of a crack very quickly regarding c (when it is barely greater than 5), with also a very sharp decrease in the void volume fraction. In presence of cracks only, the medium is stiffer than what it is in presence of voids/cavities.However, the cracked medium is still s ofter than the matrix.This proves once more the good predictive capability of the developed code, which also captures the effect on the elastic properties of cracks.

Fig. 10 Flattened nano-void

Fig. 11 Normalized in-plane bulk modulus and void volume fraction versus the flattening coefficient

It should be noted that in order to obtain the value of the effective modulus corresponding to the limit case of flattened nano-voids (almost nano-cracks), it was necessary to use successive refined meshes. The use of Heaviside enrichment in XFEM (in Eq. (20)) during the treatment of the problem with cracks only, makes it possible to obtain the value of the effective elastic modulus without densifying the mesh, as was the case for flattened nano-voids.

The results in Fig. 11 concerning the in-plane bulk modulus versus the flattening coefficient, have been obtained without considering the surface energy, i.e. for K=0. Ins the following, the surface energy will be taken into account.

Fig. 12 Normalized effective elastic stiffness versus the flattening coefficient c for different values of K′s

4.4 With surface energy

The effect of flattening is now studied for the previous conditions (f′= 0.4, EM= 70 GPa, νM= 0.32), by considering the surface energy contributions corresponding to>0 and<0, [Eqs. (32) and (33)]. The obtained variations of the normalized effective elastic stiffness versus the flattening coefficient c are displayed in Fig. 12.

The graph of Fig. 12 shows that the normalized effective elastic stiffness inevitably converges towards the result of a “simple” crack, even with taking into account of a surface energy. However, the convergence is not the same as in the case without surface energy=0). For>0,the increase is faster with c, and even exceeds the effective modulus corresponding to a crack for some flattening ratios less than 20, before decreasing trend towards the value corresponding to a crack, when c is precisely close to 20. In the case<0, the increase is slower, and the effective elastic modulus stays low longer regarding c, before tending to the crack value, later for c around 30. These results show for which flattening coefficient the surface energy, as it was described by Laplace–Young equations, may be important, and when its effect vanishes.

4.4.1 Case of multiple nano?cracks

In this last application, we use the developed code to estimate the effective bulk modulus of a REV under the same conditions as in the previous applications (f′ = 0.4,EM= 70 GPa, νM= 0.32). The REV contains either a single crack or multiple cracks and , in all cases, a fixed total length of cracks, which is equal to that of the single crack(this situation corresponds to a fixed crack density). The cases of cracks oriented successively at 0° and 90°, and randomly with regard to x-axis of the referential frame are considered, as in displayed Figs. 13, 14, 15. It should be mentioned that when the number of crack is relatively high and consequently the cracks are small (as in the case of 20 cracks, the last frame in Figs. 13, 14, 15, it was necessary to densify the meshing (at least three elements per crack) in order to ensure convergence.

Fig. 13 Nano-cracks oriented at θ = 0° parallel to x-axis

Fig. 14 Nano-cracks oriented at θ = 90° perpendicular to x-axis

The results of the various simulations related to multiple nano-cracks are reported in Fig. 16. The y all show that, with respect to the effective modulus, the increase in the number of cracks for the same crack total length is favorable to the medium stiffness. The value of this latter stabilizes for a small number of small nano-cracks (less than 10). This calculation indicates that the damage felt by the material is greater when it is circumscribed to a smaller number of cracks. The effective stiffness corresponding to the cracks oriented at 90° with respect to x-axis, is higher than that obtained in the case of horizontal cracks. This tendency is not due to the loading conditions as the material stiffness should not depend on the loading conditions, but it is purely related to the orientation of the cracks with regards to the frame axes.

By considering randomly oriented cracks, and after carrying out several realizations (5, 10, 20, and 30) to estimate the effective modulus, it is observed that for this configuration, relatively closer to the reality of structures damaged by nano-cracks, the calculation also shows an increase in effective stiffness with the number of cracks, which remains between those of 0° and 90° orientations. In addition, the effective modulus stabilizes at a value close to the most unfavorable case, horizontal cracks, when the number of cracks increases.

Fig. 15 Nano-cracks randomly oriented

Fig. 16 Normalized effective bulk modulus versus the number cracks for different crack orientations

5 Concluding rem arks

In this work, a periodic numerical homogenization tool in 2D has been developed based on XFEM and level-set. This numerical simulation tool enables to study the elastic properties of a medium containing single (or multiple) elliptical void(s)/cavity(ies) with a shape ratio, ranging from 1(cylinder/circle in 2D) to almost 0 (plane crack/line in 2D),in a relatively general context, i.e. without any restriction with respect to interactions, the size of heterogeneities, their numbers, their orientations, and by taking into account the interface/surface energy or force, in order to account for the nanoscopic nature of the heterogeneities.

Among the findings of this work after the performed simulations by using the developed simulation code, some of the salient features are as follows.

(1) The existence, for the adopted interface/free surface Laplace–Young model, of a the oretical size of nanovoid/nano-cavity below which, the porous matrix becomes stiffer than the matrix itself (free from voids/cavities).

(2) The numerical demonstration that by flattening void/cavity, its effect on the effective elasticity of a medium coincides with that of a crack of a comparable size,with the evanescence of the free surface energy contribution, when the void becomes crack-like.

(3) The numerical evidence that by multiplying cracks (for the same total length, i.e. for a fixed crack density), the cracked medium effective stiffness, at least for the bulk modulus, deteriorates less as it increases by multiplying smaller cracks for the same crack total length.

(4) The numerical evidence that horizontal cracks, with regards to the referential frame, are less favorable than vertical cracks for the effective elastic properties, at least for the in-plane bulk modulus.

(5) The numerical demonstration that randomly oriented cracks give an effective modulus lying between those of the horizontal and vertical cracks with respect to the referential frame.

(6) The numerical evidence that by increasing the number of randomly oriented cracks, the elastic properties tend towards the most unfavorable case that of horizontal cracks, at least for the in-plane effective bulk modulus.

In fact, the developed numerical simulation tool is even more general, since it allows simulating other boundary conditions for the homogenization problem, such as the assumed average strain or stress. It also enables to study the effect of the void shape, to consider several heterogeneities in the same time, like inclusions of different shapes that are s ofter or stiffer than the matrix, by also considering the interface/surface effect, without restriction on the interaction between all these heterogeneities.