- Email: [email protected]
The eigenstrain method for small defects in a lattice
J. Phys. Chem. Solids Vol. 52, No. 8. pp. 1019-1030, Printed in Great Britain.
1991
THE EIGENSTRAIN
0022-3697191 $3.00 + 0.00 0 1991 Pergamon Press plc
METHOD FOR SMALL DEFECTS IN A LATTICE
K. C. KING and T. MURA Theoretical and Applied Mechanics, Northwestern University, IL 60201, U.S.A. (Received 19 November 1990; accepted 21 March 1991)
Abstract-The eigenstrain (eigendistortion) method is used in the lattice theory to analyze the disturbance caused by the defects in a lattice. The equilibrium equations for the atoms are derived through the use of lattice theory, and the harmonic assumption for the potential energy is used for the derivation. The eigenstrain (eigendistortion) is used to model the defects, which could have a complicated shape. The results are compared with the results obtained from the elasticity theory. A simple cubic model is used to model the structure of the lattice. The Green function is obtained from the equilibrium equation by using the discretized Fourier transformation method. Therefore, the disturbance caused by the defects in the infinite lattice can be analyzed by the summation of the Green function with the strength of the eigendistortion. The problem for the defects in the thin film is solved. The thin film is coated on a substratum and bound by a free surface on the other end. The substratum is assumed to be rigid. The result is compared with the result for the same defect in the infinite lattice. The comparison indicates that the effect of the boundaries of the thin film is important. Keywords: Lattice, thin film, eigenstrain, defects, prismatic dislocation loop, inclusion.
1. INTRODUCTION
The precise description of the distortion of the surrounding atoms near defects in solid crystals is important for the calculation of the energy of defects [l], the discussion of X-ray diffraction effects [2, 31 and also the change in electrical conductivity of metals [4]. Such calculations have been done by Kanzaki [5] using the concept of the Green function, which is to calculate the displacement of the atom L caused by a unit body force acting on the atom P. Therefore, the displacement of a ‘normal’ atom caused by a substitutional or interstitial defect in a perfect lattice can be obtained by solving the equations formulated by the use of the Green function [S]. This method has been improved by Tewary [6], who introduced the ‘defect Green function’ and this defect Green function can be applied to the lattice structure with the crack surface inside it. Several works have been done by Thomson et al. using the defect Green function, such as, lattice distortion due to gas interstitial in bee metals [7], three-dimensional crack with a kink [8] and lattice theory of fracture and crack creep [9]. However, the defect problem like dislocation cannot be solved from these methods in the above, since the dislocation is prescribed by the displacement jump (Burgers vector) on the dislocation plane instead of the force acting on each atom. Therefore, eigenstrain is used in the lattice theory to overcome this kind of difficulty. Eigenstrain is a general name to describe misfit strain, thermal expansion strain,
phase transformation strain, plastic strain, etc. For example, the misfit induced by an interstitial atom can be estimated by a fundamental electronic theory [lo, 111; the misfit for the substitutional atom can be determined by the size difference between the atoms; the misfit for the dislocation is the Burgers vector [12]. The idea here is to use this misfit to calculate the distortion field without the use of the Green function in the above. This is so called ‘eigenstrain (eigendistortion) method’ which is widely used in continuum mechanics. The application of this idea in lattice theory was first introduced by Mura to calculate the moving dislocation [13], and it was again used to calculate the octahedral defect in a bee lattice and the growth of nitrogen defect cluster by Sato et al. In this work, we further improved the eigenstrain method in lattice theory by introducing the ‘eigenstrain Green function’. This eigenstrain Green function is used to calculate the distortion of atom G caused by a unit eigenstrain (misfit) between two atoms !’ and e*, therefore, the total distortion field can be obtained by the summation of the eigenstrain Green function with the strength of eigenstrains. This method can solve the dislocation defects as well as the interstitial, substitutional defects. Furthermore, the defects with complicated shape can be systemically solved with the help of the eigenstrain Green function, and this is the benefit of this method. However, in order to simplify our mathematical task considerably, a simple cubic lattice is chosen to illustrate this method although no significant difficulty can be seen here as the other lattice
1019
K. C. KING and T. MURA
1020
structures such as bee or fee are chosen. The importance of this method is that it can be further extended to model cracks and the critical criterion by calculating the equivalent eigenstrains; it is under preparation now. We also investigate the free boundary effects for the defects in lattice. Traction free boundary effects are used for the free boundary in this analysis. The equilibrium equation is formulated for a finite layer lattice domain (thin film). The thin film is coated on a substrate and the substrate is assumed to be rigid, which implies a rigid boundary condition on the other side of the thin film. The same model for the atomic structure is used for the thin film and no special surface lattice structure is assumed. In the following, lattice theory [14] and the harmonic assumption of the potential energy [ 151are used to formulate the equilibrium equations for the atoms [13]. A simple cubic lattice with the consideration of nearest and second nearest atoms interaction is assumed for the mathematical formulation. A small defect domain with cubical shape and ellipsoidal shape is calculated; the prismatic dislocation loop is also calculated in this work. The results are compared with the results from elasticity calculations [l&20]. The equations for the defects in the thin film are then formulated. The prismatic dislocation loop defect is used for the calculation. The result is compared with the result for the same defect in the infinite lattice. The comparison indicates that the effect of the boundaries of the thin film is significant.
2. THEORETICAL
DERIVATION
2.1. Equation of motion The harmonic approximation of the potential energy is introduced first in this section [ 151.Then the equation of motion will be derived by using the Hamilton principle. Finally, the equilibrium equation can be obtained from the equation of motion [13]. Let Y(LG’) be the potential energy of a pair of atoms L and e’, where G and d’ are the label of atoms at perfect lattice sites and !!’ is the vector measured from atom G to atom L’. For a perfect lattice structure, the total potential energy U of the system can be determined by the summation of the contributions of all distinct pairs as U =;I
Y(ee’). 1.f’
VY (et’)
a=~~,Y(pl.)+~~(u(C)-u(C’)) (.I’ 1.J + ; c [(u(d) - u(P)) f,c?
V12Y(68' ) + O(u’).
(3)
The gradient of potential energy, VY(c/‘), is known as the interactive force between the two atoms / and L’. Therefore, the summation of interaction force, Zr VY(!J’), is zero, since there is no net force on any atom in equilibrium condition. It follows that the linear term in eqn (3) is zero. Also, the higher order terms are neglected since small displacement is assumed in this work. Then the total potential energy can be written as follows: U = Ueq+ Uham
(4)
u’q =; 1, Y(LG’) /.I
(5)
Uhitrm= ; E, [(u(L) - u(L’)) v]*Y(ee’). 1.1
(6)
where
The equilibrium potential energy U'q is the total potential energy for the system without any defects and it can be set to zero as a reference state because the absolute potential is not important in this work. Then, the only term left for the total potential energy is Uharm. This is called the harmonic approximation. The harmonic potential energy can be written in a more general form as u=u
harm= ;;,
@,(dd’)u,(d)u,(d’),
(7)
where @,(ee’) = 6,,, 1 Y.,,(dd”) - Y.,,(&/‘). I”
(8)
The double indices appearing in the equations follow the summation convention from 1 to 3, and the commas in the subscripts denote the derivatives. The total kinetic energy T of the system is
(1)
For a distorted lattice structure caused by dislocations or inclusions, this potential energy can be calculated by the summation of the potential energy of each pair of atoms measured at the equilibrium position as u = ; 1 Y (ee’ + u(d) - I@‘)), !,L’
where u(8) is the displacement vector of atom / measured from its original lattice site. For small displacements, the total potential energy U can be expressed by Taylor’s series as
(2)
T = ;
1 Mti, (f)&(f), I
(9)
where M is the mass for each atom. The equation of motion for each atom e can be derived by using the Hamilton principle S[{,:(U-T)df]=O
(10)
Eigenstrain method for small defects in a lattice
1021
The equilibrium equation can be obtained from eqn (19) by neglecting the dynamic term as Miii(e)
=
-1 @ij(ee’)uj(e’),
(11)
f’
where 6 is the symbol for variation. It can be seen from eqn (11) that @ij(!~‘) are the force constants between the atoms L’and c”, and these force constants have the following properties ~~ij(ec’)=~~ij(~L’)=o
(12)
@,(Ld’) = Gji(ee’) = @,(e’t),
(13)
which come from the rigid body motion and the symmetric conditions [13, IS]. Then eqn (11) can be written as
1 o,(ee’)u,(e’) C’
where s”(/*L) is the Kronecker delta and it has the following properties: 8(/*k) = 1 as e* = / otherwise $((e*e) = 0. For a single unit eigendistortion I$(/*/‘) applied between atoms d* and L”, eqn (20) can be rewritten in the following ;@,(ee’)u$‘)
M&(d) = -1
@ij(Ll’)[Uj(f’)
-
= ~,(e*d’)z,:.(e*e’)x,(e~e*)~(‘(e*e),
(14)
Uj(f)].
C’
We assume the difference of the displacement vectors between e and !’ is U,(C) - Uj(/) = /I,,(ee’)x,(/‘c), where /Iy(L’e’) are the elastic distortions. tion A’,(~‘z!) is defined as
= 1 ~ij(ee’)8t(ee’)x,(c’e) I’
(1% The nota-
&(e’e) = X,(/‘) -x,(e),
(21) where the superscript E denotes the displacement caused by a unit eigendistortion. The right-hand term in eqn (21) can be interpreted as a fictitious force which is equivalent to the effective force caused by the eigendistortions. Thus, as a unit body force I,(/*) in the i-direction is applied on the atom c!*, eqn (20) can be expressed as follows
(16) ; @,,(e/‘,UF(e’, = &(d*)s(e*d),
where X,(e) is the position coordinate component k of the atom C! in a perfect lattice site, then the equation of motion (14) can be written as MYi,
= -c
@ij(de’)pkj(de’)X#‘e).
(17)
(’ However, if material is involved in any inelastic deformation, then eqn (17) is no longer valid because /Ikj(eL”) are not compatible. It is necessary to reformulate eqn (15) as u,(P) - u,(C) = [fikj(LL’) + B$(d/‘)]X#‘L),
(18)
where p$.(/L”) are the eigendistortions. The physical meaning of the eigendistortions can be interpreted as the stress-free transformation strains. For example, the differences of lattice parameters induced by the phase transformation are the eigendistortions. The definition here is similar to the definition in Micromechanics of Defects in Solids [12] that the total strain is the sum of elastic strain and eigenstrain. Substituting eqn (18) into eqn (17) and using the properties of the force constants in eqn (12), then the equations of motion become M&(e) = -c
@,(ee’)uj(e’) C’ + c ~ij(ee’)B~(e~‘)x,(~‘~). C’
(19)
(22)
where the superscript G denotes the displacement caused by a unit body force. Equations (21) and (22) are the basic equations for the following derivation. One is used for the case of distributed eigendistortions and the other one is used for the case of distributed body forces. The Green function for a unit eigendistortion can be derived from eqn (21) and the Green function for a unit body force can be derived from eqn (22). It must be pointed out here the reason to obtain two Green functions in this work. The eigenstrains (eigendistortions) used in the continuum mechanics are continuous functions and they are defined continuously in the defect domain. Therefore, the effect of the eigenstrains are just like body forces and the solution can be obtained by the integration of the Green function (obtained from unit body force) with the strength of the eigenstrain. However, the eigendistortions used in the lattice theory are not continuous functions and they are defined between any two atoms inside the defect domain. It follows that the fictitious force caused by the eigendistortion will act on two atoms simultaneously. Therefore a new Green function has to be defined to calculate the effect of the fictitious forces caused by the eigendistortions. The details of the derivation and the application of these Green functions will be discussed after we set up the model for the atomic structure.
1022
K. C. KING and T. MURA
2.2. Model for atomic structure The force constants @,,(e/‘) can be obtained from the potential energy Y(ee’) through the use of eqn (8). However, it is not convenient for the present work. A relation between force constants and elastic constants is more important as we want to compare the results between lattice theory calculation and the elasticity calculation. Thus, the simple cubic model proposed by Mindlin [21] is used. These relation equations will be shown at the end of this section, A simple cubic lattice with nearest and secondnearest atoms interaction is the model chosen for the present work. In this model each atom is interactive with 18 surrounding atoms, which are six nearest atoms and 12 second-nearest atoms. Let tl be the central force constant (i.e. force per unit relative displacement) between any nearest atoms, p be the central force constant between any second-nearest atoms and y be the angular force constant between any three, successive, non-collinear atoms (see Fig. 1). The variable b in Fig. 1 is the lattice parameter. The equation of motion for this model is given as M@l.” = #;+
@,,(f,m,n;fIfr l,m,n)= %,U,m,n;hm
-tl
+ 1,n)
=@,,(I,m,n;I,m,n
*l)=
-47
@,,(I,m,n;I+l,mfl,n) =@,,(i,m,n;t+
I,m,n + 1)
=@,,(l,m,n;I-
I,m rf 1,n)
=@,,(Z,m,n;I-l,m,n+l)=-j3
O,,(l,m,n;I,m
+ l,n)=
f-h(bm,n;l,m,n
-u
It 1)
=&(l,m,n;I+l,m,n)= ~~~(l,m~n;l,m
-47
+ I,n -t 1)
=@,,(I,m,n;I,m - 1,n + 1) =@,,(l,m,n;l+
l,m ILt:1,n)
=@,,(I,m,n;I-
1,m Ifr l,n)=
@33(1,m,n;i,m,n _t I)= -u
I.m,n+ U:- 1.m.n _ 2u;m.s)
/+i.m+in + $ l.m- 18+ *i+ I.m- 1.n + a@, +U~-l.m+l,n+ll~+I,nin~f+~~-l.m.~-l
@,,(l,m,n;Irt
1,m,n)=@33(I,m,n;I,m
_t I,n)
= -4y %(I,m,n;f+l,m,n
+U;+i.m.n-l +U~-Im.?l+I _ &@.‘7)
t- 1)
+ (p + Y)(U;+ i.m+1.N+ &- I.ffl_1.n
=@33(1,m,n;t-l,m,n
_ U/+ I.rn- I.”_ U/-. 1.m-tI.n
=@&(l,m,n;l,m + 1,n Lt: I)
2
-p
* 1)
2
=@,,(I,m,n;l,m
+~:+I,m.n+I+U~--I.m.n-I+~I+I.m.n+I
- 1,n + 1)= -/3
3 +
&
l,m.n - 1 _
+4~(U~m+l.n
+24
yl
Up
+
-
@z-
1.m.n - i _
La+
@:-
@12(lrm,n;[+ 1,m + 1,n)
1.m.n + I)
=@,z(l,m,n;f-
U:n.n+l
1_ 4u$m‘“)*
(23)
Here the subscript denotes the component of displacement vector and the superscripts denote the position coordinates of the atom in a perfect lattice site. The remaining two equations of motion can be obtained by cyclical permutation of subscripts and superscript. If we drop the term with eigendistortions @$$‘8“) in eqn (19), then we can determine the force constants by comparing eqns (19) and (23) as follows @,,(1,m, n; 1,m, n) = @22tl, m, n; 1,m, n>
=@,,(f,m,n;i,m,n)=2or+8~+16y
Fig. 1. Schematic illustration of the spring constants a, & . . . . .. y used in the simple culxc model.
I,m + I,n)= -(b +y)
@13(lrm,n;i+ Lm,n t- I) =@,,(Z,m,n;t-l,m,n+l)=-(/3+y) @,3(I,m,n;I,m+l,nf
1)
=@23(Irm,n;l,m - l,n + 1)= -(p
+y)
clrI2 = 6, @13= @3r
(24)
%3 = Q13,.
The relations between ~1,p and y in this model and the elastic constants were derived by Mindlin [21] Ithrough finite difference method as follows
G4=;
(B +2Yh
Eigenstrain method for small defects in a lattice where C, 1, C,*, C, are the Voigt elastic constants for cubic symmetric material. This relation (25) is used to compare the solution obtained from lattice theory calculation with that obtained from the elasticity calculation. The solution for defects embedded in an infinite lattice and thin layer are studied in the following sections.
1023
Solving eqn (28) for QF, we can obtain
l Nij(5)zi(e*) ew[-iW, + 452+ rtdl, (30) cY=$ o(r) where D and Nij are the determinant and cofactors of the matrix K. The displacement u?(e) can be obtained by substituting Is,” into eqn (26) as follows
2.3. Green function Two Green functions, which are explained in the previous section, are derived in this section and they are used to calculate the displacement of the atoms caused by the defects. The Green functions obtained in this section are valid for an infinite cubic symmetric lattice. The discretized Fourier transformation method is used in the present work [ 131.The discretized form of the displacement is
UP)
1 = (271)3 -*
x exp{W -P)<, + Cm- 4% (31) From the equations above, the Green function for a unit body force can be solved as
uj(e) =
1 G,(ed*) = m
x exp[i(E, + 4
+ G)l dt, G dt,,
--II
(26)
x exp{W -PK + (m - q)L
where b is the lattice parameter, and the discretized form for the Kronecker delta is
x ew{W -PK + (m - q)t2
where (I, m, n) and (p, q, r) are the position coordinates of the atoms L+and e* at the perfect lattice site, respectively. Substituting eqn (26) and eqn (27) into eqn (22), then eqn (22) becomes K/jr77 = i Z#*)exp[-
i(P5, + q
(32) The physical meaning of the Green function G,,(LL*) is the displacement of atom L(Z,m, n) in the j direction when the unit force in the i direction is applied on the atom /*(p, q, r). However this Green function is inadequate for the eigenstrain method. A similar function Ej&!*L’) (Green function for a unit eigendistortion) is derived for practical purposes. The physical meaning of this function E,,(/!*!‘) is to describe the displacement of atom L’(Z,m, n) in the j-direction when the unit eigendistortion Zg(t*/‘) is applied between atom /*(p, q, r) and one of its surrounding atoms 8’. This function Elks(L/*d’) can be derived in the same way as the Green function is derived. First, substituting eqn (26) and eqn (27) into eqn (21), eqn (21) becomes
(28)
where the coefficients Ku are x exp] - Z(P5, + q5r +
K,, = Za(cos 5, - 1) + 4fi(cos 5, cos c2 + cos 5, cos {, - 2) + 8y (cos 52 + cos 5, - 2) Kz2 = 2cr(cos & - 1) + 4B(cos r, cos
&)I. (33)
Then, solve for 17,”and substitute ri,! into eqn (26) to obtain the displacement u,?(e) (which is caused by the effect of unit eigendistortion) as
+ cos {r cos & - 2) + 8y(cos r, + cos {, - 2) K,, = 2a(cos r, - 1) + 4/?(cos 5, cos &
+ cos & cos C3- 2) + 8y (cos {, + cos & - 2) KL2= K2, = -4(fi
1 Yw) = Qn jjj @,(e*e’)x,(re*) -n
+ y)sine, sin & X
K13 = Kj2 = -4(/? + y)sin & sin & K,, = K,, = -4(/I
+ y)sin 5, sin &.
(2%
Ni’(5gt(;*“)
exp{i[(Z-p)r,
+ (m -q)& (34)
K. C. KING and T. MURA
1024
From this equation, the Green function Ejks(dL*d’) is defined as z
1 = @-#
E,,(ee*e’)
N,.(5) @,,(e*e’)x,(c’e*) i D(5)
x ev{W -_pK + (m - q)ti + (n - rkll dt!, di%d53.
(35)
This function has the following property ; E&(ee*e’)
= 0
Ejks(/d*/‘) = 0 (for L’ $ A),
(36) (37)
where the domain A represents the 18 atoms surrounding the atom L*. Equation (36) can be proved through the use of eqn (13) and eqn (16). Equation (37) comes from the fact that @is(/*e’) = 0 in eqn (35) for d’# A. These Green functions are valid for the infinite cubic symmetric lattice problems. The application of EjkJ(,,(G/*L’)is investigated in the following section. 2.4. Small defects in the infinite lattice With the help of the Green functions derived in the previous section, it is easier to treat defect problems in the lattice. Since eigendistortions can simulate defects well, this work is emphasized on the use of the Green function for a unit eigendistortion. The method of using the other Green function (the one for unit body force) is the same and does not require further explanation in this work. For a given defect domain R prescribed by the eigendistortions, the displacements of the atoms can be obtained by the summation of Ejks(Ld*d’) over R as a,(l)=
z
C Pk:(~*~‘)E,,,(~b*~‘). /’ /*En /‘E(RnA)
problem. Thus, a different method is used in this section. The thin film is assumed to be coated on a substrate and the substrate is considered to be very rigid when compared to the thin film. It follows that the displacements at the interface are zero in this analysis. The other boundary is assumed to be a traction free boundary. The displacement field caused by the defects inside the thin film is solved in the following. Let the x1 direction be the direction normal to the plane of the thin film. The discretized Fourier transform is applied to the x, and xj directions. Therefore, the discretized form of the displacement is
(38)
The summation domain n depends on the position of the atom !*. Only the atoms inside R are necessary for the summation because the eigendistortions b$$(e*/‘) are zero outside the domain R. Furthermore, if the eigendistortions are constant inside the domain R, the summation domain R (for the atom d*) can be changed into the domain aR, which is the boundary of R. This comes from the property of the function Ejks(/t*ta’) expressed in eqn (36). Thus, for the case of a constant eigendistortion distribution, only the eigendistortions on the inner surface of the defects are necessary for the computation. 2.5. Small defects in a thin$lm As mentioned before, the Green functions derived in the previous section are no longer valid for this
x exdi(mt;, + n&)1dt2 G
(3%
and for the Kronecker delta is
b(e*e) E 6,,6,6,” =
n ss
&s,, = --I
--II
x expIi[(m - q)& + (n - r)&l} d&d&,
(40)
where (I, m, n) and (p, q, r) are the position coordinates of the atoms L and L* at the perfect lattice site, respectively. Equation (39) can be integrated to obtain the solution of displacement field after tij(f, &, t3) is solved. The unknown variables 1?,(1,&, 5,) can be obtained by solving a system of linear equations, which come from the equilibrium equations of each layer. There are three equations (three components of the displacements) for each layer and these linear equations for each layer have to be solved simultaneously. Since each layer is only linked to its neighbor layers (upper and lower), a banded matrix will be formed after assembling all the equations together. The banded matrix can save computer time as well as working space. The system of equations is separated into three parts according to the position of the layers as: (1) equations for the layers inside the thin film (case 1); (2) boundary equations for the free boundary (case 2) and (3) boundary equations for the rigid boundary (case 3). Case 1. The layers 1 are assumed to be located inside the thin layer. The linear equations for these layers inside the thin film can be obtained by substituting eqns (39) and (40) into eqn (20) as
where
-(/+I), ug+ 1) ++ I)), (42) u1
Eigenstrain method for small defects in a lattice the notation
tij in eqn (42) denotes iij(l, &., &), and
C1(U~Lm.~ + U$L.m.n
1025 _
2U;L.m.n)
p(ufLmt
+
Ia+
I
+u~L.m-I.n-I+U~L.m+In-I+UfL.m-I.ntl +
u$h!-
+
(fj
IAm+ +
1.~ +
u$(L-
I),m-
lflf
I+ Up-
r)(u;L*“+
IJI _
,5u~~m~)
In-
I
_~~L.mtlr-I_~~L.m-l,n+t)
x; = 6,JC c @,(C*e’)/?$(t!*e’)X,(L’L*) R t’ x ev[--q(b
+ d31
py
+Y)(U:(L-I).m-I,n_Ulf(L-I),m+I.n)
.+. +
(u f Lmn
t
I +
u ; L,m.n -
I
+ uf,“- Ihn _ gr4;L.m.n) TY@fL.“fI.” _ ufL.W Ln)= 0 q(r&,m.n~ 1+ wjiL.m.n+ 1_ &pwl) (43)
+~(UfL~mtl.ntl+~~L.m-l~n-l fur
H,, = H,, = 2p cos & + 4y
+L.m+I.n-I+~~L.m-1.~+l+~~,L-I~.mr+l
+U:,L-I),m.n-I_6U:L,m.n
H,, = H,, = 28 cos & + 4y H,, = Hz7 = -H,,
= -Hz,
H,, = H3, = -H,,
= - H3, = i[2(B + y)sin c,]
+(p
= i[2@ f y)sin &]
+y)(U~L.nt+I.ntI+U~L,m-~."-I _U~L.m-lCti
_ti;L,m+l.n-l T(fl
H,, = H,? = u + 2~(cos t2 + cos &) Hz6 = H,5 = -4(p
H,,=
)
+Y)(*~(L-lf.mn+l
) _M~,L-if.nr,n-l)
+4QY(U~L.m+I.n+U~L.m-I.n+~~,L-I~.m.n
+ y)sin &sin & _ 3UfL.m.n)-+ y(UfL.m.n+
I _ UjtL.m.n-
I) =
0,
(45)
-(2a: +SP + 16y)+8y(cos~,+cos~,)
&=-(2or+8p+16y) + (2a cos & + 48 cos & cos & + 8y cos &) & = -(2a
f 88 + 16~)
+ (2~7cos & + 4pcos t2 cos & + 8y cos &)
where only the equation for I= +L is used for the present work because this layer is the boundary layer with traction free assumption. The ~uilib~um equations for the free boundary are obtained by substituting eqn (39) into eqn (45) as
WtLl,j x6)
H,, = H,, = H,, = Hz4 = H,, = H,2 = Hj4 = H3* = 0.
(44
Case 2. The layer +L is assumed to be a traction free boundary. The boundary equations come from the equilibrium equations for the atoms on the free boundary. For each boundary atom, the equilibrium equations are derived through the summation of the forces which come from the interaction of its neighborhood atoms inside the thin film. The summation should be zero because of the traction free condition. The general equations for the free boundary at I= t_L are expressed in the following
cL(u ? (f - I)JW_ U~L,“L)
’ [u+Ll(6x
[U+L]r = f&j”- ‘1,&-I),
(46)
@-“,
rjf, j& n$> (47)
H:,L = 2p cos & + 4y H&L = 2fl cos & + 4y HAL = - HhL = i[2y sin &] H:SL = - HGL = i[2y sin <,I
HhL = HAL = -i[Z(B + y)sin ~$1 HAL = HGL = --i[2(/? + g)sin &] HAL = d + 2j?(cos & + cos &) HgL = HGL = -4(p
+U~,L-I~,m--I,n+U~~L-I~.mtI.n_4~~L.m.n)
HhL=
T(jj +y)(U~(L-I).mtI.n_U~,L-I),m-I..
HGL = -(2cr + 68 + 12~)
+Uf,L-If.m.n+l_~f,L-l~.m.n-I
H$
+U~L.m.n-l_4ufL.m,n)~y~U~L.m-ln_U~L.m+1.n
t I) = 0
-(c( +4j
+ y)sin C2sin &
+8y)+4y(cos&+cos’r,)
+ (2or cos r* + 48 cos tz cos & f 8y cos &)
I
f2y(ufL‘“+‘“+u~tJn-In+u~t.m,nt1
UfL.m.n
0,
where
+~(UTI~--I).m+I.n+Uf(L-l).m-I,n
+ u : L,m,n - I _
I) =
= -(2cr + 68 + 12~)
+ (2or cos (3 + 48 cos {* cos & + 8y cos &) HAL = HGL = 0.
(48)
1026
K. c.
KING
and T. MURA
Case 3. The boundary layer I = 0 is assumed to be the interface layer between the thin film and substrate. The displacement of the atoms on this fayer is zero and this cannot be formed as the boundary equations. Therefore the boundary equations for the rigid boundary layer come from the equilibrium equations for the layer at 1 = 1. The equilibrium equations for the layer I = 1 can be obtained in a similar way as those obtained in case 1.
[t”]T=f~f,a:,a:,1?:,~l,~:)
Hi,.=
(50)
(for i = 1,3, j = 1,6)
Hi+3j
Finally, we assemble together as
(51)
eqns (41), (46) and (51)
(52)
Fig. 2. Three types of defects: (a) e~I~~~~da1inclusion with the principal axis %I,, 2a,, k,; (b) cuboidal inclusion with the characteristic dimension 2a, ,2a,, 20, ; and (c) circular prismatic dislocation loop with the radius R and Burgers vector B.
where
[oj’=(IT~,zz:,tz:
,..,)
ii;,u”:,t$ -+L a1
,...I(
-+L -+q. ,u2 ,743
(53)
We solve for $(I, &. &) from the equation above, then substitute it into eqn (39) to obtain the solution for the displacement of the atoms.
3. NUMERICAL
CALCULATIONS
3.1. Small defects in an infinite lattice For elasticity calculations, two types of inchtsions, which are ellipsoidal inclusion and cuboidal inclusion, and prismatic dislocation loop are chosen as the examples in the present work. The coordinate system used for each case is shown in Fig. 2(a), (b) and (c). The elastic constants used for the ellipsoidal inclusion are: Cl1 = 3.398, C,, = 0.586, CM= 0.99 (unit = lOi’% m-3 and those used for the other two
cases are isotropic with Poisson’s ratio v = 0.25 and the Young’s modulus E = 40 (unit = lOroN me2). The characteristic dimensions of the ellipsoidal are a, = 5b, a2 = 4b, a, = 36 (see Fig. 2(a)>. The eigendistortions inside this inclusion are constant and the non-zero components are arbitrarily given as /3:, = 0.5, 8% = 0.4, fib = 0.3. The characteristic dimensions of the cuboidal inclusion are: aI = Sb, a2 = 4b, a3 = 36 (see Fig. 2(b)) and the eigenstrains inside the cuboidal inclusion are the same as those used for the ellipsoidal inclusion. The prismatic dislocation loop has the radius R = 5b and the Burgers vector is normal to the plane of the circle with the magnitude 1st = b (see Fig. 2(c)). This problem can be simulated by putting the constant eigendistortion @?I= 1 inside a disk-shaped inclusion with the radius R = 5b and the thickness t = b. The results from
Eigenstrain method for small defects in a lattice . LATTICE THEORY SOLUTION
-
1027
ELASTICITY SOLUTION
1.5
1.5
l-
1
e s
e s
. 0.5-
.
0.5
.
.
. .
0
12
3
4
5
6
7
9
0 fi 0
910
1
2
.
. 0~ 0
345678
.
1
4
3
2
X&d/b
X,(0/b
5
X,(n)/6
fb) @I (a) Fig. 3. Solution for the ellipsoidal inclusion problem. Displacement distribution of ui along the xi direction (i = 1,3) are shown in (a), (b) and (c), respectively. l
LAlTICE THEORY SOLUTION
12
3
ELASTICITY SOLUTION
c/J
ch
ch 0
-
4
5
6
7
9
910
0
1
2
3
UW
4
5
6
7
8
0
1
2
Xh)lb
4
3
5
6
x,(.)/b
(W
(a)
fc) Fig. 4. Solution for the cuboidal inclusion problem. Displacement distribution of ui along the xi direction (i = 1,3) are shown in (a), (b) and (c), respectively. elasticity calculation [ 17, 18,201 are illustrated by the solid lines in Figs 3, 4 and 5.
For the lattice theory calculation, the defect domain in the lattice has to be determined first. In order to compare with the results from continuum mechanics, the domain is determined by choosing all
the atoms inside and on the boundary of the inclusions to have the best shape approximation. The eigendistortions used inside the lattice are the same as those used in the elasticity calculation, and the eigendistortions are put between each pair of interactive atoms inside the defect domain. The 0.1 -
ELASTlCilY
SCLUTICN
. LAlTlCE THEORY SCLUTICN AT 1.1 . IAlTlCE
THEORY SCiJTlCN
. IAlTlCE
THEORY SCLUTICN AT 1.2
AT l-3
0.05.
.O.l, 0
. 1
2
, 3
, 4
. 5
x&m)/6
q(m)/6
fal
@I
, S
. 7
, 5
9
I
Fig. 5. Solution for circular prismatic loop problem. Displacement of the atoms in the axial direction and radial direction along radial direction are shown in (a) and (b), respectively. PCS %2,8--H
1028
K. C. KING
equivalent force constants a, fi, y can be obtained from eqn (25) as: a = 2.226, /? = 0.293, y = 0.101 (unit = b x 10”N m-*) for the case of the ellipsoidal inclusion, and a = 4, #I = 2, y = 0 (unit = b x IO’% m*) for the other two cases. Since the eigendistortions are constant for the present study, the summation domain R in eqn (38) can be changed into c?CI(boundary of the domain a) as mentioned before. Hence, it is necessary to identify all the atoms on the boundary of the defects. The method used to identify the boundary atoms is to identify all the atoms inside Q first. Then, those atoms with any one of the 18 interactive surrounding atoms outside of R are identified as boundary atoms. For example, the atoms which are located on the surface of the cuboidal inclusion are just the boundary atoms. The boundary atoms, which simulate the boundary of the ellipsoidal inclusion on the x2-x, plane and the disk-shaped inclusion, are illustrated by the solid circles in Fig. 6(a) and Fig. 6(b).
and T. MURA In order to overcome the difficulty of high frequency in eqn (35) this integration is integrated over several sections with equal spatial dimension. The number of sections is according to the number of periods of the kernel in each direction inside the integral domain R. Gauss-Legendre quadrature is used to carry out the integration for each part. Then we substitute the numerical value of Ejks(ftf*l’) into eqn (38) to carry out the summation. The results of the displacement field for each case are obtained and shown in Figs 3, 4 and 5. In these figures, the solid circles represent the displacement of the atoms. The lattice structure after deformation is also illustrated in Fig. 7(a) for the case of prismatic dislocation loop. The vertical axis is the symmetry axis of the prismatic dislocation loop and the horizontal axis is the radial direction of the prismatic dislocation loop. The equilibrium position of the atoms after deformation is illustrated by solid circles and the cross points of the grid lines in this figure represent the position of the
6
0000000000000 0000000000000
(a)
a
1=8
7
1=7
6
I=6 1=5
35 ?4
1=4 3
I=3
2
1=2
1
I=1
0
l=O
0 000000
000000
0 0000000000000
12
3
4
5
6
7
8
910
*kW
0
X*(m,
(b) Fig. 6. Illustration of the boundary atoms (solid circles): (a) x,-x) plane of the ellipsoidal inclusion and (b) disk-shaped inclusion.
@I
Fig. 7. The equilibrium position of atoms as a prismatic dislocation loop embedded in a lattice: (a) dislocation loop in a infinite domain and (b) dislocation loop in a thin film. The shaded area represents the position of the dislocation loop.
Eigenstrain method for small defects in a lattice atoms in perfect lattice site. The eigendistortion, which is used to simulate the dislocation loop, is placed in the shaded region. 3.2. Small defects in a thin film A prismatic dislocation loop embedded in the thin film with nine layers of atoms is calculated. The layer labeled “0” is assumed to be the rigid boundary (displacement is zero) and the layer labeled “8” is assumed to be the free boundary. As explained before, the dislocation loop can be simulated by an disk-shaped inclusion; this disk-shaped inclusion is arbitrarily inserted between layer 4 and layer 5. The elastic constants, radius of the dislocation loop and the Burgers vector are the same as those used for the dislocation loop embedded in the infinite lattice. The boundary atoms can be identified in the same way as before. The result is obtained by solving for ii,(l) from eqn (52) first, then substitute t7,(() into eqn (39) and integrate numerically. The same integration technique is used here. The position of atoms after deformation is illustrated in Fig. 7(b). 4. DISCUSSION The eigenstrain Green function has been derived and used to calculate the distortion field caused by several kinds of defects with complicated shape. The most efficient way to solve the problem of defects in lattice is to construct the eigenstrain Green function table through the use of eqn (35). These results are calculated in the real space instead of the reciprocal space, and they can be used to obtain the distortion field directly through the use of eqn (38). The lattice model used in this work is just an illustration and the eigenstrain Green function can be derived in a similar way for different kinds of lattice structure. The distortion field obtained from the lattice calculation are compared to the results obtained from the elasticity calculation. The comparison indicates a good agreement as the shapes of the defects in lattice structure are similar to those in elasticity calculation. For example, the cubical shape defects in the simple cubic lattice can be matched exactly by the cubical inclusion and the distortion fields are very close to those calculated from the elasticity calculation (see Fig. 4). However, the distortion fields are different as the defect shapes prescribed by the continuous model are approximated for the other two cases, ellipsoidal and prismatic dislocation loop. They are shown in Figs 3 and 5. For example, the strain (slope of the displacement) calculated from the lattice theory is larger in the middle of the defect and smaller near the boundary when compared with the elasticity solution -(Fig. 3: ellipsoidal inclusion); the strain fields near the edge of dislocation loop are larger than those expected from the elasticity calculation, and this can be seen in Fig. S(b) that the displacement of atoms in the x2 direction at level 1 with the distance of 4b and
1029
66 from the center are significantly large. This difference decays for atoms away from the defects. This suggests that the exact shape of the defects is important for the displacement of the atoms in the vicinity of the defects. This implies that the stress concentration near the defect could be significantly different from the stress concentration calculated from the elasticity calculation. Especially for those at the kink edge of the defect domain, which is very true for the real lattice structural of the materials. For very small defects (e.g. point defects), the stress and strain fields are of great difference as compared to the elasticity calculation. It follows that the elastic strain energy can be different. This difference could be important for the nucleation of precipitates in solid, since the nucleus is sufficiently small. The results of prismatic dislocation loop embedded a thin film are significantly different from those of the same defect in an infinite lattice (see Fig. 7(a) and (b)). The displacements of atoms are no longer symmetrical to the plane of the dislocation loop which is true for this defect in an infinite lattice. The displacement of the atoms near the edge of dislocation in the radial direction is relatively large as compared to the case in infinite domain. The atoms between the defects and the substrate cannot move too much because of the rigid boundary condition. However, as more layers are put in between the dislocation loop and the rigid boundary, similar results are obtained. We conclude that the defect in a thin film with a rigid boundary can simulate the defect in a half-space lattice as sufficient layers are put in between the defect and rigid boundary. The disturbance of the defect is more significant for the atoms between the defect and the free surface. The difference of distortion fields suggests the importance of surface effect for the defects near the surface. Not only are the stress fields strongly influenced, but the elastic strain energies are of great difference. Furthermore, the elastic strain energy could change significantly as the dislocation loop moves toward the surface. Therefore, from the results of this work we are encouraged to calculate the stable criterion for a dislocation loop near the surface in the future. 5. CONCLUSION On the basis of this work we can conclude the following.
(1) Two Green functions have been derived. One is used for the case of distributed body force and the other one is used for the case of distributed eigenstrains. (2) The eigenstrain Green function is important and useful for defects with complicated shape, and the distortion fields can be calculated straightforward as the table of eigenstrain Green function is constructed.
K. C. KING and T. MURA
1030
(3) The distortion fields of dislocation type defects in
a lattice can be solved by this method directly, and this is not done before through the use of Kansaki force. (4) For small defects, the elasticity solution will give significant error especially near the boundary of the defect, and this local effect is important for This the analysis of stress concentration. difficulty can be resolved by using the lattice theory and eigenstrain method. (5) The effect of a defect in a thin film is calculated in this work, and the results are significantly different when compared with the same defect in an infinite lattice. This method is also applicable for different models of atomic structure. The potential of this method is that it can be further extended to calculate the elastic strain energy of the defect in a solid and the crack in a lattice can also be analyzed with the calculated equivalent strains. Acknowledgemenr-This research was supported by the U.S. ONR NOO014-80-K-0102 and also partially by the U.S.
AR0 DAAL03-89-K-0019.
REFERENCES I. Sato A. and Watanabe Y., J. Phys. Chem. Solids 49, 529 (1988).
2. Kuriyama M., J. Phy. Sot. Jap. 23, 1369 (1967). 3. Huang K., Proc. Roy. Sot. A190, 102 (1947). A. W. and Gorman R. L., Phys. Rev. 102, 4. Overhauser 676 (1956). 5. Kanzaki H., J. Phys. Chem. Solids 2, 24 (1957). 6. Tewary V. K., AERE HarweN Rep. TP, 388 (1969). 7. Tewary V. K., J. Phys. F: Metal Phys. 3, 1515 (1973). 8. Thomson R., Tewary V. K. and Masuda-Jindo K., J. Mar. Res. 2, 619 (1987). 9. Hsieh C. and Thomson R., J. appl. Phys. 44, 2051 (1973). 10. Slater J. C., Quantum Theory of Molecules and Solid. McGraw-Hill, New York (1974). 11. Morinaga M., Yukawa N., Adachi H., Mura T. and Sato A., Ann. Mig. Tokyo, Conf. Abstr. Japan Inst. Metals 272 (1986). 12. Mura T., Micromechanics of Defects in Solids. Martinus Nijhoff (1987). 13. Mura T., Eigensrrains in Lattice Theory l-503-519. University of Waterloo Press, Quebec (1977). 14. Born M. and Huang K., Dynamic Theory of Crystal Lattices. Oxford (1954). 15. Ashcroft N. W. and Mermin N. D., Solid State Physics. Cornell University, New York (1976). 16. Chiu Y. P., J. appl. Mech. 44, 587 (1977). 17. Kinoshita N. and Mura T., Phys. Stat. Sol. 5, 759 (1971). 18. Mura T. and Cheng P. C., J. appl. Mech. 44,591 (1977). 19. Eshelby J. D., Proc. Roy. Sot. A241, 376 (1957). 20. Kroupa F., Czech, J. Phys. BlO, 284 (1960). 21. Mindlin R. D., J. Elas. 2, 217 (1972).
1991
THE EIGENSTRAIN
0022-3697191 $3.00 + 0.00 0 1991 Pergamon Press plc
METHOD FOR SMALL DEFECTS IN A LATTICE
K. C. KING and T. MURA Theoretical and Applied Mechanics, Northwestern University, IL 60201, U.S.A. (Received 19 November 1990; accepted 21 March 1991)
Abstract-The eigenstrain (eigendistortion) method is used in the lattice theory to analyze the disturbance caused by the defects in a lattice. The equilibrium equations for the atoms are derived through the use of lattice theory, and the harmonic assumption for the potential energy is used for the derivation. The eigenstrain (eigendistortion) is used to model the defects, which could have a complicated shape. The results are compared with the results obtained from the elasticity theory. A simple cubic model is used to model the structure of the lattice. The Green function is obtained from the equilibrium equation by using the discretized Fourier transformation method. Therefore, the disturbance caused by the defects in the infinite lattice can be analyzed by the summation of the Green function with the strength of the eigendistortion. The problem for the defects in the thin film is solved. The thin film is coated on a substratum and bound by a free surface on the other end. The substratum is assumed to be rigid. The result is compared with the result for the same defect in the infinite lattice. The comparison indicates that the effect of the boundaries of the thin film is important. Keywords: Lattice, thin film, eigenstrain, defects, prismatic dislocation loop, inclusion.
1. INTRODUCTION
The precise description of the distortion of the surrounding atoms near defects in solid crystals is important for the calculation of the energy of defects [l], the discussion of X-ray diffraction effects [2, 31 and also the change in electrical conductivity of metals [4]. Such calculations have been done by Kanzaki [5] using the concept of the Green function, which is to calculate the displacement of the atom L caused by a unit body force acting on the atom P. Therefore, the displacement of a ‘normal’ atom caused by a substitutional or interstitial defect in a perfect lattice can be obtained by solving the equations formulated by the use of the Green function [S]. This method has been improved by Tewary [6], who introduced the ‘defect Green function’ and this defect Green function can be applied to the lattice structure with the crack surface inside it. Several works have been done by Thomson et al. using the defect Green function, such as, lattice distortion due to gas interstitial in bee metals [7], three-dimensional crack with a kink [8] and lattice theory of fracture and crack creep [9]. However, the defect problem like dislocation cannot be solved from these methods in the above, since the dislocation is prescribed by the displacement jump (Burgers vector) on the dislocation plane instead of the force acting on each atom. Therefore, eigenstrain is used in the lattice theory to overcome this kind of difficulty. Eigenstrain is a general name to describe misfit strain, thermal expansion strain,
phase transformation strain, plastic strain, etc. For example, the misfit induced by an interstitial atom can be estimated by a fundamental electronic theory [lo, 111; the misfit for the substitutional atom can be determined by the size difference between the atoms; the misfit for the dislocation is the Burgers vector [12]. The idea here is to use this misfit to calculate the distortion field without the use of the Green function in the above. This is so called ‘eigenstrain (eigendistortion) method’ which is widely used in continuum mechanics. The application of this idea in lattice theory was first introduced by Mura to calculate the moving dislocation [13], and it was again used to calculate the octahedral defect in a bee lattice and the growth of nitrogen defect cluster by Sato et al. In this work, we further improved the eigenstrain method in lattice theory by introducing the ‘eigenstrain Green function’. This eigenstrain Green function is used to calculate the distortion of atom G caused by a unit eigenstrain (misfit) between two atoms !’ and e*, therefore, the total distortion field can be obtained by the summation of the eigenstrain Green function with the strength of eigenstrains. This method can solve the dislocation defects as well as the interstitial, substitutional defects. Furthermore, the defects with complicated shape can be systemically solved with the help of the eigenstrain Green function, and this is the benefit of this method. However, in order to simplify our mathematical task considerably, a simple cubic lattice is chosen to illustrate this method although no significant difficulty can be seen here as the other lattice
1019
K. C. KING and T. MURA
1020
structures such as bee or fee are chosen. The importance of this method is that it can be further extended to model cracks and the critical criterion by calculating the equivalent eigenstrains; it is under preparation now. We also investigate the free boundary effects for the defects in lattice. Traction free boundary effects are used for the free boundary in this analysis. The equilibrium equation is formulated for a finite layer lattice domain (thin film). The thin film is coated on a substrate and the substrate is assumed to be rigid, which implies a rigid boundary condition on the other side of the thin film. The same model for the atomic structure is used for the thin film and no special surface lattice structure is assumed. In the following, lattice theory [14] and the harmonic assumption of the potential energy [ 151are used to formulate the equilibrium equations for the atoms [13]. A simple cubic lattice with the consideration of nearest and second nearest atoms interaction is assumed for the mathematical formulation. A small defect domain with cubical shape and ellipsoidal shape is calculated; the prismatic dislocation loop is also calculated in this work. The results are compared with the results from elasticity calculations [l&20]. The equations for the defects in the thin film are then formulated. The prismatic dislocation loop defect is used for the calculation. The result is compared with the result for the same defect in the infinite lattice. The comparison indicates that the effect of the boundaries of the thin film is significant.
2. THEORETICAL
DERIVATION
2.1. Equation of motion The harmonic approximation of the potential energy is introduced first in this section [ 151.Then the equation of motion will be derived by using the Hamilton principle. Finally, the equilibrium equation can be obtained from the equation of motion [13]. Let Y(LG’) be the potential energy of a pair of atoms L and e’, where G and d’ are the label of atoms at perfect lattice sites and !!’ is the vector measured from atom G to atom L’. For a perfect lattice structure, the total potential energy U of the system can be determined by the summation of the contributions of all distinct pairs as U =;I
Y(ee’). 1.f’
VY (et’)
a=~~,Y(pl.)+~~(u(C)-u(C’)) (.I’ 1.J + ; c [(u(d) - u(P)) f,c?
V12Y(68' ) + O(u’).
(3)
The gradient of potential energy, VY(c/‘), is known as the interactive force between the two atoms / and L’. Therefore, the summation of interaction force, Zr VY(!J’), is zero, since there is no net force on any atom in equilibrium condition. It follows that the linear term in eqn (3) is zero. Also, the higher order terms are neglected since small displacement is assumed in this work. Then the total potential energy can be written as follows: U = Ueq+ Uham
(4)
u’q =; 1, Y(LG’) /.I
(5)
Uhitrm= ; E, [(u(L) - u(L’)) v]*Y(ee’). 1.1
(6)
where
The equilibrium potential energy U'q is the total potential energy for the system without any defects and it can be set to zero as a reference state because the absolute potential is not important in this work. Then, the only term left for the total potential energy is Uharm. This is called the harmonic approximation. The harmonic potential energy can be written in a more general form as u=u
harm= ;;,
@,(dd’)u,(d)u,(d’),
(7)
where @,(ee’) = 6,,, 1 Y.,,(dd”) - Y.,,(&/‘). I”
(8)
The double indices appearing in the equations follow the summation convention from 1 to 3, and the commas in the subscripts denote the derivatives. The total kinetic energy T of the system is
(1)
For a distorted lattice structure caused by dislocations or inclusions, this potential energy can be calculated by the summation of the potential energy of each pair of atoms measured at the equilibrium position as u = ; 1 Y (ee’ + u(d) - I@‘)), !,L’
where u(8) is the displacement vector of atom / measured from its original lattice site. For small displacements, the total potential energy U can be expressed by Taylor’s series as
(2)
T = ;
1 Mti, (f)&(f), I
(9)
where M is the mass for each atom. The equation of motion for each atom e can be derived by using the Hamilton principle S[{,:(U-T)df]=O
(10)
Eigenstrain method for small defects in a lattice
1021
The equilibrium equation can be obtained from eqn (19) by neglecting the dynamic term as Miii(e)
=
-1 @ij(ee’)uj(e’),
(11)
f’
where 6 is the symbol for variation. It can be seen from eqn (11) that @ij(!~‘) are the force constants between the atoms L’and c”, and these force constants have the following properties ~~ij(ec’)=~~ij(~L’)=o
(12)
@,(Ld’) = Gji(ee’) = @,(e’t),
(13)
which come from the rigid body motion and the symmetric conditions [13, IS]. Then eqn (11) can be written as
1 o,(ee’)u,(e’) C’
where s”(/*L) is the Kronecker delta and it has the following properties: 8(/*k) = 1 as e* = / otherwise $((e*e) = 0. For a single unit eigendistortion I$(/*/‘) applied between atoms d* and L”, eqn (20) can be rewritten in the following ;@,(ee’)u$‘)
M&(d) = -1
@ij(Ll’)[Uj(f’)
-
= ~,(e*d’)z,:.(e*e’)x,(e~e*)~(‘(e*e),
(14)
Uj(f)].
C’
We assume the difference of the displacement vectors between e and !’ is U,(C) - Uj(/) = /I,,(ee’)x,(/‘c), where /Iy(L’e’) are the elastic distortions. tion A’,(~‘z!) is defined as
= 1 ~ij(ee’)8t(ee’)x,(c’e) I’
(1% The nota-
&(e’e) = X,(/‘) -x,(e),
(21) where the superscript E denotes the displacement caused by a unit eigendistortion. The right-hand term in eqn (21) can be interpreted as a fictitious force which is equivalent to the effective force caused by the eigendistortions. Thus, as a unit body force I,(/*) in the i-direction is applied on the atom c!*, eqn (20) can be expressed as follows
(16) ; @,,(e/‘,UF(e’, = &(d*)s(e*d),
where X,(e) is the position coordinate component k of the atom C! in a perfect lattice site, then the equation of motion (14) can be written as MYi,
= -c
@ij(de’)pkj(de’)X#‘e).
(17)
(’ However, if material is involved in any inelastic deformation, then eqn (17) is no longer valid because /Ikj(eL”) are not compatible. It is necessary to reformulate eqn (15) as u,(P) - u,(C) = [fikj(LL’) + B$(d/‘)]X#‘L),
(18)
where p$.(/L”) are the eigendistortions. The physical meaning of the eigendistortions can be interpreted as the stress-free transformation strains. For example, the differences of lattice parameters induced by the phase transformation are the eigendistortions. The definition here is similar to the definition in Micromechanics of Defects in Solids [12] that the total strain is the sum of elastic strain and eigenstrain. Substituting eqn (18) into eqn (17) and using the properties of the force constants in eqn (12), then the equations of motion become M&(e) = -c
@,(ee’)uj(e’) C’ + c ~ij(ee’)B~(e~‘)x,(~‘~). C’
(19)
(22)
where the superscript G denotes the displacement caused by a unit body force. Equations (21) and (22) are the basic equations for the following derivation. One is used for the case of distributed eigendistortions and the other one is used for the case of distributed body forces. The Green function for a unit eigendistortion can be derived from eqn (21) and the Green function for a unit body force can be derived from eqn (22). It must be pointed out here the reason to obtain two Green functions in this work. The eigenstrains (eigendistortions) used in the continuum mechanics are continuous functions and they are defined continuously in the defect domain. Therefore, the effect of the eigenstrains are just like body forces and the solution can be obtained by the integration of the Green function (obtained from unit body force) with the strength of the eigenstrain. However, the eigendistortions used in the lattice theory are not continuous functions and they are defined between any two atoms inside the defect domain. It follows that the fictitious force caused by the eigendistortion will act on two atoms simultaneously. Therefore a new Green function has to be defined to calculate the effect of the fictitious forces caused by the eigendistortions. The details of the derivation and the application of these Green functions will be discussed after we set up the model for the atomic structure.
1022
K. C. KING and T. MURA
2.2. Model for atomic structure The force constants @,,(e/‘) can be obtained from the potential energy Y(ee’) through the use of eqn (8). However, it is not convenient for the present work. A relation between force constants and elastic constants is more important as we want to compare the results between lattice theory calculation and the elasticity calculation. Thus, the simple cubic model proposed by Mindlin [21] is used. These relation equations will be shown at the end of this section, A simple cubic lattice with nearest and secondnearest atoms interaction is the model chosen for the present work. In this model each atom is interactive with 18 surrounding atoms, which are six nearest atoms and 12 second-nearest atoms. Let tl be the central force constant (i.e. force per unit relative displacement) between any nearest atoms, p be the central force constant between any second-nearest atoms and y be the angular force constant between any three, successive, non-collinear atoms (see Fig. 1). The variable b in Fig. 1 is the lattice parameter. The equation of motion for this model is given as M@l.” = #;+
@,,(f,m,n;fIfr l,m,n)= %,U,m,n;hm
-tl
+ 1,n)
=@,,(I,m,n;I,m,n
*l)=
-47
@,,(I,m,n;I+l,mfl,n) =@,,(i,m,n;t+
I,m,n + 1)
=@,,(l,m,n;I-
I,m rf 1,n)
=@,,(Z,m,n;I-l,m,n+l)=-j3
O,,(l,m,n;I,m
+ l,n)=
f-h(bm,n;l,m,n
-u
It 1)
=&(l,m,n;I+l,m,n)= ~~~(l,m~n;l,m
-47
+ I,n -t 1)
=@,,(I,m,n;I,m - 1,n + 1) =@,,(l,m,n;l+
l,m ILt:1,n)
=@,,(I,m,n;I-
1,m Ifr l,n)=
@33(1,m,n;i,m,n _t I)= -u
I.m,n+ U:- 1.m.n _ 2u;m.s)
/+i.m+in + $ l.m- 18+ *i+ I.m- 1.n + a@, +U~-l.m+l,n+ll~+I,nin~f+~~-l.m.~-l
@,,(l,m,n;Irt
1,m,n)=@33(I,m,n;I,m
_t I,n)
= -4y %(I,m,n;f+l,m,n
+U;+i.m.n-l +U~-Im.?l+I _ &@.‘7)
t- 1)
+ (p + Y)(U;+ i.m+1.N+ &- I.ffl_1.n
=@33(1,m,n;t-l,m,n
_ U/+ I.rn- I.”_ U/-. 1.m-tI.n
=@&(l,m,n;l,m + 1,n Lt: I)
2
-p
* 1)
2
=@,,(I,m,n;l,m
+~:+I,m.n+I+U~--I.m.n-I+~I+I.m.n+I
- 1,n + 1)= -/3
3 +
&
l,m.n - 1 _
+4~(U~m+l.n
+24
yl
Up
+
-
@z-
1.m.n - i _
La+
@:-
@12(lrm,n;[+ 1,m + 1,n)
1.m.n + I)
=@,z(l,m,n;f-
U:n.n+l
1_ 4u$m‘“)*
(23)
Here the subscript denotes the component of displacement vector and the superscripts denote the position coordinates of the atom in a perfect lattice site. The remaining two equations of motion can be obtained by cyclical permutation of subscripts and superscript. If we drop the term with eigendistortions @$$‘8“) in eqn (19), then we can determine the force constants by comparing eqns (19) and (23) as follows @,,(1,m, n; 1,m, n) = @22tl, m, n; 1,m, n>
=@,,(f,m,n;i,m,n)=2or+8~+16y
Fig. 1. Schematic illustration of the spring constants a, & . . . . .. y used in the simple culxc model.
I,m + I,n)= -(b +y)
@13(lrm,n;i+ Lm,n t- I) =@,,(Z,m,n;t-l,m,n+l)=-(/3+y) @,3(I,m,n;I,m+l,nf
1)
=@23(Irm,n;l,m - l,n + 1)= -(p
+y)
clrI2 = 6, @13= @3r
(24)
%3 = Q13,.
The relations between ~1,p and y in this model and the elastic constants were derived by Mindlin [21] Ithrough finite difference method as follows
G4=;
(B +2Yh
Eigenstrain method for small defects in a lattice where C, 1, C,*, C, are the Voigt elastic constants for cubic symmetric material. This relation (25) is used to compare the solution obtained from lattice theory calculation with that obtained from the elasticity calculation. The solution for defects embedded in an infinite lattice and thin layer are studied in the following sections.
1023
Solving eqn (28) for QF, we can obtain
l Nij(5)zi(e*) ew[-iW, + 452+ rtdl, (30) cY=$ o(r) where D and Nij are the determinant and cofactors of the matrix K. The displacement u?(e) can be obtained by substituting Is,” into eqn (26) as follows
2.3. Green function Two Green functions, which are explained in the previous section, are derived in this section and they are used to calculate the displacement of the atoms caused by the defects. The Green functions obtained in this section are valid for an infinite cubic symmetric lattice. The discretized Fourier transformation method is used in the present work [ 131.The discretized form of the displacement is
UP)
1 = (271)3 -*
x exp{W -P)<, + Cm- 4% (31) From the equations above, the Green function for a unit body force can be solved as
uj(e) =
1 G,(ed*) = m
x exp[i(E, + 4
+ G)l dt, G dt,,
--II
(26)
x exp{W -PK + (m - q)L
where b is the lattice parameter, and the discretized form for the Kronecker delta is
x ew{W -PK + (m - q)t2
where (I, m, n) and (p, q, r) are the position coordinates of the atoms L+and e* at the perfect lattice site, respectively. Substituting eqn (26) and eqn (27) into eqn (22), then eqn (22) becomes K/jr77 = i Z#*)exp[-
i(P5, + q
(32) The physical meaning of the Green function G,,(LL*) is the displacement of atom L(Z,m, n) in the j direction when the unit force in the i direction is applied on the atom /*(p, q, r). However this Green function is inadequate for the eigenstrain method. A similar function Ej&!*L’) (Green function for a unit eigendistortion) is derived for practical purposes. The physical meaning of this function E,,(/!*!‘) is to describe the displacement of atom L’(Z,m, n) in the j-direction when the unit eigendistortion Zg(t*/‘) is applied between atom /*(p, q, r) and one of its surrounding atoms 8’. This function Elks(L/*d’) can be derived in the same way as the Green function is derived. First, substituting eqn (26) and eqn (27) into eqn (21), eqn (21) becomes
(28)
where the coefficients Ku are x exp] - Z(P5, + q5r +
K,, = Za(cos 5, - 1) + 4fi(cos 5, cos c2 + cos 5, cos {, - 2) + 8y (cos 52 + cos 5, - 2) Kz2 = 2cr(cos & - 1) + 4B(cos r, cos
&)I. (33)
Then, solve for 17,”and substitute ri,! into eqn (26) to obtain the displacement u,?(e) (which is caused by the effect of unit eigendistortion) as
+ cos {r cos & - 2) + 8y(cos r, + cos {, - 2) K,, = 2a(cos r, - 1) + 4/?(cos 5, cos &
+ cos & cos C3- 2) + 8y (cos {, + cos & - 2) KL2= K2, = -4(fi
1 Yw) = Qn jjj @,(e*e’)x,(re*) -n
+ y)sine, sin & X
K13 = Kj2 = -4(/? + y)sin & sin & K,, = K,, = -4(/I
+ y)sin 5, sin &.
(2%
Ni’(5gt(;*“)
exp{i[(Z-p)r,
+ (m -q)& (34)
K. C. KING and T. MURA
1024
From this equation, the Green function Ejks(dL*d’) is defined as z
1 = @-#
E,,(ee*e’)
N,.(5) @,,(e*e’)x,(c’e*) i D(5)
x ev{W -_pK + (m - q)ti + (n - rkll dt!, di%d53.
(35)
This function has the following property ; E&(ee*e’)
= 0
Ejks(/d*/‘) = 0 (for L’ $ A),
(36) (37)
where the domain A represents the 18 atoms surrounding the atom L*. Equation (36) can be proved through the use of eqn (13) and eqn (16). Equation (37) comes from the fact that @is(/*e’) = 0 in eqn (35) for d’# A. These Green functions are valid for the infinite cubic symmetric lattice problems. The application of EjkJ(,,(G/*L’)is investigated in the following section. 2.4. Small defects in the infinite lattice With the help of the Green functions derived in the previous section, it is easier to treat defect problems in the lattice. Since eigendistortions can simulate defects well, this work is emphasized on the use of the Green function for a unit eigendistortion. The method of using the other Green function (the one for unit body force) is the same and does not require further explanation in this work. For a given defect domain R prescribed by the eigendistortions, the displacements of the atoms can be obtained by the summation of Ejks(Ld*d’) over R as a,(l)=
z
C Pk:(~*~‘)E,,,(~b*~‘). /’ /*En /‘E(RnA)
problem. Thus, a different method is used in this section. The thin film is assumed to be coated on a substrate and the substrate is considered to be very rigid when compared to the thin film. It follows that the displacements at the interface are zero in this analysis. The other boundary is assumed to be a traction free boundary. The displacement field caused by the defects inside the thin film is solved in the following. Let the x1 direction be the direction normal to the plane of the thin film. The discretized Fourier transform is applied to the x, and xj directions. Therefore, the discretized form of the displacement is
(38)
The summation domain n depends on the position of the atom !*. Only the atoms inside R are necessary for the summation because the eigendistortions b$$(e*/‘) are zero outside the domain R. Furthermore, if the eigendistortions are constant inside the domain R, the summation domain R (for the atom d*) can be changed into the domain aR, which is the boundary of R. This comes from the property of the function Ejks(/t*ta’) expressed in eqn (36). Thus, for the case of a constant eigendistortion distribution, only the eigendistortions on the inner surface of the defects are necessary for the computation. 2.5. Small defects in a thin$lm As mentioned before, the Green functions derived in the previous section are no longer valid for this
x exdi(mt;, + n&)1dt2 G
(3%
and for the Kronecker delta is
b(e*e) E 6,,6,6,” =
n ss
&s,, = --I
--II
x expIi[(m - q)& + (n - r)&l} d&d&,
(40)
where (I, m, n) and (p, q, r) are the position coordinates of the atoms L and L* at the perfect lattice site, respectively. Equation (39) can be integrated to obtain the solution of displacement field after tij(f, &, t3) is solved. The unknown variables 1?,(1,&, 5,) can be obtained by solving a system of linear equations, which come from the equilibrium equations of each layer. There are three equations (three components of the displacements) for each layer and these linear equations for each layer have to be solved simultaneously. Since each layer is only linked to its neighbor layers (upper and lower), a banded matrix will be formed after assembling all the equations together. The banded matrix can save computer time as well as working space. The system of equations is separated into three parts according to the position of the layers as: (1) equations for the layers inside the thin film (case 1); (2) boundary equations for the free boundary (case 2) and (3) boundary equations for the rigid boundary (case 3). Case 1. The layers 1 are assumed to be located inside the thin layer. The linear equations for these layers inside the thin film can be obtained by substituting eqns (39) and (40) into eqn (20) as
where
-(/+I), ug+ 1) ++ I)), (42) u1
Eigenstrain method for small defects in a lattice the notation
tij in eqn (42) denotes iij(l, &., &), and
C1(U~Lm.~ + U$L.m.n
1025 _
2U;L.m.n)
p(ufLmt
+
Ia+
I
+u~L.m-I.n-I+U~L.m+In-I+UfL.m-I.ntl +
u$h!-
+
(fj
IAm+ +
1.~ +
u$(L-
I),m-
lflf
I+ Up-
r)(u;L*“+
IJI _
,5u~~m~)
In-
I
_~~L.mtlr-I_~~L.m-l,n+t)
x; = 6,JC c @,(C*e’)/?$(t!*e’)X,(L’L*) R t’ x ev[--q(b
+ d31
py
+Y)(U:(L-I).m-I,n_Ulf(L-I),m+I.n)
.+. +
(u f Lmn
t
I +
u ; L,m.n -
I
+ uf,“- Ihn _ gr4;L.m.n) TY@fL.“fI.” _ ufL.W Ln)= 0 q(r&,m.n~ 1+ wjiL.m.n+ 1_ &pwl) (43)
+~(UfL~mtl.ntl+~~L.m-l~n-l fur
H,, = H,, = 2p cos & + 4y
+L.m+I.n-I+~~L.m-1.~+l+~~,L-I~.mr+l
+U:,L-I),m.n-I_6U:L,m.n
H,, = H,, = 28 cos & + 4y H,, = Hz7 = -H,,
= -Hz,
H,, = H3, = -H,,
= - H3, = i[2(B + y)sin c,]
+(p
= i[2@ f y)sin &]
+y)(U~L.nt+I.ntI+U~L,m-~."-I _U~L.m-lCti
_ti;L,m+l.n-l T(fl
H,, = H,? = u + 2~(cos t2 + cos &) Hz6 = H,5 = -4(p
H,,=
)
+Y)(*~(L-lf.mn+l
) _M~,L-if.nr,n-l)
+4QY(U~L.m+I.n+U~L.m-I.n+~~,L-I~.m.n
+ y)sin &sin & _ 3UfL.m.n)-+ y(UfL.m.n+
I _ UjtL.m.n-
I) =
0,
(45)
-(2a: +SP + 16y)+8y(cos~,+cos~,)
&=-(2or+8p+16y) + (2a cos & + 48 cos & cos & + 8y cos &) & = -(2a
f 88 + 16~)
+ (2~7cos & + 4pcos t2 cos & + 8y cos &)
where only the equation for I= +L is used for the present work because this layer is the boundary layer with traction free assumption. The ~uilib~um equations for the free boundary are obtained by substituting eqn (39) into eqn (45) as
WtLl,j x6)
H,, = H,, = H,, = Hz4 = H,, = H,2 = Hj4 = H3* = 0.
(44
Case 2. The layer +L is assumed to be a traction free boundary. The boundary equations come from the equilibrium equations for the atoms on the free boundary. For each boundary atom, the equilibrium equations are derived through the summation of the forces which come from the interaction of its neighborhood atoms inside the thin film. The summation should be zero because of the traction free condition. The general equations for the free boundary at I= t_L are expressed in the following
cL(u ? (f - I)JW_ U~L,“L)
’ [u+Ll(6x
[U+L]r = f&j”- ‘1,&-I),
(46)
@-“,
rjf, j& n$> (47)
H:,L = 2p cos & + 4y H&L = 2fl cos & + 4y HAL = - HhL = i[2y sin &] H:SL = - HGL = i[2y sin <,I
HhL = HAL = -i[Z(B + y)sin ~$1 HAL = HGL = --i[2(/? + g)sin &] HAL = d + 2j?(cos & + cos &) HgL = HGL = -4(p
+U~,L-I~,m--I,n+U~~L-I~.mtI.n_4~~L.m.n)
HhL=
T(jj +y)(U~(L-I).mtI.n_U~,L-I),m-I..
HGL = -(2cr + 68 + 12~)
+Uf,L-If.m.n+l_~f,L-l~.m.n-I
H$
+U~L.m.n-l_4ufL.m,n)~y~U~L.m-ln_U~L.m+1.n
t I) = 0
-(c( +4j
+ y)sin C2sin &
+8y)+4y(cos&+cos’r,)
+ (2or cos r* + 48 cos tz cos & f 8y cos &)
I
f2y(ufL‘“+‘“+u~tJn-In+u~t.m,nt1
UfL.m.n
0,
where
+~(UTI~--I).m+I.n+Uf(L-l).m-I,n
+ u : L,m,n - I _
I) =
= -(2cr + 68 + 12~)
+ (2or cos (3 + 48 cos {* cos & + 8y cos &) HAL = HGL = 0.
(48)
1026
K. c.
KING
and T. MURA
Case 3. The boundary layer I = 0 is assumed to be the interface layer between the thin film and substrate. The displacement of the atoms on this fayer is zero and this cannot be formed as the boundary equations. Therefore the boundary equations for the rigid boundary layer come from the equilibrium equations for the layer at 1 = 1. The equilibrium equations for the layer I = 1 can be obtained in a similar way as those obtained in case 1.
[t”]T=f~f,a:,a:,1?:,~l,~:)
Hi,.=
(50)
(for i = 1,3, j = 1,6)
Hi+3j
Finally, we assemble together as
(51)
eqns (41), (46) and (51)
(52)
Fig. 2. Three types of defects: (a) e~I~~~~da1inclusion with the principal axis %I,, 2a,, k,; (b) cuboidal inclusion with the characteristic dimension 2a, ,2a,, 20, ; and (c) circular prismatic dislocation loop with the radius R and Burgers vector B.
where
[oj’=(IT~,zz:,tz:
,..,)
ii;,u”:,t$ -+L a1
,...I(
-+L -+q. ,u2 ,743
(53)
We solve for $(I, &. &) from the equation above, then substitute it into eqn (39) to obtain the solution for the displacement of the atoms.
3. NUMERICAL
CALCULATIONS
3.1. Small defects in an infinite lattice For elasticity calculations, two types of inchtsions, which are ellipsoidal inclusion and cuboidal inclusion, and prismatic dislocation loop are chosen as the examples in the present work. The coordinate system used for each case is shown in Fig. 2(a), (b) and (c). The elastic constants used for the ellipsoidal inclusion are: Cl1 = 3.398, C,, = 0.586, CM= 0.99 (unit = lOi’% m-3 and those used for the other two
cases are isotropic with Poisson’s ratio v = 0.25 and the Young’s modulus E = 40 (unit = lOroN me2). The characteristic dimensions of the ellipsoidal are a, = 5b, a2 = 4b, a, = 36 (see Fig. 2(a)>. The eigendistortions inside this inclusion are constant and the non-zero components are arbitrarily given as /3:, = 0.5, 8% = 0.4, fib = 0.3. The characteristic dimensions of the cuboidal inclusion are: aI = Sb, a2 = 4b, a3 = 36 (see Fig. 2(b)) and the eigenstrains inside the cuboidal inclusion are the same as those used for the ellipsoidal inclusion. The prismatic dislocation loop has the radius R = 5b and the Burgers vector is normal to the plane of the circle with the magnitude 1st = b (see Fig. 2(c)). This problem can be simulated by putting the constant eigendistortion @?I= 1 inside a disk-shaped inclusion with the radius R = 5b and the thickness t = b. The results from
Eigenstrain method for small defects in a lattice . LATTICE THEORY SOLUTION
-
1027
ELASTICITY SOLUTION
1.5
1.5
l-
1
e s
e s
. 0.5-
.
0.5
.
.
. .
0
12
3
4
5
6
7
9
0 fi 0
910
1
2
.
. 0~ 0
345678
.
1
4
3
2
X&d/b
X,(0/b
5
X,(n)/6
fb) @I (a) Fig. 3. Solution for the ellipsoidal inclusion problem. Displacement distribution of ui along the xi direction (i = 1,3) are shown in (a), (b) and (c), respectively. l
LAlTICE THEORY SOLUTION
12
3
ELASTICITY SOLUTION
c/J
ch
ch 0
-
4
5
6
7
9
910
0
1
2
3
UW
4
5
6
7
8
0
1
2
Xh)lb
4
3
5
6
x,(.)/b
(W
(a)
fc) Fig. 4. Solution for the cuboidal inclusion problem. Displacement distribution of ui along the xi direction (i = 1,3) are shown in (a), (b) and (c), respectively. elasticity calculation [ 17, 18,201 are illustrated by the solid lines in Figs 3, 4 and 5.
For the lattice theory calculation, the defect domain in the lattice has to be determined first. In order to compare with the results from continuum mechanics, the domain is determined by choosing all
the atoms inside and on the boundary of the inclusions to have the best shape approximation. The eigendistortions used inside the lattice are the same as those used in the elasticity calculation, and the eigendistortions are put between each pair of interactive atoms inside the defect domain. The 0.1 -
ELASTlCilY
SCLUTICN
. LAlTlCE THEORY SCLUTICN AT 1.1 . IAlTlCE
THEORY SCiJTlCN
. IAlTlCE
THEORY SCLUTICN AT 1.2
AT l-3
0.05.
.O.l, 0
. 1
2
, 3
, 4
. 5
x&m)/6
q(m)/6
fal
@I
, S
. 7
, 5
9
I
Fig. 5. Solution for circular prismatic loop problem. Displacement of the atoms in the axial direction and radial direction along radial direction are shown in (a) and (b), respectively. PCS %2,8--H
1028
K. C. KING
equivalent force constants a, fi, y can be obtained from eqn (25) as: a = 2.226, /? = 0.293, y = 0.101 (unit = b x 10”N m-*) for the case of the ellipsoidal inclusion, and a = 4, #I = 2, y = 0 (unit = b x IO’% m*) for the other two cases. Since the eigendistortions are constant for the present study, the summation domain R in eqn (38) can be changed into c?CI(boundary of the domain a) as mentioned before. Hence, it is necessary to identify all the atoms on the boundary of the defects. The method used to identify the boundary atoms is to identify all the atoms inside Q first. Then, those atoms with any one of the 18 interactive surrounding atoms outside of R are identified as boundary atoms. For example, the atoms which are located on the surface of the cuboidal inclusion are just the boundary atoms. The boundary atoms, which simulate the boundary of the ellipsoidal inclusion on the x2-x, plane and the disk-shaped inclusion, are illustrated by the solid circles in Fig. 6(a) and Fig. 6(b).
and T. MURA In order to overcome the difficulty of high frequency in eqn (35) this integration is integrated over several sections with equal spatial dimension. The number of sections is according to the number of periods of the kernel in each direction inside the integral domain R. Gauss-Legendre quadrature is used to carry out the integration for each part. Then we substitute the numerical value of Ejks(ftf*l’) into eqn (38) to carry out the summation. The results of the displacement field for each case are obtained and shown in Figs 3, 4 and 5. In these figures, the solid circles represent the displacement of the atoms. The lattice structure after deformation is also illustrated in Fig. 7(a) for the case of prismatic dislocation loop. The vertical axis is the symmetry axis of the prismatic dislocation loop and the horizontal axis is the radial direction of the prismatic dislocation loop. The equilibrium position of the atoms after deformation is illustrated by solid circles and the cross points of the grid lines in this figure represent the position of the
6
0000000000000 0000000000000
(a)
a
1=8
7
1=7
6
I=6 1=5
35 ?4
1=4 3
I=3
2
1=2
1
I=1
0
l=O
0 000000
000000
0 0000000000000
12
3
4
5
6
7
8
910
*kW
0
X*(m,
(b) Fig. 6. Illustration of the boundary atoms (solid circles): (a) x,-x) plane of the ellipsoidal inclusion and (b) disk-shaped inclusion.
@I
Fig. 7. The equilibrium position of atoms as a prismatic dislocation loop embedded in a lattice: (a) dislocation loop in a infinite domain and (b) dislocation loop in a thin film. The shaded area represents the position of the dislocation loop.
Eigenstrain method for small defects in a lattice atoms in perfect lattice site. The eigendistortion, which is used to simulate the dislocation loop, is placed in the shaded region. 3.2. Small defects in a thin film A prismatic dislocation loop embedded in the thin film with nine layers of atoms is calculated. The layer labeled “0” is assumed to be the rigid boundary (displacement is zero) and the layer labeled “8” is assumed to be the free boundary. As explained before, the dislocation loop can be simulated by an disk-shaped inclusion; this disk-shaped inclusion is arbitrarily inserted between layer 4 and layer 5. The elastic constants, radius of the dislocation loop and the Burgers vector are the same as those used for the dislocation loop embedded in the infinite lattice. The boundary atoms can be identified in the same way as before. The result is obtained by solving for ii,(l) from eqn (52) first, then substitute t7,(() into eqn (39) and integrate numerically. The same integration technique is used here. The position of atoms after deformation is illustrated in Fig. 7(b). 4. DISCUSSION The eigenstrain Green function has been derived and used to calculate the distortion field caused by several kinds of defects with complicated shape. The most efficient way to solve the problem of defects in lattice is to construct the eigenstrain Green function table through the use of eqn (35). These results are calculated in the real space instead of the reciprocal space, and they can be used to obtain the distortion field directly through the use of eqn (38). The lattice model used in this work is just an illustration and the eigenstrain Green function can be derived in a similar way for different kinds of lattice structure. The distortion field obtained from the lattice calculation are compared to the results obtained from the elasticity calculation. The comparison indicates a good agreement as the shapes of the defects in lattice structure are similar to those in elasticity calculation. For example, the cubical shape defects in the simple cubic lattice can be matched exactly by the cubical inclusion and the distortion fields are very close to those calculated from the elasticity calculation (see Fig. 4). However, the distortion fields are different as the defect shapes prescribed by the continuous model are approximated for the other two cases, ellipsoidal and prismatic dislocation loop. They are shown in Figs 3 and 5. For example, the strain (slope of the displacement) calculated from the lattice theory is larger in the middle of the defect and smaller near the boundary when compared with the elasticity solution -(Fig. 3: ellipsoidal inclusion); the strain fields near the edge of dislocation loop are larger than those expected from the elasticity calculation, and this can be seen in Fig. S(b) that the displacement of atoms in the x2 direction at level 1 with the distance of 4b and
1029
66 from the center are significantly large. This difference decays for atoms away from the defects. This suggests that the exact shape of the defects is important for the displacement of the atoms in the vicinity of the defects. This implies that the stress concentration near the defect could be significantly different from the stress concentration calculated from the elasticity calculation. Especially for those at the kink edge of the defect domain, which is very true for the real lattice structural of the materials. For very small defects (e.g. point defects), the stress and strain fields are of great difference as compared to the elasticity calculation. It follows that the elastic strain energy can be different. This difference could be important for the nucleation of precipitates in solid, since the nucleus is sufficiently small. The results of prismatic dislocation loop embedded a thin film are significantly different from those of the same defect in an infinite lattice (see Fig. 7(a) and (b)). The displacements of atoms are no longer symmetrical to the plane of the dislocation loop which is true for this defect in an infinite lattice. The displacement of the atoms near the edge of dislocation in the radial direction is relatively large as compared to the case in infinite domain. The atoms between the defects and the substrate cannot move too much because of the rigid boundary condition. However, as more layers are put in between the dislocation loop and the rigid boundary, similar results are obtained. We conclude that the defect in a thin film with a rigid boundary can simulate the defect in a half-space lattice as sufficient layers are put in between the defect and rigid boundary. The disturbance of the defect is more significant for the atoms between the defect and the free surface. The difference of distortion fields suggests the importance of surface effect for the defects near the surface. Not only are the stress fields strongly influenced, but the elastic strain energies are of great difference. Furthermore, the elastic strain energy could change significantly as the dislocation loop moves toward the surface. Therefore, from the results of this work we are encouraged to calculate the stable criterion for a dislocation loop near the surface in the future. 5. CONCLUSION On the basis of this work we can conclude the following.
(1) Two Green functions have been derived. One is used for the case of distributed body force and the other one is used for the case of distributed eigenstrains. (2) The eigenstrain Green function is important and useful for defects with complicated shape, and the distortion fields can be calculated straightforward as the table of eigenstrain Green function is constructed.
K. C. KING and T. MURA
1030
(3) The distortion fields of dislocation type defects in
a lattice can be solved by this method directly, and this is not done before through the use of Kansaki force. (4) For small defects, the elasticity solution will give significant error especially near the boundary of the defect, and this local effect is important for This the analysis of stress concentration. difficulty can be resolved by using the lattice theory and eigenstrain method. (5) The effect of a defect in a thin film is calculated in this work, and the results are significantly different when compared with the same defect in an infinite lattice. This method is also applicable for different models of atomic structure. The potential of this method is that it can be further extended to calculate the elastic strain energy of the defect in a solid and the crack in a lattice can also be analyzed with the calculated equivalent strains. Acknowledgemenr-This research was supported by the U.S. ONR NOO014-80-K-0102 and also partially by the U.S.
AR0 DAAL03-89-K-0019.
REFERENCES I. Sato A. and Watanabe Y., J. Phys. Chem. Solids 49, 529 (1988).
2. Kuriyama M., J. Phy. Sot. Jap. 23, 1369 (1967). 3. Huang K., Proc. Roy. Sot. A190, 102 (1947). A. W. and Gorman R. L., Phys. Rev. 102, 4. Overhauser 676 (1956). 5. Kanzaki H., J. Phys. Chem. Solids 2, 24 (1957). 6. Tewary V. K., AERE HarweN Rep. TP, 388 (1969). 7. Tewary V. K., J. Phys. F: Metal Phys. 3, 1515 (1973). 8. Thomson R., Tewary V. K. and Masuda-Jindo K., J. Mar. Res. 2, 619 (1987). 9. Hsieh C. and Thomson R., J. appl. Phys. 44, 2051 (1973). 10. Slater J. C., Quantum Theory of Molecules and Solid. McGraw-Hill, New York (1974). 11. Morinaga M., Yukawa N., Adachi H., Mura T. and Sato A., Ann. Mig. Tokyo, Conf. Abstr. Japan Inst. Metals 272 (1986). 12. Mura T., Micromechanics of Defects in Solids. Martinus Nijhoff (1987). 13. Mura T., Eigensrrains in Lattice Theory l-503-519. University of Waterloo Press, Quebec (1977). 14. Born M. and Huang K., Dynamic Theory of Crystal Lattices. Oxford (1954). 15. Ashcroft N. W. and Mermin N. D., Solid State Physics. Cornell University, New York (1976). 16. Chiu Y. P., J. appl. Mech. 44, 587 (1977). 17. Kinoshita N. and Mura T., Phys. Stat. Sol. 5, 759 (1971). 18. Mura T. and Cheng P. C., J. appl. Mech. 44,591 (1977). 19. Eshelby J. D., Proc. Roy. Sot. A241, 376 (1957). 20. Kroupa F., Czech, J. Phys. BlO, 284 (1960). 21. Mindlin R. D., J. Elas. 2, 217 (1972).