of composites interphases: Its application in analysis of effective properties Lurie solution for spherical particle and spring layer model of Mechanics of Materials

ContentslistsavailableatScienceDirect

Mechanics of Materials

journalhomepage:www.elsevier.com/locate/mechmat

Lurie solution for spherical particle and spring layer model of interphases: Its application in analysis of effective properties of composites

Lidiia Nazarenko

^a,∗

, Swantje Bargmann

^b,c

, Henryk Stolarski

^d

aInstitute of Mechanics, Otto von Guericke University Magdeburg, Germany

bInstitute of Continuum Mechanics and Material Mechanics, Hamburg University of Technology, Germany

cInstitute of Materials Research, Helmholtz–Zentrum Geesthacht, Germany

dDepartment of Civil, Environmental and Geo- Engineering, University of Minnesota, USA

a r t i c l e i n f o

Article history:

Received 9 June 2015 Revised 8 January 2016 Available online 8 February 2016 Keywords:

Spherical equivalent inhomogeneity Spring layer model

Interphase Random composites Effective properties

a b s t r a c t

Anewapproachtothedeterminationofequivalentinhomogeneityforsphericalparticles andthespringlayermodeloftheirinterphaseswiththematrixmaterialisdeveloped.To validatethisapproachtheeffectivepropertiesofrandomcompositescontainingspherical inhomogeneitiessurroundedbyaninterphasematerialofconstantthicknessareevaluated.

Thepropertiesofequivalentinhomogeneity,incorporating onlypropertiesoftheoriginal inhomogeneityand itsinterphase, aredeterminedemployinganewapproach basedon theexactLurie’ssolutionforspheres.Thisconstitutesthecentralaspectoftheproposed approachbeing incontrast with some existingdefinitions of equivalentinhomogeneity whose properties dependent alsoonthe properties of the matrix.With the equivalent inhomogeneityspecified as proposedhere,the effectiveproperties ofthematerialwith interphasescan befoundusinganymethodapplicabletoanalysisofthematerialswith perfectinterfaces(i.e.,withoutinterphases)andanypropertiesofthematrix.Inthiswork, themethodofconditionalmomentsisemployedtothisend.Thechoiceofthatmethodis motivatedbythemethod’ssolidformalfoundations,itspotentialapplicabilitytoinhomo- geneitiesotherthanspheresand toanisotropic materials.Theresultingeffectiveproper- tiesofmaterialswithrandomlydistributedsphericalparticlesarepresentedintheclosed- formandareinexcellentagreementwithvaluesreportedintechnicalliterature,whichare basedonbothformallyexactandapproximatemethods.

© 2016TheAuthors.PublishedbyElsevierLtd.

ThisisanopenaccessarticleundertheCCBY-NC-NDlicense (http://creativecommons.org/licenses/by-nc-nd/4.0/).

1. Introduction

Connectionbetweendissimilar materials is always ac- companied by the presence of a layer, called interphase, whosepropertiesare differentthan thoseofthe adjacent

∗ Corresponding author. Tel.: +49 3916751843; fax: +49 3916712863.

E-mail address: lnazarenko@yandex.ru , lidiia.nazarenko@ovgu.de (L.

Nazarenko).

bulk materials. Some problems involving such a layer(or many of them), e.g., single spherical inhomogeneity sur- rounded by layers of a different material, can be solved analytically usingaformallyexact approach(Lurie,2005).

However, most problems of that kind are too difficult (or too demanding) to lend themselves to exact analyt- ical treatment. Consequently, over the past 40 years or so (Benveniste1985; Hashin1961, 1990,1991;Liptonand Vernescu 1995;Luo andWeng1987; MalandBose 1974;

Walpole 1978;Zhongetal.1997)variousapproacheshave http://dx.doi.org/10.1016/j.mechmat.2016.01.011

0167-6636/© 2016 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license ( http://creativecommons.org/licenses/by-nc-nd/4.0/ ).

been proposed to approximately capture theinfluence of interphaselayersonthe quantitiesofinterest,suchaslo- calstress andstrainfieldsoreffectivepropertiesofmulti- constituentmaterials(composites).

Theinterphasesarewheresomeofthemostimportant, complicated,andinterestingphenomenaincompositema- terials oftenoccur. Incomposites witha highinterphase- to-volume ratio, this may strongly affect the overall be- haviorofthematerialandrequiresappropriate modelsto capturesuchbehavior analytically.Insome situations,the interphase is well described by a mathematical surface, called interface. If both the displacements and tractions (stress vectors) can be assumed to be continuous across the interface, it is commonly called a “perfect interface”.

There exist situations, however, when it is more appro- priate to usean “imperfect interface” model, i.e., an inter- faceacrosswhichdisplacements,ortractions,orboth,suf- fer suitablydefinedjumps (cf.GurtinandMurdoch1975;

Gurtin et al. 1998; Benveniste and Miloh 2001; Hashin 1991,2002a;GuandHe2011;Guetal.2014).

In general, the interphase can be modeled by an in- terface(either perfectorimperfect)iftheratioofthe in- terphase’s thickness to a characteristic dimension of the composite (typically the size of the embedded inhomo- geneities) is “sufficiently” small (Benveniste and Miloh 2001; Hashin 1991, 2002a; Gu andHe 2011; Dong et al.

2014). The word “sufficiently” is typically understood in the asymptotic sense meaning that, treating this ratio as theonlysmallparameterinthecomposite,theinterphase can be adequately described by the leading term of the asymptoticexpansionoftherelevantequations.Depending onhowotherdataoftheproblemarerelatedtothatratio variousinterfaceconditionscanbederived(Benvenisteand Miloh2001;seealsoHashin2002a;RubinandBenveniste 2004;Dongetal.2014).

Presence of interphases, or imperfect interfaces, sig- nificantly complicates evaluationof those properties, and the literature in that area is overwhelmingly numerical (Andrianovetal.2007;AchenbachandZhu1989;Sangani and Mo 1997; Garboczi and Berryman 2001; McBride et al.2012,amongothers).Evenapproachesbasedonsimpli- fiedmethods,suchasself-consistentschemes,oftendonot leadtoclosed-formsolutionsandrequirenumericalcalcu- lations,( Hashin1991,2002a).Treatment ofperfectinter- faces, althoughcomplex too,ismoreadvanced,ascan be gleanedfromseveralgoodbooks(Christensen1991;Mura 1987;Nemat-NasserandHori,1999;Torquato2002among others). Effective(or homogenized)propertiesofcompos- iteswithperfectinterfacesareoftengivenby closed-form formulas, which is a very attractive feature. Due to the complexityoftheproblem,thoseformulasaretypicallyap- proximate;only forcomposites withregular arrangement of particles or fibers formally exact (typically numerical) resultshavebeenobtained(Andrianovetal.2007;Sangani and Mo 1997;Garboczi andBerryman 2001,among oth- ers).

Aninterfacemodelthathaslongbeenusedinanalysis of composites is the so called spring layer model. In that model tractions are continuous while displacements are allowed to experience a jump across the interface (Duan etal.2007a,b;Hashin1990,1991;Sangani andMo1997).

Thus, the continuum interphaseis replaced by a layer of normalandtangentialsprings.Thisisanadequateapprox- imationofrealityiftheinterphaseisthinandifit issuf- ficientlycompliantincomparisonwiththestiffnessofthe surroundingmaterials. Aslightlymodifiedversion ofthat modelisalsoemployed inthiswork,butitisusedinthe context ofan entirelynew approach toevaluate effective propertiesofcompositeswithsphericalinhomogeneities.

Adifficultythat isassociatedwithallinterphasemod- els, including the spring layer model, is determination of the parameters needed to describe them. Analyses of Hashin(2002a)andBenvenisteandMiloh(2001)providea rationalebehindhowtheparametersdescribingvariousin- terfacemodelsshouldberelatedtopropertiesoftheinter- phasetreatedaselasticcontinuum.Whiletheoreticallyin- sightful,thisrationaleisoflimitedpracticalutilityasdirect experimentaldeterminationofanyinterphaseparameters, whetherrelatedtoitscontinuumdescriptionorotherwise, isimpossible.Thatcanonlybedone viaaninverseanaly- sis,Linetal.(2005);Wangetal.(2008),asdonebyHashin andMonteiro(2002b). Theresults ofthat work will sub- sequentlybeused forvalidationofthe methodologypro- posedherein.

A rare, quite complicated but “formally exact” three- dimensional numerical solution (employing series expan- sionintermsofsphericalharmonics)foreffectiveproper- tiesofacompositematerial containingspherical particles withspringlayer interfaceswasobtainedby Sanganiand Mo(1997).They tabulatedthe results forvarious volume fractionsandvariouspropertiesoftheparticles,aswellas forvariouspropertiesofthespringlayer.Theseresultsare invaluableforcomparisonsandwillbe usedforthatpur- posehere.

In case of composites with interphases or imperfect interfaces and inhomogeneities in the form of particles or fibers, one viable and attractive way of evaluating theireffectivepropertiesistoreplacethe inhomogeneity- interphasesystemby an “equivalentinhomogeneity” with suitablyadjustedpropertiesthatrepresentsboththeorigi- nalinhomogeneity andthe interphaseorimperfect inter- face. In that way, the problem of the material with im- perfect interfaces is replaced by the problem with per- fectinterface,butwithchangedpropertiesoftheinhomo- geneities. Consequently, using equivalent inhomogeneity, all theexisting closed-formresultsfor two-phasemateri- alscanbeutilizedtoobtainthepropertiesofthecompos- iteswithimperfectinterfacesorinterphases.Thisapproach has been pursued in several prior publications (Duan et al., 2007a; Hashin2002a; Shen andLi2005; Sevostianov andKachanov 2007; Nazarenkoetal.2015a,b) (described inthesequel)andanewversionofitispresentedherein.

Thenewversionisanalternativetotheenergy-equivalent inhomogeneity, originally presented by Nazarenko et al., (2015) in the context of the Gurtin-Murdoch (1975) in- terphasemodel.Inthisworkthe“spring-layer” interphase modelisconsideredinstead(cf.Hashin1990,1991;Achen- bach and Zhu 1989; Sangani and Mo 1997; Duan et al 2007a,b;Hashin2002a;Guetal.2014amongothers).The resulting properties of the equivalent inhomogeneity are used in conjunction with some previously obtained re- sultsfor two-phase composites with perfectinterfaces to

determine the effective properties of random composites withspring-layerinterfaces.

Totheauthors’bestknowledgethepossibilityofusing theequivalentpropertiesoftheinhomogeneity-interphase (or imperfect interface) system is first mentioned by Hashin, (1991). The authordiscusses an extension ofthis techniqueto problemsinvolving multiplelayers of differ- entinterphasesHashin(2002a).However,Hashin’sdiscus- sionwasrestrictedtothe effectivebulkmodulus, asonly forthat casehewasable to develop closed-formexpres- sion; a numerical approach wasnecessary forevaluating theeffectiveshearmodulus.

Prediction of effective moduli of multiphase compos- itesbasedonthenotionofequivalentinhomogeneitycom- binedwiththegeneralizedself-consistentschemewasalso presentedin the two-partpaper ofDuan et al.(2007a,b) andGuetal.(2014).Sphericalparticlesorcylindricalfibers withvariousinterface effects orinterphaseswere consid- ered,includingthespringlayerinterface modeldealtwith herein. In both of the above contributions the unknown stiffnesstensoroftheequivalentinhomogeneitywasdeter- minedusingtheEshelby’sformulaforthechangeofelastic energycausedbyinsertionofaninhomogeneityinaninfi- nitematrix(Eshelby,1957).Theyassumethatsuchchange dueto embedding the equivalentinhomogeneity is equal tothechange causedbyembedding theoriginal inhomo- geneity together with its interface. This approach yields propertiesofthe equivalentparticles/fibersdependent on thepropertiesofthematrixmaterialwhichisnon-physical and it is in sharp contrast with the approach advocated byHashin(2002a)andwiththeapproachproposedinthe presentwork.Thisisthemainreasonbehindthenewfor- mulationoftheproblempresentedinthiswork.Detailsof theprocedure proposed hereby whichthe original inho- mogeneityanditsinterphaseisreplacedby anequivalent inhomogeneity,aswell as comparisonwiththe pertinent resultsobtained in thepast, are presented inthe second andthirdsectionsofthepresentarticle.

Insummary,itis notedthat theexisting solutionsex- plicitlyaccountingforpresence ofthespringlayer model ofinterphasewere obtainedin thetrue closed-formonly fortheeffectivebulkmodulus (Hashin1990,1991,2002a;

Duanetal.2007a,Guetal.2014).Other effectiveproper- ties(shear moduli)were given by a sequence of compli- catedformulas that still had to be evaluated numerically (Hashin 1990, 1991, 2002a; Duan et al. 2007a, Gu et al.

2014)orwereobtainedbynumericalmethods(Andrianov etal.2007;AchenbachandZhu1989;GarbocziandBerry- man2001;SanganiandMo1997).Thedefinitionofequiv- alent inhomogeneity proposed here is significantly sim- plerand should enhance predictive capabilitiesofall ap- proximate methods, irrespective of the material property of interest (bulk and shear). It is more direct than all ofthe severalexisting definitions, such asthat based on the Mori–Tanaka method (Shen and Li 2005), Hashin–

Shtrikman bounds (Sevostianov and Kachanov 2007) or Eshelby(1957)formulaofthe energychange(Duanetal.

2007a;Guetal.2014).Useofsomeofthoseapproachesis likelytointroduceasignificanterroralreadyatthelevelof determiningthepropertiesoftheequivalentinhomogene- ity. If such an equivalent inhomogeneity is subsequently

usedtodeterminetheeffectivepropertiesofthecomposite using additional approximate schemes, an essential error canresultforcompositeswithhighcontrastinthecompo- nentpropertiesandforhighvolumefractionofparticles.

The equivalent properties of the inhomogeneity- interphase systemin this work are determined based on the Lurie’s solution (Lurie 2005) for spheres. In contrast with that of Duan et al. (2007a) and Gu et al. (2014), they depend only on the properties of the original in- homogeneity and of the spring layer parameters (not on the properties of the matrix). Other than being non- physical, dependence of

μ

^{eqdeveloped by Duan et al.}

(2007a) and Gu et al. (2014) on the properties of the matrix seems to make it applicable only in combination with self-consistent method (which is employed in these papers). In contrast,

μ

^{eqpresented in our work, which is}

independentofthepropertiesofthematrix,maybeused in conjunction with any method (old or new) developed to evaluatetheeffectivepropertiesofcompositeswithout interphases, and with any matrix material. One small featureofthemethoddiscussedhereinisthat,asopposed to the previous springlayermodels whose thicknesswas vanishingly small, a finite thickness of the spring layer is still retained.That feature is thereason for whichthe term “spring layer interphase” is used throughout this work. Finite interphase thickness has been included in other models,such asCosserat model ofRubin and Ben- veniste (2004) andDonget al.(2014),butto the authors knowledge, notinthespringlayer models.Still,themain focusofthisarticleisonthenewdefinitionofequivalent inhomogeneity and its validation through comparison of the effective propertiesbased on that definition andthe best analytical and experimental results available in the literature.

The basic assumptions and formulas behind the def- inition of the equivalent homogeneous representation of the inhomogeneity-interphase system are presented in Section 2. In Section 3 the probability-based technique, calledthe methodofconditional moments,isbriefly out- linedand– intandemwiththeresultsobtainedinSection 2 – applied to evaluate effective properties of random composites withspring layer interphases.To validate the proposed approach some representative results are pre- sented and compared to selected existing developments.

Section 4 containsan overall discussion of the approach, presents someconclusionsandevaluates potential forfu- tureextensionsandapplicationsofthe approach.Thepa- peralsoincludesanappendixtowhichseveraldetailsper- tinent tothe developmentincludedinSection 2are rele- gated.

2. Thepropertiesoftheequivalentinhomogeneity 2.1. Ageneralviewoftheapproach

The equivalent inhomogeneity is devised to represent the system consisting of the original inhomogeneity and surroundinginterphaseofthicknessh.Theinterphaseand theinhomogeneityareassumedelasticandhavetheirown distinct properties.In the subsequent developments,it is assumedthat theinterphaseis adequatelyrepresentedby

normal and tangentially oriented linear springs (spring layer).

The concept of the energy-equivalent inhomogeneity is essentially equivalent to a two-stage homogenization.

In the first stage, using energy equivalence, the individ- ual inhomogeneity and its interphase (spring layer) is replaced by an effective (“homogenized”) inhomogeneity which combines the properties of both. In the second stage, the effective inhomogeneity is perfectly interfaced withthematrix(nointerfacialjumpsofanykind)andef- fectivepropertiesofthecompositematerialevaluated.Any homogenizationapproach(numericaloranalytical)canbe usedinthesecondstage.

Thefirst stageoftheprocess issimplyalow-levelho- mogenization step (or sub-homogenization), and like any homogenizationprocedure usedindetermining theeffec- tive propertiesofasystem, requiresthedisplacementson the boundaryofthe system(consisting ofmatrixandin- terphase, in this case) to be consistent with an average strain tensor of equivalent inhomogeneity ɛeq. At equi- librium these displacements cause attendant strain fields withintheoriginalinhomogeneityanddisplacementjumps acrosstheinterphase,bothofwhichdependontheequiv- alentstrainsɛeq.Tosolveforthosequantitiesitisassumed that the original inhomogeneity undergoes a deformation described by a strain tensor ɛ associated with the Lurie solutionforspheres (Lurie2005). Asa result,some stress components(tractions)attheboundarybetweentheinter- phaseandthematrixcanbeevaluatedasafunctionofthe properties of inhomogeneity, the properties of the inter- phaseandthestrainɛeq.Theycanbeusedtoevaluatethe average stresses in the inhomogeneity-interphase system which leadsto relations foreffective propertiesofequiv- alent inhomogeneity Ceq. It is emphasized, however, that iftheoverall(equivalent)propertiesofthecompositema- terial are sought, displacements at the matrix/interphase boundary are not necessarily constant – they depend on thehomogenizationtechniqueusedforthispurpose.

2.2. Spring-layermodeloftheinterphase

Thespring-layermodeloftheinterfacehaswidelybeen used to describe the so called “soft interphase” (Hashin 1990, 1991; Achenbach and Zhu 1989; Sangani and Mo 1997; Duan et al. 2007a,b). In this model interface trac- tions remain continuous across the interface while dis- placements suffer a jump. With respect to a curvilinear orthogonal coordinate systemon thesphere’s surfacethe mathematicaldescriptionofthosepropertiesis

σ

^dSS·n=[

σ

^2dS2−

σ

^{1dS1]·n=0, K·}u_S=

σ

^1·n.

(2.1) The vector n is unit and normal to the interface be- tween inhomogeneity andmatrix. It is assumed that the normalnpointsawayfromtheinhomogeneity.Thedouble squarebracketsindicatethejumpoffieldquantitiesacross theinterface,andthesuperscripts1and2indicatethatthe appropriate quantities are evaluated on the inhomogene- ityormatrixside ofthespringlayer. K=Knn⊗n+Kss⊗ s+Ktt⊗t isa second-order tensor. Kn,Ks andKt arethe springlayerstiffnessparametersinnormalandtangential

Fig. 1. Schematic illustration of inhomogeneity with interphase.

directions,respectively,andsandtrepresenttwoorthog- onalunitvectorsintheplanetangenttotheinterface.

The surfaceelement dS introduced in Eq. (2.1)₁ is re- latedtothefactthatthefinitethicknessofthespringlayer isretainedinthepresentdevelopment.Givenaninfinites- imalareaontheinhomogeneitysideofthelayer,bounded byacontourГ,thearea onthematrixside isoutlined by the curve formed by intersections of the lines normal to theinhomogeneity along Гwiththematrix.Such modifi- cation ofEq.(2.1) guarantees that the equilibriumof the springlayeroffinitethicknessismaintained.

By Hashin(2002a), it wasshown that thin andcom- pliant interphase can accurately be modeled by a spring layer.The sameconclusion wasalsodrawnby Benveniste andMiloh(2001),andinbothcasesitwasestablishedthat

Kn=

λ

i+2

μ

i

h , Kt=Ks=

μ

i

h , (2.2)

where

λ

iand

μ

iareLameconstantsoftheinterphase.This idea isadopted inthe presentmanuscriptandproperties of inhomogeneity/interphase system shown in Fig. 1 are determinedusinginterfacemodeldescribedbyEq.(2.1).

In this case, the displacement jump is considered as the difference between the displacement on the inter- phase/matrixsurfaceu₂ andtheone ontheinhomogene- ity/interfacesurfaceu₁

u=

^{u=u2−u1, (2.3a)}

where

u2=

ε

^{eq·[r1+hn]. (2.3b)}

r₁ is theradius oftheinhomogeneity,andhis thethick- ness of the interphase. The displacement vector u₁ will be specified on the basis of Lurie solution for spheres in Sections 2.3 and 2.4. In view of the fact that the in- homogeneity is isotropic and the spring layer properties are constant (with Kt=Ks), the properties of the equiv- alent inhomogeneity will also be isotropic. Consequently, thebulkmodulusandtheshearmoduluscanbeevaluated separately, which makes the analysis simpler and more efficient.

2.3.Modulifortheequivalentinhomogeneity

2.3.1. Bulkmodulus

Inordertoevaluatetheeffectivebulkmoduluswecon- siderhydrostaticaveragedeformationfortheinhomogene- ity/springlayersystem

ε

^eq=

β

0 0

β

0 0

0

β

0

, ^(2.4)

with corresponding displacements on the inter- phase/matrix surface. In the spherical coordinate system (r,

θ

^,

ϕ

^{),andforaspherical}inhomogeneityofradiusr₁=R, thisimpliesthat

u₂=

ε

^eq·r

|

^r=R+h ⇒ u_2[r]=

β

^[R+h],

u2[θ^]=u2[ϕ^]=0. ^(2.5)

Due to complete rotationalsymmetry ofthe problem, thedeformationwithintheoriginalinhomogeneityisu_r= u_r(^r)^{, u}_θ=0 andu_ϕ=0. Consequently, the displacement fieldvanishingatitscentercanbewrittenas

ur=Fr

R, u_θ=0, u_ϕ=0. (2.6)

TheconstantFwillbedeterminedconsideringEqs(2.1), (2.2),(2.3),(2.5)andnotingthat

u_1[r]

|

r=R =F. ^(2.7)

This,together withEq. (2.5) permitsto determine the displacementjumpsinEq.(2.1)

^ur=u_2[r]−u_1[r]=

β

^[R+h]−F,

^u_θ=

^uϕ=0. (2.8) Inordertodeterminestresseswithintheoriginalinho- mogeneity,alsoenteringthisequation,itisnotedthat

ε

rr=

∂

^ur

∂

^{r =FR,}

ε

θθ =1 r

∂

^uθ

∂θ

^{+urr= FR,}

ε

ϕϕ=u_θ

r cot

θ

^+u_r^{r= F}_R^{, (2.9)}

andall remaining strain components vanish.Thus, asex- pected,thestrainstatewithintheinhomogeneityispurely volumetric, and the resulting stress state purely hydro- static.

Inparticular,theradialstresswithintheinhomogeneity is

σ

^rr=2

μ

1

ε

^rr+

λ

1tr

ε

⁼²

μ

1

F R+3

λ

1

F R=3K₁F

R, (2.10) where

λ

1,

μ

1andK₁areLameconstantsandbulkmodulus oftheinhomogeneity.

The obtained results for displacement jumps and stressesarequitenaturallyusedtodeterminetheconstant F andthe effective bulk modulus. To this end, the inter- phaseconditionsofEqs.(2.1)and(2.2)arenowwrittenin thefollowingform

σ

^rr

|

^r=R =Kn

^ur. (2.11)

Accountingfor(2.8)and(2.10)theaboveconditionis 3K1

F

R =Kn[

β (

^R+h

)

^{−F]. (2.12)}

Introducing the dimensionless parameter

δ

= ^hR and normalized springstiffnessinnormal directionkn, where kn=RKn,wedeterminetheconstantFfromEq.(2.12) F= kn

β

^R[1+

δ

^]

3K1+kn . ^(2.13)

The surfaceelementdS₁inEq.(2.1)₁ is proportionalto R²whereasdS₂isproportionalto[R+h]².Thus,combining Eqs.(2.1)₁,(2.10)and(2.13),oneobtains

σ

2= 1

[1+

δ

^]2

σ

^rr

|

r=R = 3K₁kn

β

[1+

δ

^][3K1+kn], ^(2.14) Consequently,thebulk modulusofequivalentinhomo- geneityreads

Keq≡

σ

2

3

β

^{= [1}+

δ

^]K[13kKn1+kn]. ^(2.15) Remark. Inthelimitingcase,ifthethicknessoftheinter- phasehisassumedtobenegligiblysmallincomparisonto the particleradius r(i.e.,if

δ

→0), thebulk modulus Keq

obtainedhereis identicalwiththat ofHashin(1991),de- termined by composite assembly approachand withthat of Duanetal., (2007a) computedonthe basis ofEshelby solution.

2.3.2. Shearmodulus

In order to evaluate effective shear modulus we con- sider a homogeneous deviatoric average deformation for the inhomogeneity/spring layer system. One possibilityis toassume

ε

^eq=

₋

β

^{0 0}

0 −

β

⁰

0 0 2

β

. ^(2.16)

Inthiscase,thedisplacementsontheinterphase/matrix surfaceread

u2=

ε

^eq·r

|

r=R+h =

₋

β

^x2

−

β

^y2

2

β

^z2

=

_u

2[x]

u_2[y] u2[z]

, (2.17)

withradialandtangentialcomponents u_2[r]=[u_2[x]sin

θ

^+u2[z]cos

θ

^]

=

β

^r[2cos2

θ

−sin²

θ

^]

|

r=R+h

=

β

^[R+h][3cos2

θ

^−1],

u_2[θ]=[u_2[x]cos

θ

^−u2[z]sin

θ

^]

=−3

β

^rcos

θ

^sin

θ|

r=R+h=−3

β

^[R+h]cos

θ

^sin

θ

^,

u_2[ϕ]=−[u_2[x]sin

ϕ

^+u2[z]cos

ϕ

^]

=

β

^r[cos

ϕ

^sin

ϕ

−cos

ϕ

^sin

ϕ

^]

|

r=R+h=0. ^(2.18) Displacementsontheinhomogeneity’ssurfaceforaho- mogeneous deviatoric deformation (accounting that the displacement atthe particlecenter iszero) canbe deter- minedfromLurie’ssolution(Lurie2005).Forsphericalin- homogeneitywithradiusR

ur=

12

ν

1Ar² R² +2B

r

₃

2cos²

θ

−1 2

, u_θ =−

(

⁷⁻⁴

ν

¹

)

^A_R^r22 +B

3rcos

θ

^sin

θ

^{, u}ϕ=0. (2.19)

Onthesurfaceofinhomogeneity (r=R) displacements are

u1[r]

|

r=R=[12

ν

1A+2B]R

₃

2cos²

θ

−1 2

,

u_1[θ]

|

^r=R=−[

(

⁷⁻⁴

ν

¹

)

^A+2B]3Rcos

θ

^sin

θ

^{, (2.20)}

and they result in the following jumps across the inter- phase

^ur=u2[r]−u1[r]=[2

β (

^R+h

)

⁻

(

¹²

ν

1A+2B

)

^R]

×

₃

2cos²

θ

−1 2

,

^u_θ =u_2[θ]−u_1[θ]=3[−

β (

^R+h

)

+

( (

⁷⁻⁴

ν

¹

)

^A+B

)

^R]cos

θ

^sin

θ

^{. (2.21)}

ThosejumpsallowtoevaluateonesideofEq.(2.1)₂.To evaluate the other side, the strainsin the inhomogeneity arefirstdeterminedtobe

ε

rr =

∂

^ur

∂

^{r =}

18

ν

1Ar²

R² +B 3cos²

θ

−1 , ^(2.22)

ε

θθ = 1 r

∂

^uθ

∂θ

⁺

ur

r =−3

2

(

⁷⁻⁷

ν

¹

)

^A_R^r2₂^+B

cos²

θ

+

3

(

⁷⁻⁶

ν

1

)

^A_R^r22 +2B

, ^(2.22)

ε

ϕϕ = u_θ

r cot

θ

^+u_r^{r =}

−3

(

⁷⁻¹⁰

ν

¹

)

^A_R^r2₂

cos²

θ

−

6

ν

1Ar² R² +B

, ^(2.22)

2

ε

rθ =

γ

rθ=1 r

∂

^ur

∂θ

⁺

∂

^uθ

∂

^{r −urθ}

=−6

(

⁷⁺²

ν

1

)

^A_R^r2₂ ^+B

cos

θ

^sin

θ

^{. (2.22)}

All remaining strain components vanish. The trace of thestraintensoris

tr

ε

=

ε

rr+

ε

θθ+

ε

ϕϕ=Ar²

R²[1−2

ν

1]

1−3cos²

θ

, (2.23) whichgives

σ

rr =2

μ

1

ε

rr+

λ

1tr

ε

=2

μ

1

18

ν

1Ar²

R²+B 3cos²

θ

−1

−21

λ

1Ar²

R²[1−2

ν

1]

3cos²

θ

−1 ,

σ

rθ =2

μ

¹

ε

rθ=−6

μ

¹

(

⁷⁺²

ν

¹

)

^A_R^r2₂ ^+B

cos

θ

^sin

θ

^.

(2.24)

These are the only stress components that enter Eq.

(2.1)₂.The constantsAandBarecomputedintermsof

β

fromtheinterfaceconditions(2.1)and(2.2)(seedetailsin

Eqs.(A.8)–(A.10))

⎧ ⎪

⎪ ⎪

⎨

⎪ ⎪

⎪ ⎩

A=

(

¹⁺

δ ) β

M [

(

²

μ

^1+kt

)

^{kn −[2}

μ

^1+kn]kt]

B=

(

¹⁺

δ ) β

M [6

ν

1[kn−

μ

1]kt

−[7

(

²

μ

1+kt

)

⁺⁴

ν

1

( μ

1−kt

)

^{] kn]}

,

(2.25) and

M=−[4

μ

^12[7+5

ν

^1]+2

μ

^1[

(

⁷⁻⁴

ν

¹

)

^kn+

(

⁷⁻

ν

¹

)

^kt] +[7−10

ν

^{1]knkt], (2.26)}

where kt=RKt. With the constants A andB defined, the averagestresseswithintheinhomogeneity/interphasesys- temcanbeevaluatedand

μ

^{eqcanbedeterminedfromany}

ofthefollowingequations 2

μ

^eq≡ −S¯xx

β

⁼⁻

S¯yy

β

⁼

S¯zz

2

β

^{, (2.27)}

whereS¯xx,S¯yy andS¯zz arethe average deviatoricstresses.

Consideringthat S¯xx=S¯yy=−2S¯zz only S¯xx will be evalu- ated:

S¯xx=

σ

¯^xx−1₃^tr

σ

^{= 1}₃^[

σ

^¯xx−

σ

^{¯zz]. (2.28)}

Theappropriateaveragestressformulafortheinhomo- geneity/interphasesystemis (BenvenisteandMiloh 2001)

σ

¯i j= 1 V

∂^V2

tixjdS, ^(2.29)

wheret_iisthetractionandtheintegrationisoverthein- terphase/matrixboundary

∂

^V2.Thepositionvectorsonthe interphase/inhomogeneity and interphase/matrix surfaces havethefollowingcomponentsrespectively:

x1=Rsin

θ

^cos

ϕ

^{, y1=Rsin}

θ

^sin

ϕ

^{, z1=Rcos}

θ

^,

x2=[R+h]sin

θ

^cos

ϕ

^{, y2=[R+h]sin}

θ

^sin

ϕ

^,

z2=[R+h]cos

θ

^{. (2.30)}

ConsideringEqs.(2.1)₁ and(2.29),thecomponents

σ

¯^xx and

σ

¯zzneededinEq.(2.28)aredeterminedby

σ

¯^xx=_V¹

2

∂^V2 R²

[R+h]²txx2[R+h]²sin

θ

^d

θ

^d

ϕ

= 1 V2

∂^V2

R²txx2sin

θ

^d

θ

^d

ϕ

, ^(2.31a)

σ

¯^zz= 1 V2

∂^V2

R²tzz2sin

θ

^d

θ

^d

ϕ

. ^(2.31b)

Takingintoaccountthatforr=R

σ

^rr=tr,

σ

rθ=t_θ, ^(2.32)

tx andtzhavethefollowingform tx=[

σ

rrsin

θ

+

σ

rθcos

θ

^]cos

ϕ

,

tz=[

σ

^rrcos

θ

−

σ

rθsin

θ

^]. ^(2.33) Furthermore,fromEq.(2.24)

σ

rr and

σ

rθ evaluated atr

=Rread

σ

^{rr =2}

μ

^1[18

ν

^1A+B]

3cos²

θ

⁻¹

+21·2

μ

¹

ν

^1A

3cos²

θ

⁻¹

=2

μ

1[−3

ν

1A+B]

3cos²

θ

−1 , ^(2.34a)

σ

rθ=−6

μ

1[

(

⁷+2

ν

1

)

^A+B]cos

θ

^sin

θ

. ^(2.34b) Substituting Eq.(2.34) into Eq.(2.33) thetractions on theinhomogeneity’ssurfacebecome

tx=[−6

μ

1A

(

⁷⁺⁵

ν

1

)

^cos2

θ

^sin

θ

+

μ (

⁶

ν

1A−2B

)

^sin

θ

^]cos

ϕ

, ^(2.35a) tz=2

μ

^1[−3

ν

^1A+B]

3cos²

θ

^{−1 cos}

θ

+6

μ

^1[

(

⁷⁺²

ν

¹

)

^{A+B] cos}

θ

^sin2

θ

^{. (2.35b)}

Eqs.(2.35)and(2.30)substitutedintoEq.(2.31)yield

σ

¯^xx= −2

μ

1

5[1+

δ

^{]2[21A+5B], (2.36a)}

σ

¯^zz= 4

μ

1

5[1+

δ

^{]2[21A+5B]. (2.36b)}

Giventhat

μ

^{eqresultingfromEqs.(2.27)and(2.28)is} 2

μ

^{eq =−}

σ

^¯xx−

σ

^¯zz

3

β

⁼

21

μ

¹

5

β

^[1+

δ

^]2[21A+5B]

= 2

μ

1

5

β

^[1+

δ

^{]2[21A+5B], (2.37)}

onegets

2μeq = 2μ1 5[1+δ^]

× 4μ1[7+5ν1][2kn+3k_t]+5[7−10ν1]knk_t

4μ12(⁷+5ν1)+2μ1((⁷−4ν1)^kn+(⁷−ν1)^kt)+(⁷−10ν1)^knk_t .

(2.38)

Itisnotedthatinthelimit,if

δ

→0,anequivalentin- homogeneitywithimperfectinterfaceisobtained.Further- more,

(a) ifkn→∞andkt →∞ 2

μ

^{eq= 2}

μ

¹

[1+

δ

^{], (2.39)}

whichreducesto2

μ

^eq=2

μ

1,if

δ

=0.Thisrepresentsthe perfectinterfacecondition.

(b)ifkn=kt=02

μ

^eq=0,andtheinhomogeneitybehaves likeacavity,

(c)if kn → ∞ but kt is finite (if kt=0, it is called “free sliding”,Hashin(2002a);Duanetal.(2007b))

2

μ

^{eq= 2}

μ

¹

5[1+

δ

^]

8

μ

^1[7+5

ν

^1]+5[7−10

ν

^1]kt

[2

μ

¹

(

⁷⁻⁴

ν

¹

)

⁺

(

⁷⁻¹⁰

ν

¹

)

^kt].

(2.40) Remark. Itshouldbenotedthattheconceptofequivalent inhomogeneity presentedhere isquite differentfrom the ideaemployedbyDuanetal.(2007a)andGuetal.(2014).

Inthosepaperstheequivalentshearmodulus

μ

eqwasde- termined on the basis of Eshelby solution by comparing the change in energy if the equivalent inhomogeneity is inserted in the infinite medium under far-field load and the change due to similar insertion of the original inho- mogeneitytogetherwiththeinterface,(Eshelby1957).Asa result,

μ

eqdependsonthepropertiesofthatmedium(ma- trixmaterial).Thisisnonphysical,andinsharpcontrastto theequivalentinhomogeneityapproachadvocatedby,e.g., Hashin(2002a),andisexplicitlyexcludedintheapproach presentedhere.

The equivalent inhomogeneities described above may beusedinconjunctionwithanymethodemployedineval- uationoftheeffectivepropertiesofcompositeswithoutin- terphases.ThemethodofchoiceinthisworkistheMethod of Conditional Moments (MCM) – a rigorous formal ap- proach to analysis of random composites. The MCM is based on statistical analysis where the random structure of the materialis entirelydescribed by conditional prob- abilities, which – in comparison withthe dataspecifying deterministic structures – isthe onlyadditionaldata.The main features oftheMCMandsome basicresults associ- ated withstandard composites (without interphases) will bebrieflyoutlinedinthesubsequentsection.

3. Effectivepropertiesofrandomcomposites

3.1. Thefundamentalsofthemethodofconditionalmoments

Weexaminearepresentativemacro-volumeVofalin- earelasticcompositeconsistingofamatrixwithrandomly distributed particles. It is assumed that at each point of macro-volumextheHooke’slawisvalid

σ (

^x

)

=C

(

^x

)

^:

ε (

^x

)

, (3.1)

where

σ

^(x)andɛ(x) arethestress andstraintensors,and C(x)isthematerialstiffnesstensor.Hereitisassumedthat C(x)isstatisticallyuniformrandomfunctionofcoordinates withafinitescaleofcorrelation,whoseone-point density distributionhasthefollowingform:

f

(

^C

(

^x

) )

⁼ 2

k=1

c_k

δ (

^C

(

^x

)

−C_k

)

^{, (3.2)}

where

δ

^{(•)denotestheDiracdeltafunction,c}1,c₂andC₁, C₂arevolumefractionsandstiffnesstensorsofthematrix andoftheparticles,respectively.

If the representative volume is under a uniform load, theresultingstressesandstrainsformstatisticallyuniform randomfieldssatisfyingthepropertyofergodicity.Thisal- lows usto replace the averaging over representative vol- ume by averaging over an ensemble of realizations (e.g., Gray 2009; Gnedenko 1962). In this case, the following relationship exists betweenmacroscopic fieldsofstress

σ

¯ andstrain

ε

¯^:

σ

¯ ^=C∗:

ε

^{¯, (3.3)}

whereC∗istheeffectivestiffnesstensor,andtheoverbar

•¯ denotesstatisticalaveraging.

HavingaveragedEq.(3.1)weobtain

σ

¯ ⁼²

k=1

c_kC_k:

ε

¯k,

ε

¯⁼²

k=1

c_k

ε

¯k, ^(3.4)

where

ε

¯k=

ε

(^x)

|

k^(x) istheexpectationvalueofthestrain tensor at point x,provided that this point belong to the kthcomponent.

Itisevident bycomparisonofEqs.(3.3) and(3.4) that inordertodeterminetheeffectivestiffnesstensorC∗itis sufficienttofindtherelationshipbetweenthemeanstrain in the component

ε

¯k (e.g., in the particles

ε

¯1) and the mean strain in the macroscopicvolume

ε

¯^{.In fact, if this} relationshiphastheform:

ε

¯1=A:

ε

¯, ^(3.5)

then,consideringthat

ε

¯=c₁

ε

¯1+c₂

ε

¯2andusingEqs.(3.5) and(3.4),therelationshipfortheeffectivestiffnesstensor reads

C^∗=¯C+c1C₃:

A−⁴I

. ^(3.6)

4

I is the fourth-order unit tensor, ¯C is the expectation valueofthestiffnesstensor

¯C= 2

k=1

c_kC_k, C3=C1−C2. ^(3.7)

To derive the formula for the rank four tensor A of Eq.(3.5),associatedwitharandommaterialthefollowing process is used in the MCM. First, the equilibrium equa- tions forelasticmedium andEq.(3.1)lead tothe follow- ingstochasticdifferentialequationforthefluctuationsu⁰of therandomfieldofdisplacements(Khoroshunetal.,1993;

Nazarenkoetal.2009):

di

v

Cc:sym

∇

^u0

(

^x

)

+di

v

C⁰

(

^x

)

^:

ε (

^x

)

=0,

u⁰

(

^x

) |

∞ =0, (3.8)

where C⁰(^x)=C(^x)−Cc, u⁰(^x)=u(^x)−

ε

¯·x. Although theaboveequationisvalidforanyconstanttensorC_c,the accuracyoftheMCMisenhancedifitisselectedasfollows (Khoroshunetal.,1993;Nazarenkoetal.2015):

Cc=

_{¯C, if C}

1≤C2

C⁻¹

−1

, if C1≥C2

. ^(3.9)

UsingtheGreentensor,thesolutionofEq.(3.9)iswrit- ten inform ofan integral over theinfinite region V (see Nazarenkoetal.2009,2014)

u⁰

(

^x

)

=

Vy

G

(

^x−y

)

·di

v

C⁰

(

^y

)

^:

ε (

^y

)

−

β

dVy, ^(3.10)

where the Green tensor is the solution of the following systemofdifferentialequations:

di

v ₍

Cc:

∇

^G

(

^x

) )

⁺

δ (

^x

)

²I=0, G

(

^x

) |

∞ =0. (3.11) IntegrationofEq.(3.10)bypartsanduseoflinearkine- matic relations,leadstoa stochastic integral equation for therandomstrainfield

ε (

^x

)

=

ε

¯^+K

(

^x−y

)

∗

C⁰

(

^y

)

^:

ε (

^y

)

. (3.12)

The above operator notationfor K(^x−y) ^{has the fol-} lowinginterpretation

K

(

^x−y

)

∗

ψ (

^y

)

=

V

sym

( ∇

^x

( ∇

^xG

(

^x−y

) ) )

^:

ψ (

^y

)

−

ψ

¯

dVy, ^(3.13)

where

ε

¯^and

ψ

^{¯ are mean}(expectation) valuesofɛ(x)and

ψ

^(y).

Applying the technique of conditional averaging (see Khoroshun et al., 1993; Nazarenko et al. 2009) with re- spectto the multipointconditional densities andlimiting theprocessto atwo-point approximation– whichistan- tamounttoassumingidenticalstraindistributionsinallin- homogeneitiesthe followinglinearalgebraicequation, in- volvingconditionalone-pointmomentsofthestrainfields

ε

¯ν,isobtained:

ε

¯ν=

ε

¯+ 2

k=1

K^νk:C⁰_k:

ε

¯k, k,

ν

∈

{

^{1, 2}

}

^{, (3.14)}

where

K^νk=K

(

^x−y

)

^∗pν^k

(

^x−y

)

^{. (3.15)}

The function p_νk(^x−y) ^{denotes the} probability that point x belongs to the kth component, provided point y belongstothe

ν

^thcomponent.

Clearly,Eq.(3.14)along withthe previously usedrela- tionship

ε

¯=c₁

ε

¯1+c₂

ε

¯2,constitutesasystemoftwo(ten- sorial)equationsallowing todetermineboth

ε

¯1 and

ε

¯2 in termsof

ε

¯^{.ThisdefinestensorAofEq.(3.5)andallowsto} evaluatetheeffectivepropertiesaccordingtoEq.(3.6).

3.2.Closed-formformulasforeffectivepropertiesofrandom compositeswithsphericalparticlesandperfectinterfaces

In order to specify Eq. (3.14) one must determine the two-point conditional probabilitiesp_νk(^x−y) ^which characterize shape and arrangement of the inhomo- geneities.ThisallowsevaluatingtheconvolutionK(^x−y)∗ p_νk(^x−y) ^{of Eq. (3.15) following the general formula} showninEq.(3.13).

Detailsoftheevaluationoftheeffectivestiffnesstensor forcomposites withrandomly distributedsphericalparti- cles are presented by Khoroshun etal., 1993; Nazarenko etal.,2014,2015.Followingthemethodologydescribedin thosepapers,thefollowing expressionfortensorAofEq.

(3.5)isobtained:

A=⁴I+c2

⁴

I−L:C ⁻¹:L:C3, ^(3.16) ItsuseinEq.(3.6)yields

C^∗=¯C+c1C3:

⁴

I−L:C ⁻¹:[c2L:C3], ^(3.17) where ¯C,C₃aredeterminedinEq.(3.7)and

C=c1C₂+c2C₁−Cc. (3.18) Tensor L of Eq.(3.17) is an isotropic rank four tensor (theclassicalHilltensor,seeMura1987)

L=2b⁴I+a²I⊗²I, (3.19)