Rigorous bounds on the effective moduli of composites and inhomogeneous bodies with negative-stiffness phases

# Rigorous bounds on the effective moduli of composites and inhomogeneous bodies with negative-stiffness phases

Dennis M. Kochmann Graeme W. Milton Graduate Aerospace Laboratories, California Institute of Technology, Pasadena, CA 91125, USA Department of Mathematics, The University of Utah, Salt Lake City, UT 84112, USA
###### Abstract

We review the theoretical bounds on the effective properties of linear elastic inhomogeneous solids (including composite materials) in the presence of constituents having non-positive-definite elastic moduli (so-called negative-stiffness phases). We show that for statically stable bodies the classical displacement-based variational principles for Dirichlet and Neumann boundary problems hold but that the dual variational principle for traction boundary problems does not apply. We illustrate our findings by the example of a coated spherical inclusion whose stability conditions are obtained from the variational principles. We further show that the classical Voigt upper bound on the linear elastic moduli in multi-phase inhomogeneous bodies and composites applies and that it imposes a stability condition: overall stability requires that the effective moduli do not surpass the Voigt upper bound. This particularly implies that, while the geometric constraints among constituents in a composite can stabilize negative-stiffness phases, the stabilization is insufficient to allow for extreme overall static elastic moduli (exceeding those of the constituents). Stronger bounds on the effective elastic moduli of isotropic composites can be obtained from the Hashin-Shtrikman variational inequalities, which are also shown to hold in the presence of negative stiffness.

###### keywords:
Stability, Elasticity, Composite, Negative Stiffness, Effective Properties, Bounds

## 1 Introduction

The overall or effective properties of heterogeneous solids are uniquely linked to the properties of each composite constituent, their geometric arrangement and bonding. Owing to microstructural randomness in the arrangement of composite phases, effective physical properties in most cases cannot be determined exactly. One approach is to estimate them by the aid of rigorous upper and lower bounds. The simplest such bounds were introduced by Hill (1952) and Paul (1960): the Reuss and Voigt bounds are solely based on phase volume fractions and present upper and lower bounds on the effective linear elastic moduli of multi-phase composites. For inhomogeneous bodies the analogous bounds were obtained by Nemat-Nasser and Hori (1993) and by Willis in a 1989 private communication to Nemat-Nasser and Hori. Based on variational principles and the introduction of a polarization field, Hashin and Shtrikman (1963) derived new tighter bounds for isotropic well-ordered two-phase composites with bulk moduli and shear moduli . Their bounds on the effective bulk modulus can be attained e.g. by assemblages of coated spheres (interchanging the materials in spherical inclusions and coatings yields upper and lower bounds on the effective bulk modulus). Similarly, hierarchical laminate constructions have been shown to attain the bounds on the effective shear modulus (Norris, 1985; Milton, 1986; Francfort and Murat, 1986). Therefore, the bounds of Hashin and Shtrikman (1963) are optimal and present the strongest possible restrictions on the elastic moduli of well-ordered multi-phase solids based only on volume fractions. For non-well-ordered isotropic two-phase composites (with bulk moduli and shear moduli ) the tightest known bounds on the effective bulk modulus are those of Hill (1963a), and are attained by coated-sphere assemblages, and on the effective shear modulus are those of Milton and Phan-Thien (1982), which improve upon those of Walpole (1966) and are attained in certain parameter regimes where they coincide with the Hashin-Shtrikman formulae. By including statistical microstructural information of random composites, three-point bounds were derived e.g. by Beran and Molyneux (1966) and McCoy (1970), who used classical variational principles. For two-phase composites these bounds were simplified by Milton (1981). By improving McCoy’s bounds, Milton and Phan-Thien (1982) found stronger restrictions for the effective shear modulus. We refer to (Cherkaev, 2000; Torquato, 2002; Allaire, 2002; Milton, 2002; Tartar, 2010) for comprehensive reviews of composite bounds. Alternatively, estimates of effective composite properties have been established by an effective medium strategy, which has resulted in, among others, the self-consistent method (Hill, 1965; Budiansky, 1965; Berryman, 1980) and its generalized form (Christensen and Lo, 1979), the differential (Roscoe, 1952, 1973; Norris, 1985) and Mori-Tanaka schemes (Mori and Tanaka, 1973; Benveniste, 1987). Although beyond the scope of the present investigation, we note that nonlinear variational bounds on composite properties are available as well, see e.g. (Talbot and Willis, 1985; Ponte Castaneda, 1991; Castaneda and Willis, 1999).

All of the aforementioned bounds imply that the effective linear elastic moduli of composites (in particular the Young, bulk, and shear moduli of isotropic composites) are bounded from above by the individual moduli of the constituent materials; i.e. no composite can be stiffer than its stiffest constituent. This prohibits the creation of new composites with extreme properties (where by ‘extreme’ we refer to properties which exceed those of the constituents). However, the derivation of those bounds assume that all constituent materials possess positive-definite elastic moduli (for the specific case of isotropic solids, this is equivalent to requiring Young, bulk and shear moduli to be positive). Lakes and Drugan (2002) showed that relaxing this assumption by allowing for non-positive-definite elastic moduli (so-called negative stiffness) in one of the phases in an inhomogeneous body may lead to extreme effective stiffness. Based on exact solutions for a coated-sphere two-phase solid, they showed that a two-phase inhomogeneous body can, in principle (within the validity of the elasticity model), attain unbounded effective bulk stiffness if the constituent moduli and volume fractions are appropriately tuned. Lakes and coworkers demonstrated that various other effective physical composite properties promise to reach extreme values when including a negative-stiffness phase (Lakes, 2001a, b; Wang and Lakes, 2001, 2004, 2005). Experimentally, negative stiffness has been realized by constituents undergoing microscale instabilities such as phase transitions, see e.g. (Lakes et al., 2001; Jaglinski et al., 2006; Jaglinski and Lakes, 2007; Jaglinski et al., 2007). Similarly, on a structural level the negative-stiffness effect has been realized by buckling instabilities, see e.g. (Moore et al., 2006; Lee et al., 2007; Lee and Goverdovskiy, 2012; Kashdan et al., 2012).

While negative stiffness is generally unstable in homogeneous solids with mixed or pure-traction boundary conditions (Kirchhoff, 1859), it was shown that negative-stiffness phases can be stabilized when geometrically constrained e.g. by a sufficiently stiff and thick coating or as inclusions in a stiff matrix (Drugan, 2007; Kochmann and Drugan, 2009; Kochmann, 2012; Kochmann and Drugan, 2012). Unfortunately, the thus expanded stability regime is insufficient to stabilize extreme effective static stiffness in simple two-phase solids and in isotropic two-phase composites with equal shear moduli (Wojnar and Kochmann, 2013b, a), while allowing for interesting dynamic phenomena.

Instead of investigating particular composite geometries, here we show that arbitrary linear elastic inhomogeneous bodies and multi-phase composites cannot reach extreme stiffness by the inclusion of negative-stiffness phases if they are to be statically stable. To this end, we first review the classical variational principles in Section 2 and determine their validity in the presence of negative-stiffness phases. We illustrate the applicability or non-applicability of the various variational principles in Section 3 by the example of a coated spherical inclusion. Next, in Section 4 we apply the variational principles to show that the classical Voigt upper bound applies and, most importantly, implies a stability condition: overall stability requires that the effective moduli must not surpass the Voigt upper bound. We further show that the Hashin-Shtrikman variational inequalities apply and that they yield additional upper and lower bounds on the effective elastic moduli of isotropic composites. Our results particularly demonstrate that extreme effective (static) elastic moduli exceeding those of any of the constituents are prohibited if the solid is to be statically stable. We further review three-point bounds on the effective moduli and conclude upper and lower bounds in the presence of negative-stiffness phases. Finally, Section 5 concludes our analysis.

## 2 Stability conditions for elastic solids

### 2.1 Dirichlet problem: essential boundary conditions

We consider an inhomogeneous body containing a composite material made of linear (visco)elastic constituents and experiencing a displacement field with denoting position in -dimensional space and being time. Assume is the displacement field corresponding to a solution of the elasticity equation (linear momentum balance) with essential boundary conditions on the body’s boundary and given eigenstrains within the solid. The strain energy density of the linear elastic solid with locally-varying modulus tensor is given by

 Ψ(ε)=12[ε(x)−ε0(x)]⋅C(x)[ε(x)−ε0(x)], (1)

where is the infinitesimal symmetric strain tensor, is the infinitesimal symmetric eigenstrain tensor. We will see that a necessary condition for stability is that minimizes the total elastic energy

 W=infu(x)u(x)=v(x) on ∂Ω∫ΩΨ(ε(x))dV. (2)

To show this, assume that the solution is not the minimizer of (2) and that there is some which satisfies the essential boundary conditions and for which

 ∫ΩΨ(^ε(x))dV<∫ΩΨ(~ε(x))dV. (3)

Further, assume the inhomogeneous body has mass density and for simplicity is viscoelastic. Let the body have an initial displacement field at time given by

 u(x,0)=uinitial(x)=~u(x)+η^u(x)1+η (4)

with some . Note that (4) satisfies the essential boundary conditions as well, i.e. on . We fix the displacements to be for all times by application of appropriate body forces which we remove for all , so that for the body is out of equilibrium, i.e.

 u(x,t)={uinitial(x)for t≤0,unknownfor t>0. (5)

We maintain essential boundary conditions

 u(x,t)=v(x) on ∂Ω (6)

for all times , so that no work is done on the body and either internal motions will be damped through viscosity and drive the solid into a state of stable equilibrium, or alternatively there will be no stable equilibrium. Displacements (5) result in strains for , so that the initial energy of the body is purely elastic and equal to

 Winitial=∫ΩΨ(εinitial(x))dV=∫Ω12~ε(x)−ε0(x)+η[^ε(x)−ε0(x)]1+η⋅C(x)~ε(x)−ε0(x)+η[^ε(x)−ε0(x)]1+ηdV=1(1+η)2[∫ΩΨ(~ε(x))dV+η2∫ΩΨ(^ε(x))dV+η∫Ω[^ε(x)−ε0(x)]⋅C(x)[~ε(x)−ε0(x)]dV]. (7)

Because is a solution to the equilibrium equation, we use linear momentum balance in the absence of body forces and in static equilibrium, i.e.

 div~σ=div[C(x)(~ε(x)−ε0(x))]=0in Ω. (8)

Utilizing symmetry of the infinitesimal stress tensor and using (8), we see that

Consequently,

 Winitial=1(1+η)2[(1+2η)∫ΩΨ(~ε(x))dV+η2∫ΩΨ(^ε(x))dV]=∫ΩΨ(~ε(x))dV−η2(1+η)2[∫ΩΨ(~ε(x))dV−∫ΩΨ(^ε(x))dV]. (10)

Due to assumption (3), we know that the final term in brackets is positive and therefore

 Winitial<∫ΩΨ(~ε(x))dV. (11)

Since the energy inside cannot exceed its initial value for reasons of energy conservation, the energy can never approach the value so that, if (3) holds, cannot be the solution as . Hence, if is a stable solution of linear momentum balance, it must be the minimizer of (2). We note that we did not constrain to be positive-definite at any point in our proof. Therefore, if in the presence of negative stiffness in a heterogeneous solid a stable equilibrium solution to the Dirichlet problem exists, it must be the minimizer of (2).

### 2.2 Sufficiency of the stability conditions for the Dirichlet problem

We showed above that a necessary condition for stability is that the energy is minimized. Let us demonstrate in an elastodynamic setting, ignoring viscoelasticity, that this is indeed a sufficient condition of stability if the energy still has a minimum when we perturb the elasticity tensor by subtracting a small constant tensor

 δCijkl=η2(δikδjl+δilδjk) (12)

from it, where is a small parameter. If is isotropic, this perturbation corresponds to subtracting a small constant from the shear modulus while leaving the Lamé modulus unchanged. Suppose initially at time we begin with a small perturbation of the solution which minimizes the energy and no energy enters the domain. Initially at time we could also have a small velocity . Let be the displacement field at times . We assume that at the boundary remains zero for all times and that the density is bounded below by some constant . The total energy, i.e. the sum of elastic and kinetic parts, equals the minimum energy plus a small perturbation (due to and ) and must be conserved, i.e. is constant. We want to show that and the velocity remain small in an sense for all times. Since the elastic part cannot be less than the minimum , we conclude that the kinetic part must be at most the small perturbation , implying

 ρ02∫Ω(δ˙u(x,t))2 dV≤δW, (13)

and the elastic energy is at most +, giving

 δW≥W[~u+δu(x,t)]−W[~u]=12∫∂Ωuδε(x,t)⋅C(x)δε(x,t) dV=η∫Ωδε⋅δε dV+12∫∂Ωuδε⋅[C−δC]δε dV≥η∫Ωδε⋅δε dV. (14)

Here we have used the fact that satisfies the equilibrium equation, and that

 ∫∂Ωuδε(x,t)⋅[C(x)−δC]δε(x,t) dV≥0, (15)

which is a necessary condition for a mimimum to exist when the elasticity tensor is .

Following Ericksen and Toupin (1956), we can use the fact that on (with outward unit normal ) by writing

 0=∫∂Ω(δuiδuj,j−δujδui,j)nidS=∫Ω(δuiδuj,j−δujδui,j),idV=∫Ω((δui,i)2−δuj,iδui,j)dV (16)

so that

 ∫Ωδui,jδuj,idV=∫Ω(δui,i)2dV≥0 (17)

and therefore

 ∫Ωδε⋅δε dV=12∫Ω(δui,jδui,j+δui,jδuj,i)dV≥12∫Ωδui,jδui,j dV. (18)

Finally, using Poincaré’s inequality there exists a constant such that

 ∫Ωδui,j⋅δui,j dV≥CΩ∫Ωδui⋅δui dV, (19)

which allows us to conclude that

 ∫Ωδui⋅δui dV≤2δWηCΩ. (20)

From the inequalities (13) and (20) it is evident that and remain small in an sense for all times. Presumably, if we were to add viscoelasticity the displacement would damp to zero as .

### 2.3 Neumann problem: natural boundary conditions

As shown in Fig. 1a, we consider an inhomogeneous body containing a composite material made of linear elastic constituents, which is completely embedded in and perfectly bonded to another linear elastic surrounding solid (with positive-definite and spatially constant elastic moduli ), i.e.

 C(x)={C(x),if x∈ΩC,C0,if x∈ΩS. (21)

We assume that displacements vanish on the outer surface, i.e.  on , where denotes the entire solid. We assume that eigenstrains act within the surrounding solid but not within the inhomogeneous body:

 ε0(x)={0,if x∈ΩC,ε0(x),if x∈ΩS. (22)

Let us construct the eigenstrains in the following way: consider the surrounding solid with the inhomogeneous body removed. Next, apply tractions on the inner boundary and apply zero displacements on the outer boundary, as shown in Fig. 1b. The resulting strain field that balances the applied tractions in equilibrium will be taken as our eigenstrain , so that because of linear momentum balance we have

 divσ0(x)=div[C0ε0(x)]=0in ΩS (23)

and on the inner boundary with unit normal pointing outward from .

Going back to the system consisting of inhomogeneous body and surrounding elastic medium, assume is the displacement field corresponding to a solution of the elasticity equation (linear momentum balance) in the presence of eigenstrains as determined above with essential boundary conditions on the outer boundary . From Section 2.1 we know that a necessary condition for stability is that minimizes the total elastic energy

 W=infu(x)u(x)=0 on ∂Ω∫Ω12[ε(x)−ε0(x)]⋅C(x)[ε(x)−ε0(x)]dV. (24)

For the given problem, we can expand the energy as follows:

 ∫Ω12[ε(x)−ε0(x)]⋅C(x)[ε(x)−ε0(x)]dV=∫Ω12ε(x)⋅C(x)ε(x)dV−∫ΩSε(x)⋅C0ε0(x)dV+∫ΩS12ε0⋅C0ε0(x)dV, (25)

where the last term is constant and independent of the displacement field. The second term can further be reduced by using , which yields

Using that on and inside as well as , we thus obtain

 ∫ΩSε(x)⋅C(x)ε0dV=∫∂ΩCt0(x)⋅u(x)dS. (27)

Moreover, the first term in (25) can be decomposed into strain energy stored in the inhomogeneous body and in the surrounding solid, i.e.

 ∫Ω12ε(x)⋅C(x)ε(x)dV=∫ΩC12ε(x)⋅C(x)ε(x)dV+∫ΩS12ε(x)⋅C0ε(x)dV. (28)

Altogether we thus have

 (29)

where

 W0=∫ΩS12ε0(x)⋅C0ε0(x)dV=const.>0. (30)

Recall that we determined the eigenstresses from the application of tractions on the interface and vanishing displacements on . Now, keep tractions constant and consider a scaling of the elastic moduli in the surrounding solid of the following form:

 C0=α¯¯¯¯C0⇒ε0=C−10σ0=¯¯¯¯C−10σ0α=¯ε0α (31)

with some . This gives

 ∫Ω12[ε(x)−ε0(x)]⋅C(x)[ε(x)−ε0(x)]dV=∫ΩC12ε(x)⋅C(x)ε(x)dV+α∫ΩS12ε(x)⋅¯¯¯¯C0ε(x)dV−∫∂ΩCt0(x)⋅u(x)dS+¯¯¯¯¯¯W0α. (32)

Here, we may choose arbitrarily small so that the second term vanishes (and the final term has no effect even though it grows in an unbounded manner since it is independent of the displacement field). To a good approximation when is extremely small in is the approximate minimizer of

 W=infu(x){∫ΩC12ε(x)⋅C(x)ε(x)dV−∫∂ΩCt0(x)⋅u(x)dS} (33)

and in is approximately the minimizer of

 infu(x)∫ΩS12ε(x)⋅¯¯¯¯C0ε(x)dV (34)

subject to the constraint that on the outer boundary and that, on the boundary , equals the displacement that minimizes (33) – to ensure that this approximate solution for is continuous across . Furthermore, when is infinitesimal the stress in ,

 ~σ(x)=C0[~ε(x)−ε0(x)]=α¯¯¯¯C0~ε(x)−σ0(x) (35)

approaches and the tractions on the interface (seen from ) follow as . Therefore, due to balance of tractions (i.e. ), can be identified with the traction acting on the surface , see Fig. 1c. In summary, from (24) it follows that a necessary condition of stability with traction boundary conditions is that the displacement field is a minimizer of the total potential energy (33) with begin the tractions applied to the surface . Notice that the opposite limit of letting become infinitely large recovers the Dirichlet boundary value problem.

The dual variational principle

 ¯¯¯¯¯¯W=infσ(x)σ(x)n(x)=t(x) on ∂Ω∫Ω12σ(x)⋅C−1(x)σ(x)dV, (36)

which has been used for traction boundary problems and, among others, yields the Reuss lower bound on the effective moduli of inhomogeneous solids and composites does not apply. This will be shown by the aid of an instructive example in Section 3.3.

### 2.4 Stability conditions for homogeneous isotropic linear elastic solids

Stability conditions for homogeneous solids (for with const.) can be obtained from the variational principles shown above. For the Dirichlet problem, principle (2) can be rephrased by taking variations as

 δ2W=12∫ΩCijkl(x)δui,j(x)δuk,l(x)dV≥0∀δu(x)withδu(x)=0on ∂Ω, (37)

For a homogeneous solid with spatially constant elastic moduli, we conclude that

 Cijkl∫Ωδui,j(x)δuk,l(x)dV≥0∀δu(x)withδu(x)=0on ∂Ω. (38)

In the special case of an isotropic solid with Lamé moduli and ( being the shear modulus), we have

 Cijkl=λδijδkl+μ(δikδjl+δilδjk) (39)

so that the stability condition becomes

 λ⟨δu2i,i(x)⟩+μ⟨δui,j(x)δui,j(x)+δui,j(x)δuj,i(x)⟩≥0∀δu(x)withδu(x)=0on ∂Ω, (40)

where . Introducing the infinitesimal rotation tensor gives

 ωijωij=14(ui,j−uj,i)(ui,j−uj,i)=12(ui,jui,j−ui,juj,i), (41)

which, using (17), ultimately leads to

 λ⟨δu2i,i⟩+μ⟨δui,jδui,j+δui,jδuj,i⟩=(λ+μ)⟨δu2i,i⟩+μ⟨2δωijδωij+δui,jδuj,i⟩=(λ+2μ)⟨δu2i,i⟩+2μ⟨δωijδωij⟩≥0. (42)

Because we can make a large twist inside the body and make arbitrarily large while keeping bounded, or alternatively make a large local compression and make arbitrarily large while keeping bounded, a necessary and sufficient condition of stability is given by the well-known conditions of strong ellipticity (Ericksen and Toupin, 1956; Hill, 1957)

 μ>0andλ+2μ>0. (43)

Note that conditions (43) agree with Hadamard’s (1903) necessary conditions of pointwise stability in elastic media, which ensure real-valued wave speeds (Lord Kelvin, 1888).

For the Neumann problem, the variational principle (33) holds, from which we can obtain stability conditions for homogeneous solids again by considering the second variation:

 δ2W=∫ΩCijklδui,j(x)δuk,l(x)dV≥0∀δu(x)≠0. (44)

We may decompose the displacement gradient into its volumetric and deviatoric contributions, i.e.  with and note that . For homogeneous, isotropic, linear elastic solids we have

 Cijkl⟨δui,j(x)δuk,l(x)⟩=λ⟨δu2i,i(x)⟩+μ⟨δui,j(x)δui,j(x)+δui,j(x)δuj,i(x)⟩=(λ+23μ)⟨δu2k,k(x)⟩+2μ⟨δεij(x)devδεij(x)% dev⟩, (45)

and therefore the stability condition becomes

 (λ+23μ)⟨δε2kk⟩+2μ⟨δεdevijδε%devij⟩≥0∀δu(x)≠0. (46)

Because we can take the strain to be constant throughout the body, with either vanishing volumetric strain, or vanishing deviatoric strain, this implies the necessary and sufficient conditions of stability for the Neumann problem,

 μ>0andκ=λ+23μ>0, (47)

which are the well-known conditions of positive-definiteness of the elastic modulus tensor (Kirchhoff, 1859) with denoting the bulk modulus.

In case of homogeneous anisotropic linear elastic solids, the same variational principles apply and the existence of a unique minimizer requires quasiconvexity of the total potential energy, see e.g. (Knops and Stuart, 1984). The resultant necessary and sufficient condition of stability is positive-definiteness of the elastic modulus tensor, i.e.

 ε⋅Cε>0%forallsymmetricsecond−ordertensors ε≠0, (48)

which for isotropy automatically reduces to (47). In case of pure displacement boundary conditions, the necessary and sufficient condition of stability is strong ellipticity of the elastic modulus tensor (Hadamard, 1903), i.e.

 (a⊗n)⋅C(a⊗n)>0for all vectors a,n≠0. (49)

For isotropy this reduces to (43).

## 3 Example: coated spherical inclusion

Consider a two-phase body consisting of a homogeneous, isotropic, linear elastic spherical particle (radius , elastic moduli and ) coated by and perfectly bonded to a concentric coating of outer radius and of a different homogeneous, isotropic, linear elastic material (elastic moduli and ) as schematically shown in Fig. 2a. The system was studied before to derive effective properties and stability conditions, see e.g. (Lakes and Drugan, 2002; Kochmann and Drugan, 2012; Wojnar and Kochmann, 2013b). The same example of a two-phase solid will be used here to demonstrate the applicability and inapplicability of the standard variational principles in the presence of negative-stiffness phases and their relations to the conditions of overall stability. For simplicity, we assume radial symmetry and choose the boundary conditions accordingly.

### 3.1 Dirichlet boundary value problem

We impose a radial displacement field across the entire outer surface with outward unit normal and constant , which results in radial displacements in the inclusion (superscript i) and in the coating (superscript c) of Lamé’s type, viz.

 ui(r) =uir(r)er,uir(r)=Ar (50a) uc(r) =ucr(r)er,ucr(r)=Br+Cr2 (50b)

in spherical coordinates . In static equilibrium, constants , and are determined by application of the boundary and continuity conditions

 ucr(b)=α,σirr(a)=σ% crr(a),uir(a)=ucr(a). (51)

The stress components are determined from the displacements (50) by application of the strain-displacement relation and Hooke’s law for isotropic elasticity, . The equilibrium solution is then given by (50) with

 A=3κc+4μc3κi+4μc+3(a/b)3(κc−κi)αb,B=3κi+4μc3κi+4μc+3(a/b)3(κc−κi)αb,C=3(κc−κi)a33κi+4μc+3(a/b)3(κc−κi)αb (52)

and therefore the pressure on the outer surface follows as

 σcrr(b)=3αbκc(3κi+4μc)+4(a/b)3(κi−κc)μc3κi+4μc+3(a/b)3(κc−κi). (53)

Note that we can also obtain the effective bulk modulus of the associated two-phase assemblage of coated spheres for the Dirichlet problem via (Hashin, 1962)

 (54)

with volume averages . Therefore, an infinite effective bulk modulus is predicted when (Lakes and Drugan, 2002)

 κi⇁−3a3κc+4b3μc3(b3−a3)=κi∞, (55)

i.e. when approaches the root of the denominator in (54) from below. Also, with decreasing inclusion bulk modulus, the effective bulk modulus first goes to zero when the numerator in (54) vanishes, i.e. when

 κi=−4(b3−a3)κcμ% c3b3κc+4a3μc=κi0. (56)

Simple algebraic manipulations show that for all combinations of strongly-elliptic elastic moduli (required for pointwise stability) and radii we have that .

Next, let us verify the applicability of variational principle (2) derived above for the Dirichlet problem. To this end, we introduce a space of displacement field solutions which are continuous inside the solid and satisfy the boundary condition, such that the space of solutions contains the equilibrium solution (50) with (52). For example, consider a displacement field identical to (50) whose coefficients , and are determined by enforcing

 ~ucr(b)=α,~uir(a)=ua,~ucr(a)=ua (57)

for some interface displacement , which results in

 ~A=uaa,~B=αb2−uaa2b3−a3,~C=(uab−αa)a2b2b3−a3 (58)

with an unknown . Note that the correct equilibrium solution (52) is contained herein and attained when choosing

 ua=Aa=3κc+4μc3κi+4μc+3(a/b)3(κc−κi)αab. (59)

The total energy of the two-phase body in the absence of eigenstrains is given by

 ~W=12∫Ω~ε⋅C~εdV=12∫Ω[κ(tr~ε)2+μ~εdev⋅~εdev]dV=12∫a0[κi(tr~εi)2+μi~εidev⋅~εidev]4πr2dr+12∫ba[κc(tr~εc)2+μc~εcdev⋅~εcdev]4πr2dr=4πb3−a3[3a(b3−a3)u2aκi+3(a2ua−b2α)2κ% c+4ab(bua−aα)2μc], (60)

where is the deviatoric strain tensor. According to (2), the equilibrium solution can be found by minimization:

 ~ua=argmin~W⇒∂~W∂ua=0, (61)

which yields (59), i.e. the correct equilibrium solution. To signal whether this equilibrium solution corresponds to an energy minimum, we note that

 ∂2~W∂u2a=8aπ3a3(κi−κc)−b3(3κi+4μc)a3−b3⎧⎪ ⎪ ⎪ ⎪