- Email: [email protected]
Energy dissipation rate and kinetic relations for Eshelby transformations
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Contents lists available at ScienceDirect
Journal of the Mechanics and Physics of Solids journal homepage: www.elsevier.com/locate/jmps
Energy dissipation rate and kinetic relations for Eshelby transformations Manish Vasoya a, Babak Kondori b, Ahmed Amine Benzerga a,c, Alan Needleman a,∗ a
Department of Materials Science & Engineering, Texas A&M University, College Station, TX 77843 USA Exponent, 149 Commonwealth Drive, Menlo Park, CA 94025 USA c Department of Aerospace Engineering Texas A&M University, College Station, TX 77843 USA b
a r t i c l e
i n f o
Article history: Received 31 July 2019 Revised 19 August 2019 Accepted 21 August 2019 Available online xxx Keywords: Eshelby transformation Shear transformation zone Dissipation Kinetic relations
a b s t r a c t Vasoya et al. (J. Appl. Mech., 86, 051005, 2019) gave a general expression for the dissipation associated with a single transforming Eshelby inclusion in a linear elastic solid and showed that requiring the dissipation to be non-negative provides a strong limit on the maximum value of the transformation strain magnitude. Non-negative dissipation is a necessary, but not sufficient, condition for satisfying the Clausius–Duhem inequality which requires the dissipation rate to be non-negative. Here, a general expression is given for the dissipation rate with multiple transforming Eshelby inclusions. When specialized to a single transforming Eshelby inclusion in an infinite solid, the limit on the magnitude of transformation strain for non-negative dissipation rate is one half that for non-negative dissipation. The condition for non-negative dissipation rate is expressed as the product of a configurational force, analogous to the Peach–Koehler force for dislocations and to the J-integral for cracks, times the transformation strain rate. This gives a natural framework for specifying kinetic relations that guarantee a non-negative dissipation rate for Eshelby inclusions. Simple kinetic relations are proposed for the mesoscale modeling of the shear transformation zone (STZ) mechanism of plastic deformation in amorphous solids. Numerical examples are presented to illustrate implementation in a computational framework of a kinetic relation that guarantees a non-negative dissipation rate. Possible limitations of a linear elastic Eshelby inclusion model for the STZ plasticity mechanism are discussed in light of results from atomistic analyses. © 2019 Elsevier Ltd. All rights reserved.
1. Introduction Eshelby (1957, 1959) transforming inclusions have been used to model a wide range of material behaviors, including, phase transformations, e.g. McMeeking and Evans (1982), solute interactions, e.g. Hirth et al. (2017), the stress dependence of a generalized stacking fault energy, e.g. Andric et al. (2019), and the micromechanics of the plastic deformation in polycrystals, e.g. Hill (1965). Of particular focus here is the plastic deformation of amorphous materials, such as metallic glasses, via local atomic re-arrangements commonly termed shear transformation zones (STZs), see, for example, the review of Hufnagel et al. (2016). Within the context of mesoscale continuum mechanics, STZs can be modeled as Eshelby (1957, 1959) inclusions, as pioneered by Bulatov and Argon (1994). Subsequent mesoscale studies include those ∗
Corresponding author. E-mail address: [email protected] (A. Needleman).
https://doi.org/10.1016/j.jmps.2019.103699 0022-5096/© 2019 Elsevier Ltd. All rights reserved.
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 2
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
of Homer and Schuh (2009), Sandfeld et al. (2015),Kondori et al. (2016), Budrikis et al. (2017), Kondori et al. (2018). In particular, Homer and Schuh (2009) first developed a framework for solving general quasi-static boundary value problems for solids deforming plastically by the evolution of a collection of Eshelby inclusions. Kondori et al. (2016, 2018) introduced a method for solving such boundary value problems that combines analytical solutions for the Eshelby inclusion fields with a numerical solution for enforcing boundary conditions. Continuum mechanics mesoscale formulations, as in Homer and Schuh (2009), Kondori et al. (2016), Kondori et al. (2018), are attractive in that they incorporate the discrete nature of the plastic deformation mechanism while being sufficiently computationally efficient to allow for consideration of a large number of STZs, of a variety of loading modes and rates, and of the effects of component or specimen shape. However, being a continuum formulation, a phenomenological kinetic relation needs to be specified for the evolution of the transformation strain As noted by Kondori et al. (2019), using the kinetic relation of Kondori et al. (2018), the dissipation associated with an STZ transformation can be negative. In order to understand why this can occur, Vasoya et al. (2019) carried out an analytical investigation of the dissipation associated with a single Eshelby transformation, not only for shear transformations, but also for more general Eshelby (1957, 1959) transformations. Quite generally, it was found that negative dissipation is not associated with a particular kinetic relation but is associated with the magnitude of the transformation strain. In numerical calculations of STZ evolution (and more generally of the evolution of a wide variety of processes that can be modeled in terms of Eshelby inclusions), kinetic relations for the transformation strain rate are required that ensure a nonnegative dissipation rate. Paralleling the analysis of Vasoya et al. (2019), we formulate the requirements for a non-negative dissipation rate with multiple transforming inclusions. It is found that a non-negative dissipation rate is a more stringent requirement than is the requirement of non-negative dissipation (the integral of the dissipation rate over the time interval of the transformation). The form of the results of our analysis suggests kinetic relations that guarantee a non-negative dissipation rate. A simple form of such a kinetic relation is implemented in the numerical framework of Kondori et al. (2016), Kondori et al. (2018), Kondori et al. (2019) and numerical examples are presented. Perhaps surprisingly, estimates of transformation strains associated with Eshelby (1957, 1959) inclusions from atomistic analyses are larger than allowed by the requirement of a non-negative dissipation rate, for example, Albert et al. (2016), Ju and Atzmon (2014), Dasgupta et al. (2013). Given this discrepancy, we discuss the possibility that a purely mechanical mesoscale framework may not be adequate for modeling STZ plasticity. 2. Energy rate balance Basic ideas of the analysis presented here, but not the connection with dissipation and/or dissipation rate, are in the works of Rice (1975), Gavazza (1977), Mura (1982), Markenscoff (2010), Markenscoff and Ni (2010), Gavazza and Barnett (2018) and, as presented in the Appendix of Gavazza and Barnett (2018), Eshelby. In particular, here we present the rate version of the analysis of Vasoya et al. (2019) and also extend it to allow for multiple transforming Eshelby (1957, 1959) inclusions. We consider a purely mechanical theory, neglect body forces, confine attention to quasi-static loading histories and small deformation theory. In some initial state, we presume that the body can be characterized as a linear elastic solid with a distribution of transformed and/or transforming Eshelby inclusions, and is subject to prescribed tractions Ti0 (t ) on part of its surface, ST , and prescribed displacements Ui0 (t ) on the remaining part of its surface, SU . The stress and displacement fields at time t are written as
σi j (t ) = σˆ i j (t ) +
N
σiTjK (t )
(1)
K=1
ui (t ) = uˆi (t ) +
N
uTi K (t )
(2)
K=1
where the fields with a superscript K denote the infinite body Eshelby fields for inclusion K, and N is the number of transformed and transforming inclusions. The (ˆ ) fields are added to the Eshelby fields to satisfy the boundary conditions. The stress state at a point xq due to the transformation in inclusion K is given by
σiTjK (xq ) = Li jkl klTK (xq )
(3)
where Li jkl = Lkli j is the tensor of elastic moduli and
klTK (xq ) =
−kl∗ + klcK klcK (xq )
xq inside inclusion K xq outside inclusion K
(4) c
∗ is a uniform transformation eigenstrain and K is determined using the Eshelby (1957, 1959) constraint tensor. where kl kl
The displacements, T
c ui K ( xq ),
with klK (xq ) are not.
associated with
klcK (xq )
are continuous across the inclusion surface whereas those associated
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
3
The rate of loading, T˙i0 (t ) on ST and U˙ i0 (t ) on SU , is prescribed. The energy rate balance takes the form
˙ D˙ = W˙ − with
= and
1 2
(5)
σi j i j dV
V
(6)
W˙ =
Ti u˙ i dS
S
(7)
where S = SU + ST and a superposed ( ˙ ) denotes the time derivative. Since Li jkl = Lkli j
˙ =
σi j ˙ i j dV
V
(8)
The rate fields satisfy equilibrium
σˆ˙ i j, j = 0 σ˙ iTj,K j = 0
(9)
We can write
˙ = ˙ in + ˙ out
(10)
˙ in is the integral over all inclusions and ˙ out is the integral over the material outside the inclusions. where
˙ out =
σi j
V out
ˆ˙ i j +
N
˙ iTjK
dV =
ST +SU
K=1
Ti u˙ i dS +
Sinc
Ti u˙ i dS
(11)
but the first surface integral on the right hand side is simply W˙ . Hence,
˙ out − W˙ =
Sinc
Ti u˙ i dS
(12)
At the current time, inside each inclusion
˙ iTjK = −˙ i∗Kj + ˙ icKj and
˙ =
(13)
σi j ˆ˙ i j +
in
V in
N
˙ iTjK
dV =
K=1
V in
σi j ˆ˙ i j +
N
−
˙ i∗K j
+
˙ icjK
dV
(14)
K=1
c In inclusion M, the fields ˆ˙ i j and ˙ i jK with K = M are continuous across the inclusion boundary. Also, in inclusion M,
˙ i∗Kj = 0 for K = M. In inclusion M, define N
˙ iRjM = ˆ˙ i j +
cK −˙ i∗K j + ˙ i j
(15)
K=M
Then
˙ in
M
=
VM
cM σi j ˙ iRjM + −˙ i∗M dV j + ˙ i j
Since σi j, j = 0
˙ in
M
=
SM
σi j n j u˙ Ri M dS −
SM
σi j n j u˙ ∗M i dS +
(16) SM
σi j n j u˙ ci M dS
(17)
Here, VM and SM are the volume and surface of inclusion M and n is the outer unit normal. R c R c Since the displacements ui M and ui M (and the displacement rates u˙ i M and u˙ i M ) are continuous across inclusion boundaries, the first and third surface integrals on the right hand side of Eq. (17) cancel the part of the surface integral in Eq. (12) that is on the surface of inclusion M. This is true for all inclusions, so that
˙ − W˙ = −
N K K=1 S
σi j n j u˙ ∗K i dS
(18)
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 4
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Hence,
˙ = D˙ = W˙ −
N K K=1 S
σi j n j u˙ ∗K i dS
(19)
or changing to a volume integral for an Eshelby inclusion
D˙ =
N K=1
VK
N N N TJ ∗K σi j dV ˙ i j = σˆ i j (t ) + σi j dV ˙ i∗Kj = ZiKj ˙ i∗K j K=1
VK
J=1
(20)
K=1
where ZiKj is the configurational force type quantity conjugate to ˙ i∗K given by j
ZiKj =
VK
σˆ i j (t ) +
N
σ
TJ ij
dV
(21)
J=1
2.1. Relation between the rate of the dissipation and the dissipation rate In this section we relate the rate formulation presented here to the total dissipation formulation given by Vasoya et al. (2019), who confined attention to the special case of a single Eshelby inclusion in an infinite solid subjected to a time independent stress σi∞ at infinity (since only one inclusion is considered, the superscript K is dropped). j In Vasoya et al. (2019) the total dissipation associated with a transformation strain i∗j was calculated, where the dissipation is
D = W − = Ri j i∗j with
Ri j =
1 σi∞j + σiTj Vincl 2
(22)
(23)
Hence,
D˙ = R˙ i j i∗j + Ri j ˙ i∗j =
1 σi∞j ˙ i∗j + σ˙ iTj i∗j + σiTj ˙ i∗j Vincl 2
(24)
In the inclusion
iTj = −i∗j + icj icj = Si jkl kl∗
(25)
where Sijkl is Eshelby’s constraint tensor and, in general Sijkl = Sklij . We can then write
iTj = −Ii jkl + Si jkl kl∗ = Hi jkl kl∗
(26)
where Iijkl is the fourth order identity tensor and
σiTj = Li jmn Hmnkl kl∗ = Mi jkl kl∗
(27)
Although Hijkl = Hklij , Mi jkl = Mkli j . This follows from the symmetry of Hill’s polarization tensor, Barnett and Cai (2018); Hill (1965). As a consequence of this symmetry, σ˙ iTj i∗j + σiTj ˙ i∗j = 2σiTj ˙ i∗j , so that Eq. (24) gives
D˙ =
σi∞j + σiTj ˙ i∗jVincl
(28)
Alternatively, specializing Eq. (20) to the case considered by Vasoya et al. (2019) gives Eq. (28) and integrating with respect to time gives the expression of Vasoya et al. (2019) for the dissipation
D=
1 σi∞j + σiTj i∗jVincl 2
(29)
Consistency between the direct calculation of D of Vasoya et al. (2019) and the calculation of the dissipation rate here and then integrating with respect to time, requires that Mi jkl = Mkli j . Hence, this energy based result can be regarded as another proof of the symmetry of Mijkl . 2.2. Comparison of limits set by the dissipation and by the dissipation rate We compare the limit on the magnitude of the transformation strain set by requiring the total dissipation to be nonnegative to the limit set by requiring the dissipation rate to be non-negative. As in Section 2.1, we confine attention to a single STZ and to the two plane strain cases in Vasoya et al. (2019). Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
5
2.2.1. Plane strain biaxial tension We consider a circular inclusion in plane strain having radius a and Aincl = π a2 . For an equal biaxial strain transformation, ∗ = ∗ = ∗ and ∗ = 0. The transformation stress field inside the circular inclusion is 11 22 12 T T σ11T = σ22 = σ T , σ12 =0
(30)
The stress σ T in the inclusion is given by, Mura (1982)
σT = −
E∗ 2 (1 − ν 2 )
(31)
so that from Vasoya et al. (2019)
∞ ∗ − σββ
D=
E ( ∗ )2 Aincl 2 (1 − ν 2 )
(32)
∞ = σ ∞ + σ ∞ and non-negative dissipation requires where σββ 11 22
( ∗ )2 ≤ 2 ( 1 − ν 2 )
∞ σββ
E
On the other hand from
D˙ =
∞ σββ −
∗
(33)
E∗ A ˙ ∗ (1 − ν 2 ) incl
(34)
∞ , ∗ and ˙ ∗ all positive, a non-negative dissipation rate requires With σββ
∗ ≤ (1 − ν 2 )
∞ σββ
(35)
E
The limiting value of ∗ for non-negative D˙ is one half the limiting value for non-negative D. 2.2.2. Plane strain shear ∗ = ∗ = γ ∗ and ∗ = Here, in plane strain, the circular inclusion of radius a undergoes a shear transformation with 12 21 11 ∗ = 0. The stress field induced by such a transformation is 22 T T T σ11T = σ22 = 0 , σ12 = σ21 = τT
(36)
and from Mura (1982)
τT = −
Eγ ∗ 4 (1 − ν 2 )
(37)
The dissipation D is given by
∞ ∗ D = 2σ12 γ −
E ( γ ∗ )2 Aincl 4 (1 − ν 2 )
(38)
∞ the remote shear stress and non-negative D requires with σ12
( γ ∗ )2 ≤ 4 ( 1 − ν )
∞ σ12
G
γ∗
(39)
where G = E/2(1 + ν ). ∞ , γ ∗ and γ˙ ∗ positive, gives that a non-negative dissipation rate Differentiating Eq. (38) with respect to time, with σ12 requires
γ ∗ ≤ 2 (1 − ν )
∞ σ12
(40)
G
Here also the limiting value of γ ∗ for non-negative D˙ is one half the limiting value for non-negative D. 3. Kinetic relations In STZ dynamics, thermodynamically consistent kinetics is needed for each transforming inclusion. Focusing attention on Eshelby inclusion K, the simplest relation that guarantees non-negative dissipation is a linear relation of the form
˙ i∗Kj = KiKjkl ZklK K Zkl
where is given by Eq. (21) and the transformability tensor any second order tensor b.
(41) KiKjkl
is positive semi-definite, i.e. satisfies
bi j KiKjkl bkl
≥ 0 for
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 6
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
We presume that for an amorphous solid, the relationship in Eq. (41) is isotropic. Then, KiKjkl must take the form
KiKjkl =
1 K κ1 δik δ jl + δil δ jk + κ2K δi j δkl 2
(42)
so that K ˙ i∗Kj = κ1K ZiKj + κ2K Zkk δi j
(43)
where δ ij is the Kronecker delta. Introducing the deviatoric quantities
1 3
1 3
∗K K ˙ i∗Kj = ˙ i∗Kj − ˙ kk δi j ZiKj = ZiKj − Zkk δi j
(44)
Eq. (43) can be rewritten in terms of the deviatoric and spherical parts as
3 2
∗K K ˙ i∗Kj = κSK ZiKj ˙ kk = 3κHK Zkk
(45)
with κSK = 2κ1K /3 and κHK = κ1K /3 + κ2K denoting the shear and volumetric coefficients of the transformability tensor. The dissipation rate in inclusion K can now be written as K K 2 K K 2 D˙ K = ZiKj ˙ i∗K j = κS (ZS ) + κH (ZH )
(46)
where
3 2
K (ZSK )2 = ZiKj ZiKj ZHK = Zkk
(47)
Here, for simplicity, the scalars κSK and κHK are considered to be constants. If both κSK and κHK are non-negative, then D˙ K is guaranteed to be non-negative. However, we note that a positive dissipation rate does not require that κSK and κHK are each positive, only that the sum in Eq. (46) is positive. In particular, ∗K at first increases and for STZs the shear transformation can be accompanied by a transient volume change, so that kk then decreases so that at the end of the transformation there is no volume change. Hence, over part of the time for the transformation κHK < 0. In such a case, non-negative dissipation requires
κSK (ZSK )2 ≥ −κHK (ZHK )2 3.1. Plane strain biaxial tension Here, we consider the same plane strain biaxial tension case with a single transforming inclusion as considered in Section 2.2.1. The evolution equation for ∗ is taken to be
1 B¯
˙ ∗ = Zββ Since
D˙ =
(48)
∞ σββ −
E∗ A ˙ ∗ (1 − ν 2 ) incl
it is convenient to take
1 ˙ = B ∗
Zββ
1
Aincl E/[(1 − ν 2 )]
(49)
2 ∞ 1 (1 − ν )σββ = − ∗ B E
(50)
so that B is a drag-like factor having units of time. Define ∗ max =
∞ σββ
[E/(1 − ν 2 )]
(51)
Then Eq. (56) can be written as
1 B
1 B
∗ ˙ ∗ + ∗ = max
(52)
so that the time evolution of ∗ is ∗ ∗ (t ) = max [1 − exp (−t/B )]
(53)
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
7
3.2. Plane strain shear For the plane strain shear case with a single transforming inclusion in Section 2.2.2, the dissipation rate D˙ is given by
∞ D˙ = 2σ12 −
Eγ ∗ Aincl γ˙ ∗ = Z12 γ˙ ∗ 2 (1 − ν 2 )
(54)
and we take
1 B¯
γ˙ ∗ = Z12 or
γ˙ ∗ =
(55)
1 1 Z12 B Aincl [E/2(1 − ν 2 )]
=
∞ σ12 1 −γ∗ B [E/4(1 − ν 2 )]
(56)
With ∗ γmax =
∞ σ12 [E/4(1 − ν 2 )]
(57)
Eq. (56) can be written as
1 B
1 B
∗ γ˙ ∗ + γ ∗ = γmax
(58)
so that ∗ γ ∗ (t ) = γmax [1 − exp (−t/B )]
(59)
In Eq. (59), as in Eq. (53), the maximum transformation strain is approached asymptotically as t → ∞ with characteristic time B. 4. Numerical implementation and results In order to illustrate the use of a kinetic relation guaranteeing non-negative dissipation, a few results are presented for plane strain tension of a rectangular bar with a collection of potentially active shear transformation sites. The numerical formulation follows that in Kondori et al. (2016, 2018, 2019) except for the kinetic relation and will be summarized here for completeness. Geometry changes are neglected and boundary value problems are solved by superposing the analytical solutions for STZ sites and a numerical solution for the image fields that enforce the prescribed boundary conditions. The rectangular bar of width Wm in the x1 −direction and length L in the x2 −direction is subject to the boundary conditions
σ11 = σ12 = 0 on x1 = 0, Wm u2 = U (t ) and u2 = 0 and
(60)
σ12 = 0 on x2 = L
σ12 = 0 on x2 = 0
(61a) (61b)
and u1 = 0 at one point to prevent rigid body motions. At each instant of time the stress and displacement fields are given by Eqs. (1) and (2). The evolution of STZs is determined by (i) the activation of new STZs and (ii) the accumulation of eigenstrain in active STZs. STZ activation takes place when the shear strain energy density φ at site K reaches a specified value φ K , where for an elastically isotropic solid φ is
φ=
σij σij 2G
=
σM2
(62)
3G
2 = 3σ σ /2 with the stress deviator σ = σ − σ δ /3. The value of φ K , for an where the Mises stress, σ M is given by σM ij kk i j ij ij ij STZ site is selected randomly from a Gaussian distribution with a mean value φ nuc and taken to lie within ± 2φ std , where φ std is the standard deviation of φ nuc . The transformation strain for STZ K is taken to be of the form i∗K (t ) = γ ∗K (t ) pKi j with pKi j = (sKi nKj + nKi sKj ) with sK a unit j
vector along a maximum shear stress direction at activation of STZ K and nK a unit normal to sK . The evolution of γ ∗K (t) is given by an expression of the form of Eq. (55) but for multiple STZs as in Eq. (21) and normalized as in Eq. (56) so that
γ˙ ∗K =
K 1 1 Zsn B Aincl [E/4(1 − ν 2 )]
(63)
K = Z K = Z K nK sK . The orientation of sK , nK is chosen so that Z K where B is taken to have the same value for all STZs and Zsn ns sn ij i j
is positive and therefore so is γ˙ ∗K .
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS 8
ARTICLE IN PRESS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Fig. 1. Schematic representation of dynamic STZ nucleation. (a) A potential STZ activation site with radius a in a material under applied boundary conditions. (b) The STZ in part (a) is activated and has created four new potential STZ activation sites at a distance of 2a, along the maximum shear directions. From Kondori et al. (2019).
When an STZ site activates, new potential STZ activation sites can nucleate in the vicinity of the activated site. The nucleated potential STZ activation sites are placed along the maximum shear directions at locations where the transformation induces the maximum shear stress magnitude in the material surrounding the STZ. There are a maximum of four potential STZ activation sites that can nucleate when STZ K activates. However, fewer than four new activation sites may be introduced because new sites are not permitted to overlap an existing activation site or the domain boundary. Fig. 1 from Kondori et al. (2019) illustrates this process. Within this framework, strains greater than the transformation strain can occur by reactivation of a previously activated STZ site. Except for parameters associated with kinetic relation, the parameter values used in the calculations are the same as those for the un-notched bar calculations in Kondori et al. (2019). The bar width is Wm = 64 nm, the bar length is L = 4Wm = 256 nm, the initial number of potential STZ nucleation sites is 2200 (an initial potential activation site density of 1.36 × 1017 m−2 ), the prescribed strain rate is U˙ /L = 0.1 s−1 , Young’s modulus E = 100 GPa, Poisson’s ratio ν = 0.325, nuc = 1.5 GPa and with φ the STZ radius a is 0.6 nm and φnuc = 20MJ/m3 corresponding to σM std = 0.25φnuc . With these pa = σ being the only non-zero components of σ , rameter values E/[4(1 − ν 2 )] = 27.95 GPa and for a stress state with σ12 21 ij ∞ in Eq. (57) with σ nuc = 0.866 GPa gives an estimate of γ ∗ replacing σ12 max of 0.031. 12 K reaching the value zero. However, the kinetic In some cases, the transformation governed by Eq. (63) is completed by Zsn relation in Eq. (63) can give rise to a situation as in Eq. (59) where a complete transformation is only approached asymptotically. To limit the time over which the transformation is allowed to take place, a uniform termination time tN = 4.61B is introduced for all STZs, which is the time for γ ∗ in Eq. (59) to reach 99% of its asymptotic value. The transformation is terminated at t = tNK + tN , where tNK is the time at which the transformation has begun in STZ K. The time delay for reactivation of each STZ site is taken to be tD = ρ tN , where ρ is a specified constant. Thus, in the calculations here all material characteristic times are taken to be proportional to B. To illustrate the sort of responses that are obtained when implementing the kinetic relation in Eq. (63) in a numerical calculation of STZ plasticity, results are shown for three values of the reactivation time tD , with fixed B, and for two values of B, with fixed reactivation time tD . With B fixed at 0.001, Fig. 2 shows results for various values of the reactivation time tD . Fig. 2(a) computed stress-strain curves are shown with σ ave is the average stress on x2 = L (where U(t) is prescribed). As the reactivation time decreases, the load drops more rapidly. Fig. 2(b) shows the evolution of the number of activated STZ sites. Bursts of STZ activation occur with the bursts of STZ activation associated with the drops in σ ave in Fig. 2(a). The strain interval between bursts of activation increases with increasing tD . For tD = 5 tN and 10 tN there is a clear increase in σ ave between the drops, whereas for tD = 3 tN , the activation bursts occur so close together that the value of σ ave does not increase between them. The delay in the bursts of STZ activation for tD = 10 tN leads to strains at which there are more activated STZs for tD = 10 tN than for tD = 5 tN . Fig. 3 shows the distribution of activated STZ sites along with the number of times a site has been activated at values of σ ave ≈ 0.5 GPa. The distribution of activated sites are aligned in bands. With tD = 3 tN and tD = 5 tN a significant number of STZ sites have been activated twice whereas with tD = 10 tN only very few STZ sites have been activated more than once. For all three values of tD , the activated STZ sites are located in deformation bands. The effect of varying B with tD = 5 tN on the stress-strain response and the evolution of activated STZs is shown in Fig. 4. The results for B = 0.001 are repeated from Fig. 2. Since B is the only characteristic time scale in the calculations here (all other characteristic material times are taken to be proportional to B), a governing non-dimensional quantity is (U˙ /L )B. Hence, the results in Fig. 4 can also be regarded as pertaining to an increased U˙ /L with fixed B. Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS
ARTICLE IN PRESS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
9
Fig. 2. The dependence of the stress-strain response and the evolution of STZs on the delay time for reactivation tD with B = 0.001 s. (a) Stress-strain curves. (b) Evolution of the number of activated STZs with overall strain, U/L.
Fig. 3. Distribution of active STZ sites, with n the number of times a site has been activated, for B = 0.001 and for various values of tD . (a) tD = 3 tN s at σave = 0.519. (b) tD = 5 tN s at σave = 0.514. (c) tD = 10 tN s at σave = 0.546.
In Fig. 4(a), with B increased by an order of magnitude from 0.001 s to 0.01 s, the peak stress is somewhat higher, the stress oscillations occur over a larger strain range and σ ave does not fall to zero over the strain range shown. At a given value of strain, U/L, the number of activated STZs in Fig. 4(b) is smaller for the calculation with B = 0.01 than for the calculation with B = 0.001. However, because much larger values of U/L are reached with B = 0.01 before the σ ave vanishes, the number of activated STZs with B = 0.01 eventually exceeds the number activated with B = 0.001 As described here, it is rather straightforward to implement a kinetic relation that guarantees a non-negative dissipation rate in the computational framework of Kondori et al. (2018, 2019). The results in Fig. 2 with B = 0.001 and tD = 5 TN (also shown in Fig. 4) are qualitatively similar to and quantitatively fairly close to those in Fig. 3 of Kondori et al. (2019) for an un-notched bar. 5. Discussion Phenomenological relations, as specified in Homer and Schuh (2009), Kondori et al. (2016), Kondori et al. (2018), for the transformation strain of an Eshelby (1957, 1959) inclusion do not guarantee a non-negative dissipation rate. Indeed, as noted by Kondori et al. (2019), a negative dissipation rate can be obtained in such calculations using typical values for the transformation strain. Section 3 presents a framework for developing fairly general kinetic relations that guarantee a nonnegative dissipation rate within a purely mechanical theory and simple kinetic relations are specified that are analogous Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 10
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Fig. 4. The dependence of the stress-strain response and the evolution of STZs on B for calculations with tD = 5 tN . (a) Stress-strain curves. (b) Evolution of the number of activated STZs with overall strain, U/L.
to a J − a relation for crack growth or a dislocation velocity-Peach-Koehler force relation for dislocation motion. However, any such kinetic relation limits the magnitude of the transformation strain to be of the order of a stress level divided by an elastic modulus. Much larger strains at a material point are observed in amorphous solids that deform plastically by the STZ mechanism, for example in a shear band. This can only occur within a purely mechanical framework if multiple transformations take place at the same material point. One possibility is multiple activations of the same STZ (as in the implementation in Section 4) and/or by the nucleation of additional STZ sites that overlap the initial site, possibly after sometime delay following completion of the transformation. Atomistic analyses can provide a physical basis for a kinetic relation for transformation strain evolution and atomistic analyses of STZ plasticity have given rise to a deformation mode that can be related to the strain field associated with a transforming Eshelby inclusion. However, the transformation strains so obtained can be much larger than permitted by the requirement of a non-negative dissipation rate. For example, Ju and Atzmon (2014) suggest a STZ transformation shear strain1 of 0.15 based on their analyses of simulated data. Dasgupta et al. (2013) carried out molecular dynamics (MD) analyses of a two dimensional Lennard–Jones solid and found that a shear transformation zone like deformation mode was obtained and the resulting energy associated with the transformations could be fit by the Eshelby inclusion formula with a shear transformation strain of γ ∗ ≈ 0.09. Albert et al. (2016) carried out three dimensional MD calculations and found that the deformation mode was well fit by a transformation strain of the form
= ∗
0.03 0.12 −0.14
0.12 −0.08 −0.07
−0.14 −0.07 0.06
(64)
Clearly, these fits of atomistic calculations to a transformation of an Eshelby (1957, 1959) inclusion give magnitudes of the transformation strain that are much larger than a non-negative dissipation rate (or even non-negative dissipation) allows. Presuming that the atomistic calculations do not violate the requirement of non-negative dissipation rate and the calculations of transformation strain in the atomistic analyses are the appropriate values to relate to the continuum strain measure, there are several possible reasons for the discrepancy between the atomistic analyses and our continuum analysis. For example, it is possible that due to high stresses in the atomistic calculations linear elasticity is not an appropriate constitutive description or that surface energy plays a significant role in the atomistic analyses. In addition, it is worth noting that the predicted transformation strain magnitude with a non-negative dissipation rate increases with increasing stress level (in Sections 2.2.1 and 2.2.2 the limiting strain is of the order of a stress level divided by an elastic modulus). However, transformation strains associated with such large stresses would not be relevant for mesoscale modeling of plastic deformation at much lower stress levels. Another possibility is that, within the continuum context, a purely mechanical theory does not provide an adequate mesoscale framework. A broader thermodynamic framework may be needed. 1 ∗ ∗ , in others it is identified with 212 and in others which of these In some reports, as here, the transformation shear strain is identified with 12 definitions is used is not explicitly specified.
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
11
Within the context of continuum thermodynamics and for a spatially uniform and time independent density (as appropriate for our small strain formulation), the Clausius–Duhem inequality can be written as
˙ − W˙ + s˙ T −
qi ∂ T /∂ xi ≥0 T
(65)
where T is the absolute temperature, qi are the components of the heat flux and s is the entropy. Assuming that qi and/or ∂ T/∂ xi vanish, Eq. (65) simplifies to
˙ ≥0 W˙ + s˙ T −
(66)
There are several possibilities for allowing larger transformation strains. One possibility is that rather than the material response during transformation being an elastic process, it is an elastic-plastic (or more precisely an inelastic process of p some kind), so that the total strain rate during transformation is ˙ i j = ˙ iej + ˙ i j and within a transforming inclusion
˙ + σi j ˙ p W˙ = σi j ˙ i j = σi j ˙ iej + σi j ˙ ipj = ij
(67)
then Eq. (66) becomes
σi j ˙ ipj + s˙ T ≥ 0
(68)
However, it is unlikely that inelastic deformation in the inclusion permits larger transformation strains presuming that the material response is elastic prior to transformation and the transformation reduces the stress level in the inclusion so that the process involves elastic unloading. Furthermore, Vasoya et al. (2019) showed that a reduction in stiffness in the inclusion relative to that of the material outside the inclusion reduces the transformation strain that can occur with nonnegative dissipation. p If σi j ˙ i j = 0 in Eq. (68), what is left is the entropic contribution and, since the rate of change of entropy is always positive, s˙ T gives a positive contribution and thus can contribute to allowing a larger transformation strain for an Eshelby (1957, 1959) inclusion. One source of entropy production is that arising from atomic thermal vibrations in the inclusion. Another source, that could possibly give a contribution, arises from an increase in configurational temperature and entropy due to the structural rearrangement in the transforming inclusion. Indeed, atomistic analyses and experiment indicate that entropy can play a role in dislocation nucleation at finite temperatures, see for example, Chen et al. (2015); Nguyen et al. (2011); Ryu et al. (2011). It is also worth noting that a formulation for using statistical mechanics to relate the entropy of an atomic system to mesoscale quantities for dislocation dynamics has been presented by Acharya (2011). To give an estimate of the possible effect of entropy production, we presume the dominant entropy effect is configurational and take the configurational temperature to be proportional to the transformation strain magnitude and the configurational entropy to be proportional to the configurational temperature, and focus on the plane strain shear case of Sections 2.2.2 and 3.2. With these restrictions, the Clausius-Duhem inequality can be written as ∗ + CD γ ∗ − γ ∗ ]γ˙ ∗ ≥ 0 [γmax
where the mechanical limit strain
γth∗ −lim =
γ
∗ max
1 − CD
(69) ∗ γmax
is defined in Eq. (57) and the thermodynamic limit on the transformation strain is
(70)
∗ so if γmax = 0.03, then CD = 2/3 to CD = 3/4 would allow transformation strain magnitudes in the range 0.09 to 0.12 with non-negative dissipation.
6. Conclusions 1. The analysis of Vasoya et al. (2019) has been modified to give a general expression for the dissipation rate and extended to account for multiple transforming Eshelby (1957, 1959) inclusions. 2. The requirement of non-negative dissipation rate is more stringent than the requirement of non-negative dissipation. • For a single transforming inclusion in an infinite solid with a remote uniform stress, the requirement of non-negative dissipation rate gives a limiting transformation strain that is a factor of two smaller than is given by the requirement of non-negative dissipation. 3. A framework for specifying the evolution of the transformation strain in terms of the associated configurational force was proposed. Such kinetic relations guarantee a non-negative dissipation rate. • Implementing a kinetic relation that guarantees a non-negative dissipation rate into a computational framework is rather straightforward. In a simple quasi-static initial/boundary value problem, and with an appropriate choice of parameter values, results are obtained that are qualitatively consistent with previous results. The implementation here was in the computational framework of Kondori et al. (2018, 2019) but such kinetic relations can be incorporated into other computational frameworks. 4. Atomistic analyses indicate that the magnitude of STZ transformation strains can exceed the value permitted by the requirement of non-negative dissipation rate. There are a variety of possibilities for this discrepancy, one of which is that a purely mechanical mesoscale framework may not be adequate for modeling STZ plasticity. Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS 12
ARTICLE IN PRESS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Acknowledgments We gratefully acknowledge financial support from the Weizmann Institute of Science through the Weizmann-Texas A&M program. References Acharya, A., 2011. Microcanonical entropy and mesoscale dislocation mechanics and plasticity. J. Elast. 104, 23–44. Albert, T., Tanguy, A., Boioli, F., Rodney, D., 2016. Mapping between atomistic simulations and eshelby inclusions in the shear deformation of an amorphous silicon model. Phys. Rev. E 93, 053002. Andric, P., Yin, B., Curtin, W.A., 2019. Stress-dependence of generalized stacking fault energies. J. Mech. Phys. Solids 122, 262–279. Barnett, D.M., Cai, W., 2018. Properties of the eshelby tensor and existence of the equivalent ellipsoidal inclusion solution. J. Mech. Phys. Solids 121, 71–80. Budrikis, Z., Castellanos, D.F., Sandfeld, S., Zaiser, M., Zapperi, S., 2017. Universal features of amorphous plasticity. Nature Commun. 8, 15928. Bulatov, V.V., Argon, A.S., 1994. A stochastic model for continuum elasto-plastic behavior. i. numerical approach and strain localization. Modell. Simul. Mater. Sci. Engin. 2, 167–184. Chen, L.Y., He, M., Shin, J., Richter, G., Gianola, D.S., 2015. Measuring surface dislocation nucleation in defect-scarce nanostructures. Nature Mater. 14, 707–713. Dasgupta, R., Hentschel, H.G.E., Procaccia, I., 2013. Yield strain in shear banding amorphous solids. Phys. Rev. E 87, 022810. Eshelby, J.D., 1957. The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proc. R. Soc. Lond. A 241, 376–396. Eshelby, J.D., 1959. The elastic field outside an ellipsoidal inclusion. Proc. R. Soc. Lond. A 252, 561–569. Gavazza, S.D., 1977. Forces on pure inclusions and Somigliana dislocations. Scr. Metall. 11, 979–981. Gavazza, S.D., Barnett, D.M., 2018. The pointwise Eshelby force on the interface between a transformed inclusion and its surrounding matrix. Math. Mech. Solid 23, 233–239. Hill, R., 1965. Continuum micro-mechanics of elastoplastic polycrystals. J. Mech. Phys. Solids 13, 89–101. Hirth, J.P., Barnett, D., Hoagland, R.G., 2017. Solute atmospheres at dislocations. Acta Mater. 131, 574–593. Homer, E.R., Schuh, C.A., 2009. Mesoscale modeling of amorphous metals by shear transformation zone dynamics. Acta Mater. 57, 2823–2833. Hufnagel, T.C., Schuh, C.A., Falk, M.L., 2016. Deformation of metallic glasses: recent developments in theory, simulations, and experiments. Acta Mater. 109, 375–393. Ju, J.D., Atzmon, M., 2014. A comprehensive atomistic analysis of the experimental dynamic-mechanical response of a metallic glass. Acta Mater. 74, 183–188. Kondori, B., Benzerga, A.A., Needleman, A., 2016. Discrete shear transformation zone plasticity. Extreme Mech. Lett. 9, 21–29. Kondori, B., Benzerga, A.A., Needleman, A., 2018. Discrete shear-transformation-zone plasticity modeling of notched bars. J. Mech. Phys. Solids 111, 18–42. Kondori, B., Benzerga, A.A., Needleman, A., 2019. Errata to: discrete shear-transformation-zone plasticity modeling of notched bars. J. Mech. Phys. Solids 127, 151–153. Markenscoff, X., 2010. Evolution equation of moving defects: dislocations and inclusions. Int. J. Fract. 166, 35–40. Markenscoff, X., Ni, L., 2010. The energy-release rate and “self-force” of dynamically expanding spherical and plane inclusion boundaries with dilatational eigenstrain. J. Mech. Phys. Solids 58, 1–11. McMeeking, R.M., Evans, A.G., 1982. Mechanics of transformation- toughening in brittle materials. J. Am. Ceram. Soc. 65, 242–246. Mura, T., 1982. Micromechanics of Defects in Solids. Martinus Nijhoff, The Hague. Nguyen, L.D., Baker, K.L., Warner, D.H., 2011. Atomistic predictions of dislocation nucleation with transition state theory. Phys. Rev. B. 84, 024118. Rice, J.R., 1975. Continuum Mechanics and Thermodynamics of Plasticity in Relation to Microscale Deformation Mechanisms. In: Argon, A.S. (Ed.), Chapter 2 of Constitutive Equations in Plasticity. M.I.T. Press, pp. 23–79. Ryu, S., Kang, K., Cai, W., 2011. Entropic effect on the rate of dislocation nucleation. Proc. Nat. Acad. Sci. 108, 5174–5178. Sandfeld, S., Budrikis, Z., Zapperi, S., Castellanos, D.F., 2015. Avalanches, loading and finite size effects in 2d amorphous plasticity: results from a finite element model. J. Stat. Mech. P02011. Vasoya, M., Kondori, B., Benzerga, A.A., Needleman, A., 2019. Limits on transformation strains for non-negative dissipation. J. Appl. Mech. 86, 051005.
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS
[m3Gsc;August 28, 2019;13:32]
Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Contents lists available at ScienceDirect
Journal of the Mechanics and Physics of Solids journal homepage: www.elsevier.com/locate/jmps
Energy dissipation rate and kinetic relations for Eshelby transformations Manish Vasoya a, Babak Kondori b, Ahmed Amine Benzerga a,c, Alan Needleman a,∗ a
Department of Materials Science & Engineering, Texas A&M University, College Station, TX 77843 USA Exponent, 149 Commonwealth Drive, Menlo Park, CA 94025 USA c Department of Aerospace Engineering Texas A&M University, College Station, TX 77843 USA b
a r t i c l e
i n f o
Article history: Received 31 July 2019 Revised 19 August 2019 Accepted 21 August 2019 Available online xxx Keywords: Eshelby transformation Shear transformation zone Dissipation Kinetic relations
a b s t r a c t Vasoya et al. (J. Appl. Mech., 86, 051005, 2019) gave a general expression for the dissipation associated with a single transforming Eshelby inclusion in a linear elastic solid and showed that requiring the dissipation to be non-negative provides a strong limit on the maximum value of the transformation strain magnitude. Non-negative dissipation is a necessary, but not sufficient, condition for satisfying the Clausius–Duhem inequality which requires the dissipation rate to be non-negative. Here, a general expression is given for the dissipation rate with multiple transforming Eshelby inclusions. When specialized to a single transforming Eshelby inclusion in an infinite solid, the limit on the magnitude of transformation strain for non-negative dissipation rate is one half that for non-negative dissipation. The condition for non-negative dissipation rate is expressed as the product of a configurational force, analogous to the Peach–Koehler force for dislocations and to the J-integral for cracks, times the transformation strain rate. This gives a natural framework for specifying kinetic relations that guarantee a non-negative dissipation rate for Eshelby inclusions. Simple kinetic relations are proposed for the mesoscale modeling of the shear transformation zone (STZ) mechanism of plastic deformation in amorphous solids. Numerical examples are presented to illustrate implementation in a computational framework of a kinetic relation that guarantees a non-negative dissipation rate. Possible limitations of a linear elastic Eshelby inclusion model for the STZ plasticity mechanism are discussed in light of results from atomistic analyses. © 2019 Elsevier Ltd. All rights reserved.
1. Introduction Eshelby (1957, 1959) transforming inclusions have been used to model a wide range of material behaviors, including, phase transformations, e.g. McMeeking and Evans (1982), solute interactions, e.g. Hirth et al. (2017), the stress dependence of a generalized stacking fault energy, e.g. Andric et al. (2019), and the micromechanics of the plastic deformation in polycrystals, e.g. Hill (1965). Of particular focus here is the plastic deformation of amorphous materials, such as metallic glasses, via local atomic re-arrangements commonly termed shear transformation zones (STZs), see, for example, the review of Hufnagel et al. (2016). Within the context of mesoscale continuum mechanics, STZs can be modeled as Eshelby (1957, 1959) inclusions, as pioneered by Bulatov and Argon (1994). Subsequent mesoscale studies include those ∗
Corresponding author. E-mail address: [email protected] (A. Needleman).
https://doi.org/10.1016/j.jmps.2019.103699 0022-5096/© 2019 Elsevier Ltd. All rights reserved.
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 2
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
of Homer and Schuh (2009), Sandfeld et al. (2015),Kondori et al. (2016), Budrikis et al. (2017), Kondori et al. (2018). In particular, Homer and Schuh (2009) first developed a framework for solving general quasi-static boundary value problems for solids deforming plastically by the evolution of a collection of Eshelby inclusions. Kondori et al. (2016, 2018) introduced a method for solving such boundary value problems that combines analytical solutions for the Eshelby inclusion fields with a numerical solution for enforcing boundary conditions. Continuum mechanics mesoscale formulations, as in Homer and Schuh (2009), Kondori et al. (2016), Kondori et al. (2018), are attractive in that they incorporate the discrete nature of the plastic deformation mechanism while being sufficiently computationally efficient to allow for consideration of a large number of STZs, of a variety of loading modes and rates, and of the effects of component or specimen shape. However, being a continuum formulation, a phenomenological kinetic relation needs to be specified for the evolution of the transformation strain As noted by Kondori et al. (2019), using the kinetic relation of Kondori et al. (2018), the dissipation associated with an STZ transformation can be negative. In order to understand why this can occur, Vasoya et al. (2019) carried out an analytical investigation of the dissipation associated with a single Eshelby transformation, not only for shear transformations, but also for more general Eshelby (1957, 1959) transformations. Quite generally, it was found that negative dissipation is not associated with a particular kinetic relation but is associated with the magnitude of the transformation strain. In numerical calculations of STZ evolution (and more generally of the evolution of a wide variety of processes that can be modeled in terms of Eshelby inclusions), kinetic relations for the transformation strain rate are required that ensure a nonnegative dissipation rate. Paralleling the analysis of Vasoya et al. (2019), we formulate the requirements for a non-negative dissipation rate with multiple transforming inclusions. It is found that a non-negative dissipation rate is a more stringent requirement than is the requirement of non-negative dissipation (the integral of the dissipation rate over the time interval of the transformation). The form of the results of our analysis suggests kinetic relations that guarantee a non-negative dissipation rate. A simple form of such a kinetic relation is implemented in the numerical framework of Kondori et al. (2016), Kondori et al. (2018), Kondori et al. (2019) and numerical examples are presented. Perhaps surprisingly, estimates of transformation strains associated with Eshelby (1957, 1959) inclusions from atomistic analyses are larger than allowed by the requirement of a non-negative dissipation rate, for example, Albert et al. (2016), Ju and Atzmon (2014), Dasgupta et al. (2013). Given this discrepancy, we discuss the possibility that a purely mechanical mesoscale framework may not be adequate for modeling STZ plasticity. 2. Energy rate balance Basic ideas of the analysis presented here, but not the connection with dissipation and/or dissipation rate, are in the works of Rice (1975), Gavazza (1977), Mura (1982), Markenscoff (2010), Markenscoff and Ni (2010), Gavazza and Barnett (2018) and, as presented in the Appendix of Gavazza and Barnett (2018), Eshelby. In particular, here we present the rate version of the analysis of Vasoya et al. (2019) and also extend it to allow for multiple transforming Eshelby (1957, 1959) inclusions. We consider a purely mechanical theory, neglect body forces, confine attention to quasi-static loading histories and small deformation theory. In some initial state, we presume that the body can be characterized as a linear elastic solid with a distribution of transformed and/or transforming Eshelby inclusions, and is subject to prescribed tractions Ti0 (t ) on part of its surface, ST , and prescribed displacements Ui0 (t ) on the remaining part of its surface, SU . The stress and displacement fields at time t are written as
σi j (t ) = σˆ i j (t ) +
N
σiTjK (t )
(1)
K=1
ui (t ) = uˆi (t ) +
N
uTi K (t )
(2)
K=1
where the fields with a superscript K denote the infinite body Eshelby fields for inclusion K, and N is the number of transformed and transforming inclusions. The (ˆ ) fields are added to the Eshelby fields to satisfy the boundary conditions. The stress state at a point xq due to the transformation in inclusion K is given by
σiTjK (xq ) = Li jkl klTK (xq )
(3)
where Li jkl = Lkli j is the tensor of elastic moduli and
klTK (xq ) =
−kl∗ + klcK klcK (xq )
xq inside inclusion K xq outside inclusion K
(4) c
∗ is a uniform transformation eigenstrain and K is determined using the Eshelby (1957, 1959) constraint tensor. where kl kl
The displacements, T
c ui K ( xq ),
with klK (xq ) are not.
associated with
klcK (xq )
are continuous across the inclusion surface whereas those associated
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
3
The rate of loading, T˙i0 (t ) on ST and U˙ i0 (t ) on SU , is prescribed. The energy rate balance takes the form
˙ D˙ = W˙ − with
= and
1 2
(5)
σi j i j dV
V
(6)
W˙ =
Ti u˙ i dS
S
(7)
where S = SU + ST and a superposed ( ˙ ) denotes the time derivative. Since Li jkl = Lkli j
˙ =
σi j ˙ i j dV
V
(8)
The rate fields satisfy equilibrium
σˆ˙ i j, j = 0 σ˙ iTj,K j = 0
(9)
We can write
˙ = ˙ in + ˙ out
(10)
˙ in is the integral over all inclusions and ˙ out is the integral over the material outside the inclusions. where
˙ out =
σi j
V out
ˆ˙ i j +
N
˙ iTjK
dV =
ST +SU
K=1
Ti u˙ i dS +
Sinc
Ti u˙ i dS
(11)
but the first surface integral on the right hand side is simply W˙ . Hence,
˙ out − W˙ =
Sinc
Ti u˙ i dS
(12)
At the current time, inside each inclusion
˙ iTjK = −˙ i∗Kj + ˙ icKj and
˙ =
(13)
σi j ˆ˙ i j +
in
V in
N
˙ iTjK
dV =
K=1
V in
σi j ˆ˙ i j +
N
−
˙ i∗K j
+
˙ icjK
dV
(14)
K=1
c In inclusion M, the fields ˆ˙ i j and ˙ i jK with K = M are continuous across the inclusion boundary. Also, in inclusion M,
˙ i∗Kj = 0 for K = M. In inclusion M, define N
˙ iRjM = ˆ˙ i j +
cK −˙ i∗K j + ˙ i j
(15)
K=M
Then
˙ in
M
=
VM
cM σi j ˙ iRjM + −˙ i∗M dV j + ˙ i j
Since σi j, j = 0
˙ in
M
=
SM
σi j n j u˙ Ri M dS −
SM
σi j n j u˙ ∗M i dS +
(16) SM
σi j n j u˙ ci M dS
(17)
Here, VM and SM are the volume and surface of inclusion M and n is the outer unit normal. R c R c Since the displacements ui M and ui M (and the displacement rates u˙ i M and u˙ i M ) are continuous across inclusion boundaries, the first and third surface integrals on the right hand side of Eq. (17) cancel the part of the surface integral in Eq. (12) that is on the surface of inclusion M. This is true for all inclusions, so that
˙ − W˙ = −
N K K=1 S
σi j n j u˙ ∗K i dS
(18)
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 4
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Hence,
˙ = D˙ = W˙ −
N K K=1 S
σi j n j u˙ ∗K i dS
(19)
or changing to a volume integral for an Eshelby inclusion
D˙ =
N K=1
VK
N N N TJ ∗K σi j dV ˙ i j = σˆ i j (t ) + σi j dV ˙ i∗Kj = ZiKj ˙ i∗K j K=1
VK
J=1
(20)
K=1
where ZiKj is the configurational force type quantity conjugate to ˙ i∗K given by j
ZiKj =
VK
σˆ i j (t ) +
N
σ
TJ ij
dV
(21)
J=1
2.1. Relation between the rate of the dissipation and the dissipation rate In this section we relate the rate formulation presented here to the total dissipation formulation given by Vasoya et al. (2019), who confined attention to the special case of a single Eshelby inclusion in an infinite solid subjected to a time independent stress σi∞ at infinity (since only one inclusion is considered, the superscript K is dropped). j In Vasoya et al. (2019) the total dissipation associated with a transformation strain i∗j was calculated, where the dissipation is
D = W − = Ri j i∗j with
Ri j =
1 σi∞j + σiTj Vincl 2
(22)
(23)
Hence,
D˙ = R˙ i j i∗j + Ri j ˙ i∗j =
1 σi∞j ˙ i∗j + σ˙ iTj i∗j + σiTj ˙ i∗j Vincl 2
(24)
In the inclusion
iTj = −i∗j + icj icj = Si jkl kl∗
(25)
where Sijkl is Eshelby’s constraint tensor and, in general Sijkl = Sklij . We can then write
iTj = −Ii jkl + Si jkl kl∗ = Hi jkl kl∗
(26)
where Iijkl is the fourth order identity tensor and
σiTj = Li jmn Hmnkl kl∗ = Mi jkl kl∗
(27)
Although Hijkl = Hklij , Mi jkl = Mkli j . This follows from the symmetry of Hill’s polarization tensor, Barnett and Cai (2018); Hill (1965). As a consequence of this symmetry, σ˙ iTj i∗j + σiTj ˙ i∗j = 2σiTj ˙ i∗j , so that Eq. (24) gives
D˙ =
σi∞j + σiTj ˙ i∗jVincl
(28)
Alternatively, specializing Eq. (20) to the case considered by Vasoya et al. (2019) gives Eq. (28) and integrating with respect to time gives the expression of Vasoya et al. (2019) for the dissipation
D=
1 σi∞j + σiTj i∗jVincl 2
(29)
Consistency between the direct calculation of D of Vasoya et al. (2019) and the calculation of the dissipation rate here and then integrating with respect to time, requires that Mi jkl = Mkli j . Hence, this energy based result can be regarded as another proof of the symmetry of Mijkl . 2.2. Comparison of limits set by the dissipation and by the dissipation rate We compare the limit on the magnitude of the transformation strain set by requiring the total dissipation to be nonnegative to the limit set by requiring the dissipation rate to be non-negative. As in Section 2.1, we confine attention to a single STZ and to the two plane strain cases in Vasoya et al. (2019). Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
5
2.2.1. Plane strain biaxial tension We consider a circular inclusion in plane strain having radius a and Aincl = π a2 . For an equal biaxial strain transformation, ∗ = ∗ = ∗ and ∗ = 0. The transformation stress field inside the circular inclusion is 11 22 12 T T σ11T = σ22 = σ T , σ12 =0
(30)
The stress σ T in the inclusion is given by, Mura (1982)
σT = −
E∗ 2 (1 − ν 2 )
(31)
so that from Vasoya et al. (2019)
∞ ∗ − σββ
D=
E ( ∗ )2 Aincl 2 (1 − ν 2 )
(32)
∞ = σ ∞ + σ ∞ and non-negative dissipation requires where σββ 11 22
( ∗ )2 ≤ 2 ( 1 − ν 2 )
∞ σββ
E
On the other hand from
D˙ =
∞ σββ −
∗
(33)
E∗ A ˙ ∗ (1 − ν 2 ) incl
(34)
∞ , ∗ and ˙ ∗ all positive, a non-negative dissipation rate requires With σββ
∗ ≤ (1 − ν 2 )
∞ σββ
(35)
E
The limiting value of ∗ for non-negative D˙ is one half the limiting value for non-negative D. 2.2.2. Plane strain shear ∗ = ∗ = γ ∗ and ∗ = Here, in plane strain, the circular inclusion of radius a undergoes a shear transformation with 12 21 11 ∗ = 0. The stress field induced by such a transformation is 22 T T T σ11T = σ22 = 0 , σ12 = σ21 = τT
(36)
and from Mura (1982)
τT = −
Eγ ∗ 4 (1 − ν 2 )
(37)
The dissipation D is given by
∞ ∗ D = 2σ12 γ −
E ( γ ∗ )2 Aincl 4 (1 − ν 2 )
(38)
∞ the remote shear stress and non-negative D requires with σ12
( γ ∗ )2 ≤ 4 ( 1 − ν )
∞ σ12
G
γ∗
(39)
where G = E/2(1 + ν ). ∞ , γ ∗ and γ˙ ∗ positive, gives that a non-negative dissipation rate Differentiating Eq. (38) with respect to time, with σ12 requires
γ ∗ ≤ 2 (1 − ν )
∞ σ12
(40)
G
Here also the limiting value of γ ∗ for non-negative D˙ is one half the limiting value for non-negative D. 3. Kinetic relations In STZ dynamics, thermodynamically consistent kinetics is needed for each transforming inclusion. Focusing attention on Eshelby inclusion K, the simplest relation that guarantees non-negative dissipation is a linear relation of the form
˙ i∗Kj = KiKjkl ZklK K Zkl
where is given by Eq. (21) and the transformability tensor any second order tensor b.
(41) KiKjkl
is positive semi-definite, i.e. satisfies
bi j KiKjkl bkl
≥ 0 for
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 6
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
We presume that for an amorphous solid, the relationship in Eq. (41) is isotropic. Then, KiKjkl must take the form
KiKjkl =
1 K κ1 δik δ jl + δil δ jk + κ2K δi j δkl 2
(42)
so that K ˙ i∗Kj = κ1K ZiKj + κ2K Zkk δi j
(43)
where δ ij is the Kronecker delta. Introducing the deviatoric quantities
1 3
1 3
∗K K ˙ i∗Kj = ˙ i∗Kj − ˙ kk δi j ZiKj = ZiKj − Zkk δi j
(44)
Eq. (43) can be rewritten in terms of the deviatoric and spherical parts as
3 2
∗K K ˙ i∗Kj = κSK ZiKj ˙ kk = 3κHK Zkk
(45)
with κSK = 2κ1K /3 and κHK = κ1K /3 + κ2K denoting the shear and volumetric coefficients of the transformability tensor. The dissipation rate in inclusion K can now be written as K K 2 K K 2 D˙ K = ZiKj ˙ i∗K j = κS (ZS ) + κH (ZH )
(46)
where
3 2
K (ZSK )2 = ZiKj ZiKj ZHK = Zkk
(47)
Here, for simplicity, the scalars κSK and κHK are considered to be constants. If both κSK and κHK are non-negative, then D˙ K is guaranteed to be non-negative. However, we note that a positive dissipation rate does not require that κSK and κHK are each positive, only that the sum in Eq. (46) is positive. In particular, ∗K at first increases and for STZs the shear transformation can be accompanied by a transient volume change, so that kk then decreases so that at the end of the transformation there is no volume change. Hence, over part of the time for the transformation κHK < 0. In such a case, non-negative dissipation requires
κSK (ZSK )2 ≥ −κHK (ZHK )2 3.1. Plane strain biaxial tension Here, we consider the same plane strain biaxial tension case with a single transforming inclusion as considered in Section 2.2.1. The evolution equation for ∗ is taken to be
1 B¯
˙ ∗ = Zββ Since
D˙ =
(48)
∞ σββ −
E∗ A ˙ ∗ (1 − ν 2 ) incl
it is convenient to take
1 ˙ = B ∗
Zββ
1
Aincl E/[(1 − ν 2 )]
(49)
2 ∞ 1 (1 − ν )σββ = − ∗ B E
(50)
so that B is a drag-like factor having units of time. Define ∗ max =
∞ σββ
[E/(1 − ν 2 )]
(51)
Then Eq. (56) can be written as
1 B
1 B
∗ ˙ ∗ + ∗ = max
(52)
so that the time evolution of ∗ is ∗ ∗ (t ) = max [1 − exp (−t/B )]
(53)
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
7
3.2. Plane strain shear For the plane strain shear case with a single transforming inclusion in Section 2.2.2, the dissipation rate D˙ is given by
∞ D˙ = 2σ12 −
Eγ ∗ Aincl γ˙ ∗ = Z12 γ˙ ∗ 2 (1 − ν 2 )
(54)
and we take
1 B¯
γ˙ ∗ = Z12 or
γ˙ ∗ =
(55)
1 1 Z12 B Aincl [E/2(1 − ν 2 )]
=
∞ σ12 1 −γ∗ B [E/4(1 − ν 2 )]
(56)
With ∗ γmax =
∞ σ12 [E/4(1 − ν 2 )]
(57)
Eq. (56) can be written as
1 B
1 B
∗ γ˙ ∗ + γ ∗ = γmax
(58)
so that ∗ γ ∗ (t ) = γmax [1 − exp (−t/B )]
(59)
In Eq. (59), as in Eq. (53), the maximum transformation strain is approached asymptotically as t → ∞ with characteristic time B. 4. Numerical implementation and results In order to illustrate the use of a kinetic relation guaranteeing non-negative dissipation, a few results are presented for plane strain tension of a rectangular bar with a collection of potentially active shear transformation sites. The numerical formulation follows that in Kondori et al. (2016, 2018, 2019) except for the kinetic relation and will be summarized here for completeness. Geometry changes are neglected and boundary value problems are solved by superposing the analytical solutions for STZ sites and a numerical solution for the image fields that enforce the prescribed boundary conditions. The rectangular bar of width Wm in the x1 −direction and length L in the x2 −direction is subject to the boundary conditions
σ11 = σ12 = 0 on x1 = 0, Wm u2 = U (t ) and u2 = 0 and
(60)
σ12 = 0 on x2 = L
σ12 = 0 on x2 = 0
(61a) (61b)
and u1 = 0 at one point to prevent rigid body motions. At each instant of time the stress and displacement fields are given by Eqs. (1) and (2). The evolution of STZs is determined by (i) the activation of new STZs and (ii) the accumulation of eigenstrain in active STZs. STZ activation takes place when the shear strain energy density φ at site K reaches a specified value φ K , where for an elastically isotropic solid φ is
φ=
σij σij 2G
=
σM2
(62)
3G
2 = 3σ σ /2 with the stress deviator σ = σ − σ δ /3. The value of φ K , for an where the Mises stress, σ M is given by σM ij kk i j ij ij ij STZ site is selected randomly from a Gaussian distribution with a mean value φ nuc and taken to lie within ± 2φ std , where φ std is the standard deviation of φ nuc . The transformation strain for STZ K is taken to be of the form i∗K (t ) = γ ∗K (t ) pKi j with pKi j = (sKi nKj + nKi sKj ) with sK a unit j
vector along a maximum shear stress direction at activation of STZ K and nK a unit normal to sK . The evolution of γ ∗K (t) is given by an expression of the form of Eq. (55) but for multiple STZs as in Eq. (21) and normalized as in Eq. (56) so that
γ˙ ∗K =
K 1 1 Zsn B Aincl [E/4(1 − ν 2 )]
(63)
K = Z K = Z K nK sK . The orientation of sK , nK is chosen so that Z K where B is taken to have the same value for all STZs and Zsn ns sn ij i j
is positive and therefore so is γ˙ ∗K .
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS 8
ARTICLE IN PRESS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Fig. 1. Schematic representation of dynamic STZ nucleation. (a) A potential STZ activation site with radius a in a material under applied boundary conditions. (b) The STZ in part (a) is activated and has created four new potential STZ activation sites at a distance of 2a, along the maximum shear directions. From Kondori et al. (2019).
When an STZ site activates, new potential STZ activation sites can nucleate in the vicinity of the activated site. The nucleated potential STZ activation sites are placed along the maximum shear directions at locations where the transformation induces the maximum shear stress magnitude in the material surrounding the STZ. There are a maximum of four potential STZ activation sites that can nucleate when STZ K activates. However, fewer than four new activation sites may be introduced because new sites are not permitted to overlap an existing activation site or the domain boundary. Fig. 1 from Kondori et al. (2019) illustrates this process. Within this framework, strains greater than the transformation strain can occur by reactivation of a previously activated STZ site. Except for parameters associated with kinetic relation, the parameter values used in the calculations are the same as those for the un-notched bar calculations in Kondori et al. (2019). The bar width is Wm = 64 nm, the bar length is L = 4Wm = 256 nm, the initial number of potential STZ nucleation sites is 2200 (an initial potential activation site density of 1.36 × 1017 m−2 ), the prescribed strain rate is U˙ /L = 0.1 s−1 , Young’s modulus E = 100 GPa, Poisson’s ratio ν = 0.325, nuc = 1.5 GPa and with φ the STZ radius a is 0.6 nm and φnuc = 20MJ/m3 corresponding to σM std = 0.25φnuc . With these pa = σ being the only non-zero components of σ , rameter values E/[4(1 − ν 2 )] = 27.95 GPa and for a stress state with σ12 21 ij ∞ in Eq. (57) with σ nuc = 0.866 GPa gives an estimate of γ ∗ replacing σ12 max of 0.031. 12 K reaching the value zero. However, the kinetic In some cases, the transformation governed by Eq. (63) is completed by Zsn relation in Eq. (63) can give rise to a situation as in Eq. (59) where a complete transformation is only approached asymptotically. To limit the time over which the transformation is allowed to take place, a uniform termination time tN = 4.61B is introduced for all STZs, which is the time for γ ∗ in Eq. (59) to reach 99% of its asymptotic value. The transformation is terminated at t = tNK + tN , where tNK is the time at which the transformation has begun in STZ K. The time delay for reactivation of each STZ site is taken to be tD = ρ tN , where ρ is a specified constant. Thus, in the calculations here all material characteristic times are taken to be proportional to B. To illustrate the sort of responses that are obtained when implementing the kinetic relation in Eq. (63) in a numerical calculation of STZ plasticity, results are shown for three values of the reactivation time tD , with fixed B, and for two values of B, with fixed reactivation time tD . With B fixed at 0.001, Fig. 2 shows results for various values of the reactivation time tD . Fig. 2(a) computed stress-strain curves are shown with σ ave is the average stress on x2 = L (where U(t) is prescribed). As the reactivation time decreases, the load drops more rapidly. Fig. 2(b) shows the evolution of the number of activated STZ sites. Bursts of STZ activation occur with the bursts of STZ activation associated with the drops in σ ave in Fig. 2(a). The strain interval between bursts of activation increases with increasing tD . For tD = 5 tN and 10 tN there is a clear increase in σ ave between the drops, whereas for tD = 3 tN , the activation bursts occur so close together that the value of σ ave does not increase between them. The delay in the bursts of STZ activation for tD = 10 tN leads to strains at which there are more activated STZs for tD = 10 tN than for tD = 5 tN . Fig. 3 shows the distribution of activated STZ sites along with the number of times a site has been activated at values of σ ave ≈ 0.5 GPa. The distribution of activated sites are aligned in bands. With tD = 3 tN and tD = 5 tN a significant number of STZ sites have been activated twice whereas with tD = 10 tN only very few STZ sites have been activated more than once. For all three values of tD , the activated STZ sites are located in deformation bands. The effect of varying B with tD = 5 tN on the stress-strain response and the evolution of activated STZs is shown in Fig. 4. The results for B = 0.001 are repeated from Fig. 2. Since B is the only characteristic time scale in the calculations here (all other characteristic material times are taken to be proportional to B), a governing non-dimensional quantity is (U˙ /L )B. Hence, the results in Fig. 4 can also be regarded as pertaining to an increased U˙ /L with fixed B. Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS
ARTICLE IN PRESS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
9
Fig. 2. The dependence of the stress-strain response and the evolution of STZs on the delay time for reactivation tD with B = 0.001 s. (a) Stress-strain curves. (b) Evolution of the number of activated STZs with overall strain, U/L.
Fig. 3. Distribution of active STZ sites, with n the number of times a site has been activated, for B = 0.001 and for various values of tD . (a) tD = 3 tN s at σave = 0.519. (b) tD = 5 tN s at σave = 0.514. (c) tD = 10 tN s at σave = 0.546.
In Fig. 4(a), with B increased by an order of magnitude from 0.001 s to 0.01 s, the peak stress is somewhat higher, the stress oscillations occur over a larger strain range and σ ave does not fall to zero over the strain range shown. At a given value of strain, U/L, the number of activated STZs in Fig. 4(b) is smaller for the calculation with B = 0.01 than for the calculation with B = 0.001. However, because much larger values of U/L are reached with B = 0.01 before the σ ave vanishes, the number of activated STZs with B = 0.01 eventually exceeds the number activated with B = 0.001 As described here, it is rather straightforward to implement a kinetic relation that guarantees a non-negative dissipation rate in the computational framework of Kondori et al. (2018, 2019). The results in Fig. 2 with B = 0.001 and tD = 5 TN (also shown in Fig. 4) are qualitatively similar to and quantitatively fairly close to those in Fig. 3 of Kondori et al. (2019) for an un-notched bar. 5. Discussion Phenomenological relations, as specified in Homer and Schuh (2009), Kondori et al. (2016), Kondori et al. (2018), for the transformation strain of an Eshelby (1957, 1959) inclusion do not guarantee a non-negative dissipation rate. Indeed, as noted by Kondori et al. (2019), a negative dissipation rate can be obtained in such calculations using typical values for the transformation strain. Section 3 presents a framework for developing fairly general kinetic relations that guarantee a nonnegative dissipation rate within a purely mechanical theory and simple kinetic relations are specified that are analogous Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS 10
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Fig. 4. The dependence of the stress-strain response and the evolution of STZs on B for calculations with tD = 5 tN . (a) Stress-strain curves. (b) Evolution of the number of activated STZs with overall strain, U/L.
to a J − a relation for crack growth or a dislocation velocity-Peach-Koehler force relation for dislocation motion. However, any such kinetic relation limits the magnitude of the transformation strain to be of the order of a stress level divided by an elastic modulus. Much larger strains at a material point are observed in amorphous solids that deform plastically by the STZ mechanism, for example in a shear band. This can only occur within a purely mechanical framework if multiple transformations take place at the same material point. One possibility is multiple activations of the same STZ (as in the implementation in Section 4) and/or by the nucleation of additional STZ sites that overlap the initial site, possibly after sometime delay following completion of the transformation. Atomistic analyses can provide a physical basis for a kinetic relation for transformation strain evolution and atomistic analyses of STZ plasticity have given rise to a deformation mode that can be related to the strain field associated with a transforming Eshelby inclusion. However, the transformation strains so obtained can be much larger than permitted by the requirement of a non-negative dissipation rate. For example, Ju and Atzmon (2014) suggest a STZ transformation shear strain1 of 0.15 based on their analyses of simulated data. Dasgupta et al. (2013) carried out molecular dynamics (MD) analyses of a two dimensional Lennard–Jones solid and found that a shear transformation zone like deformation mode was obtained and the resulting energy associated with the transformations could be fit by the Eshelby inclusion formula with a shear transformation strain of γ ∗ ≈ 0.09. Albert et al. (2016) carried out three dimensional MD calculations and found that the deformation mode was well fit by a transformation strain of the form
= ∗
0.03 0.12 −0.14
0.12 −0.08 −0.07
−0.14 −0.07 0.06
(64)
Clearly, these fits of atomistic calculations to a transformation of an Eshelby (1957, 1959) inclusion give magnitudes of the transformation strain that are much larger than a non-negative dissipation rate (or even non-negative dissipation) allows. Presuming that the atomistic calculations do not violate the requirement of non-negative dissipation rate and the calculations of transformation strain in the atomistic analyses are the appropriate values to relate to the continuum strain measure, there are several possible reasons for the discrepancy between the atomistic analyses and our continuum analysis. For example, it is possible that due to high stresses in the atomistic calculations linear elasticity is not an appropriate constitutive description or that surface energy plays a significant role in the atomistic analyses. In addition, it is worth noting that the predicted transformation strain magnitude with a non-negative dissipation rate increases with increasing stress level (in Sections 2.2.1 and 2.2.2 the limiting strain is of the order of a stress level divided by an elastic modulus). However, transformation strains associated with such large stresses would not be relevant for mesoscale modeling of plastic deformation at much lower stress levels. Another possibility is that, within the continuum context, a purely mechanical theory does not provide an adequate mesoscale framework. A broader thermodynamic framework may be needed. 1 ∗ ∗ , in others it is identified with 212 and in others which of these In some reports, as here, the transformation shear strain is identified with 12 definitions is used is not explicitly specified.
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
ARTICLE IN PRESS
JID: MPS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
11
Within the context of continuum thermodynamics and for a spatially uniform and time independent density (as appropriate for our small strain formulation), the Clausius–Duhem inequality can be written as
˙ − W˙ + s˙ T −
qi ∂ T /∂ xi ≥0 T
(65)
where T is the absolute temperature, qi are the components of the heat flux and s is the entropy. Assuming that qi and/or ∂ T/∂ xi vanish, Eq. (65) simplifies to
˙ ≥0 W˙ + s˙ T −
(66)
There are several possibilities for allowing larger transformation strains. One possibility is that rather than the material response during transformation being an elastic process, it is an elastic-plastic (or more precisely an inelastic process of p some kind), so that the total strain rate during transformation is ˙ i j = ˙ iej + ˙ i j and within a transforming inclusion
˙ + σi j ˙ p W˙ = σi j ˙ i j = σi j ˙ iej + σi j ˙ ipj = ij
(67)
then Eq. (66) becomes
σi j ˙ ipj + s˙ T ≥ 0
(68)
However, it is unlikely that inelastic deformation in the inclusion permits larger transformation strains presuming that the material response is elastic prior to transformation and the transformation reduces the stress level in the inclusion so that the process involves elastic unloading. Furthermore, Vasoya et al. (2019) showed that a reduction in stiffness in the inclusion relative to that of the material outside the inclusion reduces the transformation strain that can occur with nonnegative dissipation. p If σi j ˙ i j = 0 in Eq. (68), what is left is the entropic contribution and, since the rate of change of entropy is always positive, s˙ T gives a positive contribution and thus can contribute to allowing a larger transformation strain for an Eshelby (1957, 1959) inclusion. One source of entropy production is that arising from atomic thermal vibrations in the inclusion. Another source, that could possibly give a contribution, arises from an increase in configurational temperature and entropy due to the structural rearrangement in the transforming inclusion. Indeed, atomistic analyses and experiment indicate that entropy can play a role in dislocation nucleation at finite temperatures, see for example, Chen et al. (2015); Nguyen et al. (2011); Ryu et al. (2011). It is also worth noting that a formulation for using statistical mechanics to relate the entropy of an atomic system to mesoscale quantities for dislocation dynamics has been presented by Acharya (2011). To give an estimate of the possible effect of entropy production, we presume the dominant entropy effect is configurational and take the configurational temperature to be proportional to the transformation strain magnitude and the configurational entropy to be proportional to the configurational temperature, and focus on the plane strain shear case of Sections 2.2.2 and 3.2. With these restrictions, the Clausius-Duhem inequality can be written as ∗ + CD γ ∗ − γ ∗ ]γ˙ ∗ ≥ 0 [γmax
where the mechanical limit strain
γth∗ −lim =
γ
∗ max
1 − CD
(69) ∗ γmax
is defined in Eq. (57) and the thermodynamic limit on the transformation strain is
(70)
∗ so if γmax = 0.03, then CD = 2/3 to CD = 3/4 would allow transformation strain magnitudes in the range 0.09 to 0.12 with non-negative dissipation.
6. Conclusions 1. The analysis of Vasoya et al. (2019) has been modified to give a general expression for the dissipation rate and extended to account for multiple transforming Eshelby (1957, 1959) inclusions. 2. The requirement of non-negative dissipation rate is more stringent than the requirement of non-negative dissipation. • For a single transforming inclusion in an infinite solid with a remote uniform stress, the requirement of non-negative dissipation rate gives a limiting transformation strain that is a factor of two smaller than is given by the requirement of non-negative dissipation. 3. A framework for specifying the evolution of the transformation strain in terms of the associated configurational force was proposed. Such kinetic relations guarantee a non-negative dissipation rate. • Implementing a kinetic relation that guarantees a non-negative dissipation rate into a computational framework is rather straightforward. In a simple quasi-static initial/boundary value problem, and with an appropriate choice of parameter values, results are obtained that are qualitatively consistent with previous results. The implementation here was in the computational framework of Kondori et al. (2018, 2019) but such kinetic relations can be incorporated into other computational frameworks. 4. Atomistic analyses indicate that the magnitude of STZ transformation strains can exceed the value permitted by the requirement of non-negative dissipation rate. There are a variety of possibilities for this discrepancy, one of which is that a purely mechanical mesoscale framework may not be adequate for modeling STZ plasticity. Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699
JID: MPS 12
ARTICLE IN PRESS
[m3Gsc;August 28, 2019;13:32]
M. Vasoya, B. Kondori and A.A. Benzerga et al. / Journal of the Mechanics and Physics of Solids xxx (xxxx) xxx
Acknowledgments We gratefully acknowledge financial support from the Weizmann Institute of Science through the Weizmann-Texas A&M program. References Acharya, A., 2011. Microcanonical entropy and mesoscale dislocation mechanics and plasticity. J. Elast. 104, 23–44. Albert, T., Tanguy, A., Boioli, F., Rodney, D., 2016. Mapping between atomistic simulations and eshelby inclusions in the shear deformation of an amorphous silicon model. Phys. Rev. E 93, 053002. Andric, P., Yin, B., Curtin, W.A., 2019. Stress-dependence of generalized stacking fault energies. J. Mech. Phys. Solids 122, 262–279. Barnett, D.M., Cai, W., 2018. Properties of the eshelby tensor and existence of the equivalent ellipsoidal inclusion solution. J. Mech. Phys. Solids 121, 71–80. Budrikis, Z., Castellanos, D.F., Sandfeld, S., Zaiser, M., Zapperi, S., 2017. Universal features of amorphous plasticity. Nature Commun. 8, 15928. Bulatov, V.V., Argon, A.S., 1994. A stochastic model for continuum elasto-plastic behavior. i. numerical approach and strain localization. Modell. Simul. Mater. Sci. Engin. 2, 167–184. Chen, L.Y., He, M., Shin, J., Richter, G., Gianola, D.S., 2015. Measuring surface dislocation nucleation in defect-scarce nanostructures. Nature Mater. 14, 707–713. Dasgupta, R., Hentschel, H.G.E., Procaccia, I., 2013. Yield strain in shear banding amorphous solids. Phys. Rev. E 87, 022810. Eshelby, J.D., 1957. The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proc. R. Soc. Lond. A 241, 376–396. Eshelby, J.D., 1959. The elastic field outside an ellipsoidal inclusion. Proc. R. Soc. Lond. A 252, 561–569. Gavazza, S.D., 1977. Forces on pure inclusions and Somigliana dislocations. Scr. Metall. 11, 979–981. Gavazza, S.D., Barnett, D.M., 2018. The pointwise Eshelby force on the interface between a transformed inclusion and its surrounding matrix. Math. Mech. Solid 23, 233–239. Hill, R., 1965. Continuum micro-mechanics of elastoplastic polycrystals. J. Mech. Phys. Solids 13, 89–101. Hirth, J.P., Barnett, D., Hoagland, R.G., 2017. Solute atmospheres at dislocations. Acta Mater. 131, 574–593. Homer, E.R., Schuh, C.A., 2009. Mesoscale modeling of amorphous metals by shear transformation zone dynamics. Acta Mater. 57, 2823–2833. Hufnagel, T.C., Schuh, C.A., Falk, M.L., 2016. Deformation of metallic glasses: recent developments in theory, simulations, and experiments. Acta Mater. 109, 375–393. Ju, J.D., Atzmon, M., 2014. A comprehensive atomistic analysis of the experimental dynamic-mechanical response of a metallic glass. Acta Mater. 74, 183–188. Kondori, B., Benzerga, A.A., Needleman, A., 2016. Discrete shear transformation zone plasticity. Extreme Mech. Lett. 9, 21–29. Kondori, B., Benzerga, A.A., Needleman, A., 2018. Discrete shear-transformation-zone plasticity modeling of notched bars. J. Mech. Phys. Solids 111, 18–42. Kondori, B., Benzerga, A.A., Needleman, A., 2019. Errata to: discrete shear-transformation-zone plasticity modeling of notched bars. J. Mech. Phys. Solids 127, 151–153. Markenscoff, X., 2010. Evolution equation of moving defects: dislocations and inclusions. Int. J. Fract. 166, 35–40. Markenscoff, X., Ni, L., 2010. The energy-release rate and “self-force” of dynamically expanding spherical and plane inclusion boundaries with dilatational eigenstrain. J. Mech. Phys. Solids 58, 1–11. McMeeking, R.M., Evans, A.G., 1982. Mechanics of transformation- toughening in brittle materials. J. Am. Ceram. Soc. 65, 242–246. Mura, T., 1982. Micromechanics of Defects in Solids. Martinus Nijhoff, The Hague. Nguyen, L.D., Baker, K.L., Warner, D.H., 2011. Atomistic predictions of dislocation nucleation with transition state theory. Phys. Rev. B. 84, 024118. Rice, J.R., 1975. Continuum Mechanics and Thermodynamics of Plasticity in Relation to Microscale Deformation Mechanisms. In: Argon, A.S. (Ed.), Chapter 2 of Constitutive Equations in Plasticity. M.I.T. Press, pp. 23–79. Ryu, S., Kang, K., Cai, W., 2011. Entropic effect on the rate of dislocation nucleation. Proc. Nat. Acad. Sci. 108, 5174–5178. Sandfeld, S., Budrikis, Z., Zapperi, S., Castellanos, D.F., 2015. Avalanches, loading and finite size effects in 2d amorphous plasticity: results from a finite element model. J. Stat. Mech. P02011. Vasoya, M., Kondori, B., Benzerga, A.A., Needleman, A., 2019. Limits on transformation strains for non-negative dissipation. J. Appl. Mech. 86, 051005.
Please cite this article as: M. Vasoya, B. Kondori and A.A. Benzerga et al., Energy dissipation rate and kinetic relations for Eshelby transformations, Journal of the Mechanics and Physics of Solids, https://doi.org/10.1016/j.jmps.2019.103699