- Email: [email protected]
Atomistic simulation of stress-induced phase transformation and recrystallization at the crack tip in bcc iron
Acta Materialia 55 (2007) 401–407 www.actamat-journals.com
Atomistic simulation of stress-induced phase transformation and recrystallization at the crack tip in bcc iron Ya-Fang Guo b
a,*
, Yue-Sheng Wang a, Dong-Liang Zhao
b
a Institute of Engineering Mechanics, Northern JiaoTong University, Beijing 100044, China Functional Materials Division, Central Iron and Steel Research Institute, Beijing 100081, China
Received 7 July 2006; received in revised form 8 August 2006; accepted 18 August 2006 Available online 1 November 2006
Abstract The mechanisms of low-temperature deformation at the crack tip in body-centered cubic (bcc) iron have been studied by molecular dynamics simulations. Phase transformation and recrystallization are found to occur at the crack tip when a sufficiently high stress concentration exists. The mechanisms of new grain nucleation and phase transformation in particular are discussed. For {1 1 0}Æ1 1 0æ and {1 11}Æ1 1 0æ cracks, phase transformation from the bcc structure to a typical close-packed hexagonal structure is observed accompanying new grain nucleation at the crack tip. It is found to be achieved by the emission of 16 h1 1 1i partial dislocations on the {1 1 2} slip plane. For the {1 0 0}Æ1 1 0æ crack, two long twin strips appear at the crack tip along the Æ1 1 1æ slip direction. A new grain with a different orientation from the original bcc crystal is formed with the growth of the twin region. 2006 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: Molecular dynamics simulations; Fracture; Phase transformation; Twinning; Iron
1. Introduction The fracture strength of metals depends strongly on the behavior of crack propagation, especially the deformation mechanisms at the crack tip region. In a previous work, deformation twinning was widely studied and expected to be an important mechanism in metals deformed at lowtemperature and/or high strain rates. For example, twin formation in body-centered cubic (bcc) crystals accompanying partial dislocation nucleation was observed in early experimental studies [1]. Recent experimental and computer simulation results have shown that deformation twinning may play a significant role as a mechanism of plasticity around crack tips [2–5]. Moreover, experimental results have also indicated that even face-centered cubic (fcc) materials that are not normally associated with defor-
*
Corresponding author. Fax: +86 10 51682094. E-mail address: [email protected] (Y.-F. Guo).
mation twinning, such as aluminum, will twin given a sufficiently high stress concentration such as at a crack tip [6]. For certain combinations of loading mode and orientation, deformation twinning occurring at aluminum crack tips was observed by atomistic simulation, and the mechanism of the twinning process at the atomistic scale was discussed [7]. As we have known, the structural phase transformation in steel is a widely observed phenomenon, and mechanisms of the transformation from the high-temperature austenitic phase of steel in the fcc structure to the low-temperature martensitic phase in the bcc structure are well discussed in the early studies [8,9]. Moreover, the pioneering experiments [10] discovered that the ground state crystal structure of bcc iron undergoes a stress-induced martensitic phase transformation to a hexagonally close-packed (hcp) structure, and the transformation pressure is specially discussed [11–14]. Afterwards, experimental results also revealed that for a certain range of pressure and temperature, the bcc iron
1359-6454/$30.00 2006 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.actamat.2006.08.022
402
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
was found to transform into a reversible fcc structure, and other observations showed nano grain formation associated with stress-induced phase transformation in iron [15,16]. Furthermore, multimillion-atom molecular dynamics simulations were used to investigate the shock-induced phase transformation of solid iron, and the results indicated that the dynamics and orientation of the developing close-packed grains depend on the shock strength and especially on the crystallographic shock direction [17]. In addition, a similar structural transition from the bcc structure to the close-packed hexagonal (hcp) structure was observed for cracks in iron [18]. More recently, twinning and stress-induced recrystallization phenomena were observed at the crack tip of both single crystal and nanocrystalline in bcc iron, and the deformation mechanisms in the crack tip region at low-temperatures are studied [19,20]. Based on these results, we can conclude that the structural phase transformation play an important role on the plastic deformation behavior of materials at certain loading and temperature conditions, and mechanisms of the phase transformation from bcc to hcp structures should be carefully examined. In the present paper, deformation mechanisms of the crack tip in bcc iron at low-temperatures and high loading levels are analyzed for different types of crack, while the orientations of the crack plane and the crack front are different. Twinning and recrystallization are observed at the crack tip region, and phase transformation from bcc to hcp structures are found accompanying new grain nucleation. We give the details of atomistic structures of both the new phase and the grain boundary. The mechanisms of the nucleation and growth of the new grain at the crack tip are discussed.
are all along the [1 0 1] direction, while the crack planes lie on (0 1 0), ð1 0 1Þ and ð1 1 1Þ plane, respectively. We consider the plane strain problem and apply the periodic boundary condition in the direction parallel to the crack front. The outermost two layers of atoms in the model are identified as the boundary region (except for the crack faces) in atomistic simulations, which is equal to the cut-off distance of the potential. The Finnis–Sinclair N-body potential [21] was used in this work. Before atomistic simulations, we first get the pre-existing crack in the model according to the continuum theory by shifting the atoms from their positions in a perfect crystal to positions specified by the anisotropic elasticity continuum mechanics equations [22] for the desired stress intensity KI. After getting the initial configurations of the crack at different loads, the atoms at the crack tip are fully relaxed for 4000 steps of magnitude 5 · 1015 s with the fixed-displacement boundary condition, and the structure evolution of the crack tip at different loading levels is investigated. The temperature of the system is invariant throughout the simulation, which is achieved by scaling all atoms’ instantaneous velocities with the appropriate Maxwell– Boltzmann distribution at a specified temperature. The XMD program is used in our atomistic simulations [23]. According to anisotropic linear elastic continuum theory, the displacement and stress field of a straight crack for the plane-strain condition [22] is given by
2. Conditions of computer simulations
and
In our simulations three types of mode I (opening mode) cracks in bcc iron crystal are considered, including types of {1 0 0}Æ1 1 0æ, {1 1 0}Æ1 1 0æ and {1 1 1}Æ1 1 0æ cracks. As shown in Fig. 1, crack fronts of these three types of cracks
Fig. 1. Lattice structure of bcc crystal.
fug ¼ K I ðrÞ1=2 ff ðhÞg frg ¼ K I ðrÞ
ð1=2Þ
ð1Þ
fgðhÞg
ð2Þ
where T
fug ¼ fux ; uy g ;
T
ff g ¼ ffx ; fy g ; T
¼ frx ; ry ; rxy g ;
K I ¼ rðapÞ
1=2
frg
fgg ¼ fgx ; gy ; gxy g
T
ð3Þ
ð4Þ
In Eqs. (1)–(4), KI is the stress intensity factor, r and h are the polar coordinates of the atom for which the displacement is calculated, r labels the distance from the atom to the crack tip, g is a function of the angle h between r and the crack plane, r is the applied load and a is the length of the crack. In our study we use the stress intensity factor KI to express the corresponding external load r, which can be defined as ‘‘nKIC’’, where KIC is the critical stress intensity factor. We tested several systems of different sizes beforehand, and verified that the basic results of the structure evolution at the crack tip in our simulation, such as the twinning and the phase transformation, are neither affected by the presence of the fixed boundary nor insensitive to the change of the system size if the simulated atomistic region is sufficiently large. In the present simulation, we select systems with dimensions of up to 35 · 35 by 1.1 nm, which contains about 138,000 atoms in total. The initial crack length in ˚. our calculation model is about 37 A
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
3. Results 3.1. Structure evolution of the crack tip for different types of crack Fig. 2 shows the crack configuration at 5 K at the load of 1.9 KIC for different types of cracks. In Fig. 2a–c for cracks with a Æ1 1 0æ crack front, both twinning and recrystallization at the crack tip are observed. For the {1 0 0}Æ1 1 0æ crack in Fig. 2a, we can see clearly that a new grain is formed at the crack tip region accompanying the formation of two long twin strips. The structure of the new grain is the same as that of the original bcc crystal, but the orientation is different. Thus no phase transformation occurs in the process of recrystallization for the {1 0 0}Æ1 1 0æ crack. However, for the {1 1 0}Æ1 1 0æ and {1 1 1}Æ1 1 0æ cracks in Fig. 2b and c, new grains are
403
formed at the crack tip with a typical hexagonal closepacked (hcp) structure, which is different from the original structure of the bcc crystal. This means that phase transformation occurs accompanied by new grain nucleation at the crack tip. Meanwhile, we can find twin formation at the crack tip in Fig. 2b and c, but the dimension of the twin region is rather small compared with that in Fig. 2a. Furthermore, we have obtained the crack configuration at different loads, and found that, with the load increasing, both the thickness of twinned region and the dimension of the new grain increase. The results also indicate that the stress-induced phase transformation or recrystallization can only occur when the load is higher than a value of about 1.8 KIC. Furthermore, crack branches are observed in Fig. 2a–c accompanied with new grain nucleation at the crack tip.
Fig. 2. Atomistic configuration of crack at 5 K under the load of 1.9 KIC for: (a) {1 0 0}Æ1 1 0æ crack; (b) {110}Æ1 1 0æ crack; and (c) {1 1 1}Æ1 1 0æ crack.
Fig. 3. Atomistic structure of the crack tip at 5 K under the load of 2.0 KIC. (a), (b) and (c) represent a {1 1 0}Æ1 1 0æ crack after 90, 110 and 140 steps relaxation, respectively; (d), (e) and (f) represent a {1 1 1}Æ1 1 0æ crack after 40, 80 and 140 steps relaxation, respectively. Solid and open circles indicate atoms in A and B layers.
404
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
3.2. Phase transformation at the crack tip for {1 1 0}Æ1 1 0æ and {1 1 1}Æ1 1 0æ cracks Fig. 3 gives the details of the atomistic process of phase transformation at the crack tip. In Fig. 3a–c for the {1 1 0}Æ1 1 0æ crack, we observe that atoms at the crack tip move along the Æ1 1 0æ direction after 90 time steps relaxation, which is perpendicular to the crack plane. Atoms in layer A move along the ½ 1 0 1 direction, while atoms in layer B move along the opposite ½1 0 1 direction. After 140 time steps relaxation, a unit of hcp structure is formed at the crack tip. Meanwhile, the crack tip branches appear in the ½ 1 1 1 and ½1 1 1 directions due to the formation of the hcp unit. For the {1 1 1}Æ1 1 0æ crack in Fig. 3d–f, the direction of atoms movement is just the same as that in the {1 1 0}Æ1 1 0æ crack, but is not perpendicular to the crack plane. After 140 time steps relaxation, a unit of hcp structure is generated at the crack tip just above the crack plane. In Fig. 4a we give the boundary structure of the new grain at the crack tip for the {1 1 0}Æ1 1 0æ crack. M and N represent ð 1 2 1Þ slip planes, while K and L represent ð 1 2 1Þ slip planes. We observe in Fig. 4a that the slip planes of M, N, K and L can move along the slip directions of ½1 1 1 and ½1 1 1, thus 16 ½ 1 1 1 and 16 ½1 1 1 partial dislocations are generated. The grain boundary is formed by the emission of 16 ½ 1 1 1 and 16 ½1 1 1 partial dislocations on successive ð 1 2 1Þ and ð 1 2 1Þ planes, which finally causes the growth of the new grain. After 1000 time steps relaxation, the diame˚. ter of the grain is up to approximately 40 A In Fig. 4b for the {1 1 1}Æ1 1 0æ crack, a new grain nucleated above the crack plane is found. The structure of the
grain boundary in the {1 1 1}Æ1 1 0æ crack is similar to that 1 1 1 in the {1 1 0}Æ1 1 0æ crack, and consists of a row of 16 ½ partial dislocations on the ð1 2 1Þ plane. Below the crack plane, a four layer twin is observed lying along the ½1 1 1 direction with a ð1 2 1Þ twin plane, and four 16 ½1 1 1 partial dislocations are emitted from the crack tip. After 400 time steps relaxation in Fig. 4b, the length of the nucleated twin ˚. is up to approximately 60 A We compare the atomic structure of bcc and hcp crystals in Fig. 5. It is found that the phase transformation from bcc to hcp structures is caused by atoms movement of 1 ½1 0 1 and 13 ½1 0 1 along the ½1 0 1 direction and the oppo3 site ½1 0 1 direction. After phase transformation, the distance of atoms in the [0 1 0] direction decreases from ˚ (a) to 2.733 A ˚ (a 0 ), while the distance of atoms 2.8665 A ˚ (b) to in the ½1 0 1 direction increases from 4.0525 A
Fig. 5. Phase transformation of: (a) the bcc structure to (b) the hcp structure. Solid and open circles indicate atoms in A and B layers, respectively.
Fig. 4. Details of the crack tip region at 5 K under the load of 2.0 KIC: (a) {1 1 0}Æ1 1 0æ crack after 1000 steps relaxation and (b) {1 1 1}Æ1 1 0æ crack after 400 steps relaxation. Solid and open circles indicate atoms in A and B layers.
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
405
˚ (b 0 ). This mechanism at an atomic point is consis4.7336 A tent with the particular transformation path proposed by Mao et al. [24]. Furthermore, from the analysis of the boundary structure in Fig. 4, we find that the new grain nucleates and grows up by the emission of 16 ½ 1 1 1 partial dislocations. Therefore, atoms movement of 13 ½1 0 1 can be achieved by two steps, the emission of a 16 ½1 1 1 partial dislocation on the ð1 2 1Þ slip plane and the subsequent emission of a 16 ½1 1 1 partial dislocation on the ð1 2 1Þ slip plane. This process can be described as the following dislocation reaction: a a a ½1 0 1 ! ½1 1 1 þ ½1 1 1 ð5Þ 3 6 6 3.3. Grain nucleation at the crack tip for the {1 0 0}Æ1 1 0æ crack For the {1 0 0}Æ1 1 0æ crack, no phase transformation occurs at the crack tip. The stress concentration at the crack tip is released by the nucleation of the new grain with a different orientation from the original bcc crystal. In Fig. 6b we find that 16 ½1 1 1 and 16 ½1 1 1 partial dislocations are nucleated separately on the ð1 2 1Þ and ð1 2 1Þ slip planes after 60 time steps relaxation, and it causes twin formation at the crack tip. After 140 time steps relaxation in Fig. 6c, the bcc unit at the crack tip is stretched in the [0 1 0] direction due to the formation of two twin strips separately along the ½1 1 1 and ½1 1 1 directions, which causes new grain nucleation at the crack tip. The boundary between the twin strip and the new grain is a Æ1 1 0æ symmetrical tilt boundary with a rotation angle of 70.5. Subsequently, new grain grows with the growth of two twin regions, and the lengths of the new grain in ½1 1 1 and ½ 1 1 1 directions are always equal to the thickness of the twin region. After 500 time steps relaxation in Fig. 7, two five-layered twins are formed at the crack tip both above and below the crack plane, and the diameter of the ˚ . The difference of new grain is up to approximately 20 A the angle between the new grain and the original bcc crystal is 90. Meanwhile, the length of the twin is up to approxi˚. mately 50 A
Fig. 7. Details of the crack tip region at 5 K for a {1 0 0}Æ1 1 0æ crack at 1.9 KIC after 500 steps relaxation. Solid and open circles indicate atoms in A and B layers.
3.4. Stress and energy analysis of phase transformation and recrystallization According to anisotropic linear elastic continuum theory, the shear stress of a straight crack for the plane-strain condition [22] is given in Eq. (2). In Fig. 8 we calculate the values of gxy(h), varying with the change of h from 0 to 360 for three types of cracks. For {1 0 0}Æ1 1 0æ and {1 1 0}Æ1 1 0æ cracks, we find that the distributions of the shear stress are symmetric, and the maximum values are at angles of about 72 (288) and 63 (297), respectively. Meanwhile, the value of gxy(h) for the {1 1 0}Æ1 1 0æ crack is much greater than that for the {1 0 0}Æ1 1 0æ crack. Thus phase transformation occurs for the {1 1 0}Æ1 1 0æ crack, while only twinning and recrystallization exist for the {1 0 0}Æ1 1 0æ crack. For the {1 1 1}Æ1 1 1æ crack, the distribution of the shear stress are unsymmetrical, and it causes the unsymmetrical behavior of the crack tip deformation. The maximum values exist at the angle of 81, where the phase transformation occurs above the crack plane. At the angle
Fig. 6. Atomistic structure of the crack tip at 5 K for a {1 0 0}Æ1 1 0æ crack at 1.9 KICafter: (a) 20 steps; (b) 60 steps; and (c) 140 steps relaxation. Solid and open circles indicate atoms in A and B layers.
406
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
Fig. 8. Distribution of shear stress around the crack tip. The dash line represents the results for a {1 0 0}Æ1 1 0æ crack; dash-dot line for {110}Æ1 1 0æ crack; and solid line for {1 1 1}Æ1 1 0æ crack.
of 279, the value of gxy(h) is rather smaller, thus a multilayer twin is formed below the crack plane. Based on the potential we used here, we obtained the energy of the bcc and hcp structures. The results indicate that the calculated energy of atoms in the hcp structure is ˚ , which 4.17 eV/atom at a lattice constant a of 2.733 A is 0.11 eV/atom higher than that of the bcc phase, with a ˚ . Thus, for the crack with a lattice constant a of 2.8665 A Æ1 1 0æ front, it is easier for stacking fault formation and twinning to occur compared with the phase transformation, and the phase transformation can only occur at a higher loading level. Moreover, the average energy of atoms in bcc and hcp crystals varies with the change of the lattice constant a. When a is below a value of about ˚ , the atom energy of the hcp crystal becomes lower 2.74 A than that of the bcc crystal. This is why the new grain with an hcp structure can nucleate and exist at the crack tip when a high loading level is applied. 4. Discussion From Section 3, we find that twinning and recrystallization are the main deformation mechanisms at the crack tip at low-temperatures. When the crack is loaded at a higher loading level, the crack tip region is subjected to high stresses. There are two ways to release the stress concentrated at the crack tip. One is to change the crystal structure, which causes the stress-induced phase transformation at the crack tip. The other way is to change the crystal orientation, which causes new grain nucleation or twinning. After new grain nucleation (with or without phase transformation processes) and/or twinning, the high stress concentration at the crack tip region is released, and the energy is
transformed into boundary energy, phase transition energy and elastic strain energy owing to twinning and new grain formation. Experimental results have revealed that twinning is a main deformation mechanism at low-temperatures and plays a significant role in plasticity deformation around crack tips [1–3]. In addition, computer simulations are widely used to study the twining process at the crack tip in bcc iron [4,5]. Moreover, phase transformation from bcc to hcp iron under a special stress or temperature condition is a well observed phenomenon, and many experiments have shown the occurrence of a stress-induced martensitic phase transformation from bcc to hcp structures [10–14]. Furthermore, martensitic phase transformations have also been observed in creep damage and fatigue experiments, accompanied by crack initiation and growth [25]. In an early study of atomistic simulation, the pressureinduced martensitic phase transformation from bcc to hcp structures in a perfect-crystal iron was investigated using different interatomic potentials [26]. Subsequently, the same structural transitions in a bcc crystal under a uniaxial stress [27] and a shocking load [17] were obtained by atomistic simulations. The crack tip is a region that undergoes a rather high stress concentration, thus stress-induced phase transformation and recrystallization can occur, which has also been proved recently by computational simulations by using an EAM interatomic potential [20]. The Finnis–Sinclair N-body potential we applied in this work is widely used to study the deformation behavior of bcc crystals, such as dislocation nucleation and emission, twinning, fracture and plasticity [5,28,29], including the exploration of the pressure dependence of the deformation response of bcc crystals [29]. Furthermore, atomistic simulations of the stacking fault formation in bcc iron using the Finnis–Sinclair N-body potential in particular have been discussed [30], confirming that the Finnis–Sinclair potential is efficient for the calculation of deformation behavior in bcc crystal. Moreover, the structural transition from bcc to hcp iron has been obtained using different interatomic potentials both in a perfect-crystal iron and the crack tip region [20,26,27]. Thus we are sure that the phenomenon of phase transformation and recrystallization at the crack tip obtained by our simulation is credible and does not depend on the special interatomic potential we are used. 5. Conclusions In this paper, molecular dynamics simulations are used to investigate the mechanisms of low-temperature deformation at the crack tip in bcc iron, and the mechanisms of new grain formation and the phase transformation at the crack tip in particular are studied. Atomistic simulation shows that the phase transformation from bcc to hcp structures at the crack tip can be obtained by the successive emission of 16 h1 1 1i partial dislocations on the {1 1 2} slip plane, while the atoms’ distance along the [0 1 0] direction contracted. This mechanism is exactly consistent with the par-
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
ticular transformation path proposed by previous work based on the analysis of crystal structures [24], and it also agrees well with the atomistic simulation results of the martensitic phase transformation in a perfect-crystal iron [26]. As we have known, the martensitic transformation depends on the special stress and temperature conditions, but the effect of temperature is not considered in the present work. Thus, in future work, the effect of the temperature on the phase transformation should be carefully explored, especially the mechanism of transformation at the atomistic scale. Moreover, for other bcc structure materials, such as NiAl alloys with a B2-type structure, martensitic transformation has been observed [31]. It would be useful to study the deformation behavior at the crack tip in NiAl alloys by atomistic simulation, especially the mechanism of phase transformation, so that we can compare it with the mechanism of martenstic transformation at the crack tip in bcc iron. Acknowledgements The research was supported by Chinese Nature Science Foundation (Grant No. 10572019 and No. 10672016) and NJTU Science Foundation (No. 2005SM0035). References [1] [2] [3] [4] [5] [6] [7]
Ogawa K. Phil Mag 1965;11:217. M} ullner P. Mater Sci Eng A 1997;94:234. Smida T, Bosansky J. Kovove Mater 2002;40:146. Wei X, Moriarty JA. Phys Rev B 1996;54:6941. Guo YF, Wang CY, Zhao DL. Mater Sci Eng A 2003;349:29. Chen Q, Huang Y, Qiao L, Chu W. Sci China (Ser. E) 1999;42:1. Hai S, Tadmor EB. Acta Mater 2003;51:117.
407
[8] Bain EC, Dunkirk NY. The nature of martensite. Trans AIME 1924;70:25. [9] Wayman CM. The nature of martensite phase transformations nondiffusive. In: Cahn RW, Haasen P, editors. Physical metallurgy. Amsterdam: North-Holland; 1983. p. 1031 [chapter 15]. [10] Bancroft D, Peterson E, Minshall S. J Appl Phys 1956;27:291. [11] Rueff JP, Krisch M, Cai YQ, Kaprolat A, Hanfland M, Lorenzen M, et al. Phys Rev B 1999;60:14510. [12] Barker L, Hollenbach R. J Appl Phys 1974;45:4872. [13] Taylor R, Pasternak M, Jeanloz R. J Appl Phys 1991;69:6126. [14] Von Bargen N, Boehler R. High Pressure Res 1990;6:133. [15] Fujiwara H, Inomoto H, Sanada R, Ameyama K. Scripta Mater 2001;44:2039. [16] Zhang J, Guyot F. Phys Chem Miner 1999;26:419. [17] Kadau K, Germann TC, Lohmdahl PS, Holian BL. Science 2002;296:1681. [18] Nishimura K, Miyazaki N. Comput Model Eng Sci 2001;2:143. [19] Farkas D. Phil Mag 2005;85:387. [20] Latapie L, Farkas D. Modelling Simul Mater Sci Eng 2003;11: 745. [21] (a) Finnis MW, Sinclair JE. Phil Mag A 1984;50:45; (b) Finnis MW, Sinclair JE. Phil Mag A 1986;53:161. [22] Sih GC, Liebowitz H. Mathematical theories of brittle facture. In: Fracture: an advanced treatise. New York: Academic Press; 1968. [23] Riffkin J, Center for simulation, University of Connecticut, CT. Available from: http://www.ims.uconn.edu/centers/simul/. [24] Mao H, Bassett W, Takahashi T. J Appl Phys 1967;38:272. [25] Nagae Y. Mater Sci Eng A 2004;387–389:665. [26] Lee KY, Ray JR. Phys Rev B 1989;39:565. [27] Cheung KS, Harrison RJ, Yip S. J Appl Phys 1992;71:4009. [28] Chang J, Cai W, Bulatov VV, Yip S. Mater Sci Eng A 2001;309– 310:160. [29] Xu DS, Chang JP, Li J, Yang R, Li D, Yip S. Mater Sci Eng A 2004;387–389:840. [30] Machova´ A, Beltz GE, Chang M. Modelling Simul Mater Sci Eng 1999;7:949. [31] Lasˇek J, Chra´ska T, Krˇecˇek P, Bartusˇka P. Scripta Mater 1997;37:897.
Atomistic simulation of stress-induced phase transformation and recrystallization at the crack tip in bcc iron Ya-Fang Guo b
a,*
, Yue-Sheng Wang a, Dong-Liang Zhao
b
a Institute of Engineering Mechanics, Northern JiaoTong University, Beijing 100044, China Functional Materials Division, Central Iron and Steel Research Institute, Beijing 100081, China
Received 7 July 2006; received in revised form 8 August 2006; accepted 18 August 2006 Available online 1 November 2006
Abstract The mechanisms of low-temperature deformation at the crack tip in body-centered cubic (bcc) iron have been studied by molecular dynamics simulations. Phase transformation and recrystallization are found to occur at the crack tip when a sufficiently high stress concentration exists. The mechanisms of new grain nucleation and phase transformation in particular are discussed. For {1 1 0}Æ1 1 0æ and {1 11}Æ1 1 0æ cracks, phase transformation from the bcc structure to a typical close-packed hexagonal structure is observed accompanying new grain nucleation at the crack tip. It is found to be achieved by the emission of 16 h1 1 1i partial dislocations on the {1 1 2} slip plane. For the {1 0 0}Æ1 1 0æ crack, two long twin strips appear at the crack tip along the Æ1 1 1æ slip direction. A new grain with a different orientation from the original bcc crystal is formed with the growth of the twin region. 2006 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: Molecular dynamics simulations; Fracture; Phase transformation; Twinning; Iron
1. Introduction The fracture strength of metals depends strongly on the behavior of crack propagation, especially the deformation mechanisms at the crack tip region. In a previous work, deformation twinning was widely studied and expected to be an important mechanism in metals deformed at lowtemperature and/or high strain rates. For example, twin formation in body-centered cubic (bcc) crystals accompanying partial dislocation nucleation was observed in early experimental studies [1]. Recent experimental and computer simulation results have shown that deformation twinning may play a significant role as a mechanism of plasticity around crack tips [2–5]. Moreover, experimental results have also indicated that even face-centered cubic (fcc) materials that are not normally associated with defor-
*
Corresponding author. Fax: +86 10 51682094. E-mail address: [email protected] (Y.-F. Guo).
mation twinning, such as aluminum, will twin given a sufficiently high stress concentration such as at a crack tip [6]. For certain combinations of loading mode and orientation, deformation twinning occurring at aluminum crack tips was observed by atomistic simulation, and the mechanism of the twinning process at the atomistic scale was discussed [7]. As we have known, the structural phase transformation in steel is a widely observed phenomenon, and mechanisms of the transformation from the high-temperature austenitic phase of steel in the fcc structure to the low-temperature martensitic phase in the bcc structure are well discussed in the early studies [8,9]. Moreover, the pioneering experiments [10] discovered that the ground state crystal structure of bcc iron undergoes a stress-induced martensitic phase transformation to a hexagonally close-packed (hcp) structure, and the transformation pressure is specially discussed [11–14]. Afterwards, experimental results also revealed that for a certain range of pressure and temperature, the bcc iron
1359-6454/$30.00 2006 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.actamat.2006.08.022
402
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
was found to transform into a reversible fcc structure, and other observations showed nano grain formation associated with stress-induced phase transformation in iron [15,16]. Furthermore, multimillion-atom molecular dynamics simulations were used to investigate the shock-induced phase transformation of solid iron, and the results indicated that the dynamics and orientation of the developing close-packed grains depend on the shock strength and especially on the crystallographic shock direction [17]. In addition, a similar structural transition from the bcc structure to the close-packed hexagonal (hcp) structure was observed for cracks in iron [18]. More recently, twinning and stress-induced recrystallization phenomena were observed at the crack tip of both single crystal and nanocrystalline in bcc iron, and the deformation mechanisms in the crack tip region at low-temperatures are studied [19,20]. Based on these results, we can conclude that the structural phase transformation play an important role on the plastic deformation behavior of materials at certain loading and temperature conditions, and mechanisms of the phase transformation from bcc to hcp structures should be carefully examined. In the present paper, deformation mechanisms of the crack tip in bcc iron at low-temperatures and high loading levels are analyzed for different types of crack, while the orientations of the crack plane and the crack front are different. Twinning and recrystallization are observed at the crack tip region, and phase transformation from bcc to hcp structures are found accompanying new grain nucleation. We give the details of atomistic structures of both the new phase and the grain boundary. The mechanisms of the nucleation and growth of the new grain at the crack tip are discussed.
are all along the [1 0 1] direction, while the crack planes lie on (0 1 0), ð1 0 1Þ and ð1 1 1Þ plane, respectively. We consider the plane strain problem and apply the periodic boundary condition in the direction parallel to the crack front. The outermost two layers of atoms in the model are identified as the boundary region (except for the crack faces) in atomistic simulations, which is equal to the cut-off distance of the potential. The Finnis–Sinclair N-body potential [21] was used in this work. Before atomistic simulations, we first get the pre-existing crack in the model according to the continuum theory by shifting the atoms from their positions in a perfect crystal to positions specified by the anisotropic elasticity continuum mechanics equations [22] for the desired stress intensity KI. After getting the initial configurations of the crack at different loads, the atoms at the crack tip are fully relaxed for 4000 steps of magnitude 5 · 1015 s with the fixed-displacement boundary condition, and the structure evolution of the crack tip at different loading levels is investigated. The temperature of the system is invariant throughout the simulation, which is achieved by scaling all atoms’ instantaneous velocities with the appropriate Maxwell– Boltzmann distribution at a specified temperature. The XMD program is used in our atomistic simulations [23]. According to anisotropic linear elastic continuum theory, the displacement and stress field of a straight crack for the plane-strain condition [22] is given by
2. Conditions of computer simulations
and
In our simulations three types of mode I (opening mode) cracks in bcc iron crystal are considered, including types of {1 0 0}Æ1 1 0æ, {1 1 0}Æ1 1 0æ and {1 1 1}Æ1 1 0æ cracks. As shown in Fig. 1, crack fronts of these three types of cracks
Fig. 1. Lattice structure of bcc crystal.
fug ¼ K I ðrÞ1=2 ff ðhÞg frg ¼ K I ðrÞ
ð1=2Þ
ð1Þ
fgðhÞg
ð2Þ
where T
fug ¼ fux ; uy g ;
T
ff g ¼ ffx ; fy g ; T
¼ frx ; ry ; rxy g ;
K I ¼ rðapÞ
1=2
frg
fgg ¼ fgx ; gy ; gxy g
T
ð3Þ
ð4Þ
In Eqs. (1)–(4), KI is the stress intensity factor, r and h are the polar coordinates of the atom for which the displacement is calculated, r labels the distance from the atom to the crack tip, g is a function of the angle h between r and the crack plane, r is the applied load and a is the length of the crack. In our study we use the stress intensity factor KI to express the corresponding external load r, which can be defined as ‘‘nKIC’’, where KIC is the critical stress intensity factor. We tested several systems of different sizes beforehand, and verified that the basic results of the structure evolution at the crack tip in our simulation, such as the twinning and the phase transformation, are neither affected by the presence of the fixed boundary nor insensitive to the change of the system size if the simulated atomistic region is sufficiently large. In the present simulation, we select systems with dimensions of up to 35 · 35 by 1.1 nm, which contains about 138,000 atoms in total. The initial crack length in ˚. our calculation model is about 37 A
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
3. Results 3.1. Structure evolution of the crack tip for different types of crack Fig. 2 shows the crack configuration at 5 K at the load of 1.9 KIC for different types of cracks. In Fig. 2a–c for cracks with a Æ1 1 0æ crack front, both twinning and recrystallization at the crack tip are observed. For the {1 0 0}Æ1 1 0æ crack in Fig. 2a, we can see clearly that a new grain is formed at the crack tip region accompanying the formation of two long twin strips. The structure of the new grain is the same as that of the original bcc crystal, but the orientation is different. Thus no phase transformation occurs in the process of recrystallization for the {1 0 0}Æ1 1 0æ crack. However, for the {1 1 0}Æ1 1 0æ and {1 1 1}Æ1 1 0æ cracks in Fig. 2b and c, new grains are
403
formed at the crack tip with a typical hexagonal closepacked (hcp) structure, which is different from the original structure of the bcc crystal. This means that phase transformation occurs accompanied by new grain nucleation at the crack tip. Meanwhile, we can find twin formation at the crack tip in Fig. 2b and c, but the dimension of the twin region is rather small compared with that in Fig. 2a. Furthermore, we have obtained the crack configuration at different loads, and found that, with the load increasing, both the thickness of twinned region and the dimension of the new grain increase. The results also indicate that the stress-induced phase transformation or recrystallization can only occur when the load is higher than a value of about 1.8 KIC. Furthermore, crack branches are observed in Fig. 2a–c accompanied with new grain nucleation at the crack tip.
Fig. 2. Atomistic configuration of crack at 5 K under the load of 1.9 KIC for: (a) {1 0 0}Æ1 1 0æ crack; (b) {110}Æ1 1 0æ crack; and (c) {1 1 1}Æ1 1 0æ crack.
Fig. 3. Atomistic structure of the crack tip at 5 K under the load of 2.0 KIC. (a), (b) and (c) represent a {1 1 0}Æ1 1 0æ crack after 90, 110 and 140 steps relaxation, respectively; (d), (e) and (f) represent a {1 1 1}Æ1 1 0æ crack after 40, 80 and 140 steps relaxation, respectively. Solid and open circles indicate atoms in A and B layers.
404
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
3.2. Phase transformation at the crack tip for {1 1 0}Æ1 1 0æ and {1 1 1}Æ1 1 0æ cracks Fig. 3 gives the details of the atomistic process of phase transformation at the crack tip. In Fig. 3a–c for the {1 1 0}Æ1 1 0æ crack, we observe that atoms at the crack tip move along the Æ1 1 0æ direction after 90 time steps relaxation, which is perpendicular to the crack plane. Atoms in layer A move along the ½ 1 0 1 direction, while atoms in layer B move along the opposite ½1 0 1 direction. After 140 time steps relaxation, a unit of hcp structure is formed at the crack tip. Meanwhile, the crack tip branches appear in the ½ 1 1 1 and ½1 1 1 directions due to the formation of the hcp unit. For the {1 1 1}Æ1 1 0æ crack in Fig. 3d–f, the direction of atoms movement is just the same as that in the {1 1 0}Æ1 1 0æ crack, but is not perpendicular to the crack plane. After 140 time steps relaxation, a unit of hcp structure is generated at the crack tip just above the crack plane. In Fig. 4a we give the boundary structure of the new grain at the crack tip for the {1 1 0}Æ1 1 0æ crack. M and N represent ð 1 2 1Þ slip planes, while K and L represent ð 1 2 1Þ slip planes. We observe in Fig. 4a that the slip planes of M, N, K and L can move along the slip directions of ½1 1 1 and ½1 1 1, thus 16 ½ 1 1 1 and 16 ½1 1 1 partial dislocations are generated. The grain boundary is formed by the emission of 16 ½ 1 1 1 and 16 ½1 1 1 partial dislocations on successive ð 1 2 1Þ and ð 1 2 1Þ planes, which finally causes the growth of the new grain. After 1000 time steps relaxation, the diame˚. ter of the grain is up to approximately 40 A In Fig. 4b for the {1 1 1}Æ1 1 0æ crack, a new grain nucleated above the crack plane is found. The structure of the
grain boundary in the {1 1 1}Æ1 1 0æ crack is similar to that 1 1 1 in the {1 1 0}Æ1 1 0æ crack, and consists of a row of 16 ½ partial dislocations on the ð1 2 1Þ plane. Below the crack plane, a four layer twin is observed lying along the ½1 1 1 direction with a ð1 2 1Þ twin plane, and four 16 ½1 1 1 partial dislocations are emitted from the crack tip. After 400 time steps relaxation in Fig. 4b, the length of the nucleated twin ˚. is up to approximately 60 A We compare the atomic structure of bcc and hcp crystals in Fig. 5. It is found that the phase transformation from bcc to hcp structures is caused by atoms movement of 1 ½1 0 1 and 13 ½1 0 1 along the ½1 0 1 direction and the oppo3 site ½1 0 1 direction. After phase transformation, the distance of atoms in the [0 1 0] direction decreases from ˚ (a) to 2.733 A ˚ (a 0 ), while the distance of atoms 2.8665 A ˚ (b) to in the ½1 0 1 direction increases from 4.0525 A
Fig. 5. Phase transformation of: (a) the bcc structure to (b) the hcp structure. Solid and open circles indicate atoms in A and B layers, respectively.
Fig. 4. Details of the crack tip region at 5 K under the load of 2.0 KIC: (a) {1 1 0}Æ1 1 0æ crack after 1000 steps relaxation and (b) {1 1 1}Æ1 1 0æ crack after 400 steps relaxation. Solid and open circles indicate atoms in A and B layers.
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
405
˚ (b 0 ). This mechanism at an atomic point is consis4.7336 A tent with the particular transformation path proposed by Mao et al. [24]. Furthermore, from the analysis of the boundary structure in Fig. 4, we find that the new grain nucleates and grows up by the emission of 16 ½ 1 1 1 partial dislocations. Therefore, atoms movement of 13 ½1 0 1 can be achieved by two steps, the emission of a 16 ½1 1 1 partial dislocation on the ð1 2 1Þ slip plane and the subsequent emission of a 16 ½1 1 1 partial dislocation on the ð1 2 1Þ slip plane. This process can be described as the following dislocation reaction: a a a ½1 0 1 ! ½1 1 1 þ ½1 1 1 ð5Þ 3 6 6 3.3. Grain nucleation at the crack tip for the {1 0 0}Æ1 1 0æ crack For the {1 0 0}Æ1 1 0æ crack, no phase transformation occurs at the crack tip. The stress concentration at the crack tip is released by the nucleation of the new grain with a different orientation from the original bcc crystal. In Fig. 6b we find that 16 ½1 1 1 and 16 ½1 1 1 partial dislocations are nucleated separately on the ð1 2 1Þ and ð1 2 1Þ slip planes after 60 time steps relaxation, and it causes twin formation at the crack tip. After 140 time steps relaxation in Fig. 6c, the bcc unit at the crack tip is stretched in the [0 1 0] direction due to the formation of two twin strips separately along the ½1 1 1 and ½1 1 1 directions, which causes new grain nucleation at the crack tip. The boundary between the twin strip and the new grain is a Æ1 1 0æ symmetrical tilt boundary with a rotation angle of 70.5. Subsequently, new grain grows with the growth of two twin regions, and the lengths of the new grain in ½1 1 1 and ½ 1 1 1 directions are always equal to the thickness of the twin region. After 500 time steps relaxation in Fig. 7, two five-layered twins are formed at the crack tip both above and below the crack plane, and the diameter of the ˚ . The difference of new grain is up to approximately 20 A the angle between the new grain and the original bcc crystal is 90. Meanwhile, the length of the twin is up to approxi˚. mately 50 A
Fig. 7. Details of the crack tip region at 5 K for a {1 0 0}Æ1 1 0æ crack at 1.9 KIC after 500 steps relaxation. Solid and open circles indicate atoms in A and B layers.
3.4. Stress and energy analysis of phase transformation and recrystallization According to anisotropic linear elastic continuum theory, the shear stress of a straight crack for the plane-strain condition [22] is given in Eq. (2). In Fig. 8 we calculate the values of gxy(h), varying with the change of h from 0 to 360 for three types of cracks. For {1 0 0}Æ1 1 0æ and {1 1 0}Æ1 1 0æ cracks, we find that the distributions of the shear stress are symmetric, and the maximum values are at angles of about 72 (288) and 63 (297), respectively. Meanwhile, the value of gxy(h) for the {1 1 0}Æ1 1 0æ crack is much greater than that for the {1 0 0}Æ1 1 0æ crack. Thus phase transformation occurs for the {1 1 0}Æ1 1 0æ crack, while only twinning and recrystallization exist for the {1 0 0}Æ1 1 0æ crack. For the {1 1 1}Æ1 1 1æ crack, the distribution of the shear stress are unsymmetrical, and it causes the unsymmetrical behavior of the crack tip deformation. The maximum values exist at the angle of 81, where the phase transformation occurs above the crack plane. At the angle
Fig. 6. Atomistic structure of the crack tip at 5 K for a {1 0 0}Æ1 1 0æ crack at 1.9 KICafter: (a) 20 steps; (b) 60 steps; and (c) 140 steps relaxation. Solid and open circles indicate atoms in A and B layers.
406
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
Fig. 8. Distribution of shear stress around the crack tip. The dash line represents the results for a {1 0 0}Æ1 1 0æ crack; dash-dot line for {110}Æ1 1 0æ crack; and solid line for {1 1 1}Æ1 1 0æ crack.
of 279, the value of gxy(h) is rather smaller, thus a multilayer twin is formed below the crack plane. Based on the potential we used here, we obtained the energy of the bcc and hcp structures. The results indicate that the calculated energy of atoms in the hcp structure is ˚ , which 4.17 eV/atom at a lattice constant a of 2.733 A is 0.11 eV/atom higher than that of the bcc phase, with a ˚ . Thus, for the crack with a lattice constant a of 2.8665 A Æ1 1 0æ front, it is easier for stacking fault formation and twinning to occur compared with the phase transformation, and the phase transformation can only occur at a higher loading level. Moreover, the average energy of atoms in bcc and hcp crystals varies with the change of the lattice constant a. When a is below a value of about ˚ , the atom energy of the hcp crystal becomes lower 2.74 A than that of the bcc crystal. This is why the new grain with an hcp structure can nucleate and exist at the crack tip when a high loading level is applied. 4. Discussion From Section 3, we find that twinning and recrystallization are the main deformation mechanisms at the crack tip at low-temperatures. When the crack is loaded at a higher loading level, the crack tip region is subjected to high stresses. There are two ways to release the stress concentrated at the crack tip. One is to change the crystal structure, which causes the stress-induced phase transformation at the crack tip. The other way is to change the crystal orientation, which causes new grain nucleation or twinning. After new grain nucleation (with or without phase transformation processes) and/or twinning, the high stress concentration at the crack tip region is released, and the energy is
transformed into boundary energy, phase transition energy and elastic strain energy owing to twinning and new grain formation. Experimental results have revealed that twinning is a main deformation mechanism at low-temperatures and plays a significant role in plasticity deformation around crack tips [1–3]. In addition, computer simulations are widely used to study the twining process at the crack tip in bcc iron [4,5]. Moreover, phase transformation from bcc to hcp iron under a special stress or temperature condition is a well observed phenomenon, and many experiments have shown the occurrence of a stress-induced martensitic phase transformation from bcc to hcp structures [10–14]. Furthermore, martensitic phase transformations have also been observed in creep damage and fatigue experiments, accompanied by crack initiation and growth [25]. In an early study of atomistic simulation, the pressureinduced martensitic phase transformation from bcc to hcp structures in a perfect-crystal iron was investigated using different interatomic potentials [26]. Subsequently, the same structural transitions in a bcc crystal under a uniaxial stress [27] and a shocking load [17] were obtained by atomistic simulations. The crack tip is a region that undergoes a rather high stress concentration, thus stress-induced phase transformation and recrystallization can occur, which has also been proved recently by computational simulations by using an EAM interatomic potential [20]. The Finnis–Sinclair N-body potential we applied in this work is widely used to study the deformation behavior of bcc crystals, such as dislocation nucleation and emission, twinning, fracture and plasticity [5,28,29], including the exploration of the pressure dependence of the deformation response of bcc crystals [29]. Furthermore, atomistic simulations of the stacking fault formation in bcc iron using the Finnis–Sinclair N-body potential in particular have been discussed [30], confirming that the Finnis–Sinclair potential is efficient for the calculation of deformation behavior in bcc crystal. Moreover, the structural transition from bcc to hcp iron has been obtained using different interatomic potentials both in a perfect-crystal iron and the crack tip region [20,26,27]. Thus we are sure that the phenomenon of phase transformation and recrystallization at the crack tip obtained by our simulation is credible and does not depend on the special interatomic potential we are used. 5. Conclusions In this paper, molecular dynamics simulations are used to investigate the mechanisms of low-temperature deformation at the crack tip in bcc iron, and the mechanisms of new grain formation and the phase transformation at the crack tip in particular are studied. Atomistic simulation shows that the phase transformation from bcc to hcp structures at the crack tip can be obtained by the successive emission of 16 h1 1 1i partial dislocations on the {1 1 2} slip plane, while the atoms’ distance along the [0 1 0] direction contracted. This mechanism is exactly consistent with the par-
Y.-F. Guo et al. / Acta Materialia 55 (2007) 401–407
ticular transformation path proposed by previous work based on the analysis of crystal structures [24], and it also agrees well with the atomistic simulation results of the martensitic phase transformation in a perfect-crystal iron [26]. As we have known, the martensitic transformation depends on the special stress and temperature conditions, but the effect of temperature is not considered in the present work. Thus, in future work, the effect of the temperature on the phase transformation should be carefully explored, especially the mechanism of transformation at the atomistic scale. Moreover, for other bcc structure materials, such as NiAl alloys with a B2-type structure, martensitic transformation has been observed [31]. It would be useful to study the deformation behavior at the crack tip in NiAl alloys by atomistic simulation, especially the mechanism of phase transformation, so that we can compare it with the mechanism of martenstic transformation at the crack tip in bcc iron. Acknowledgements The research was supported by Chinese Nature Science Foundation (Grant No. 10572019 and No. 10672016) and NJTU Science Foundation (No. 2005SM0035). References [1] [2] [3] [4] [5] [6] [7]
Ogawa K. Phil Mag 1965;11:217. M} ullner P. Mater Sci Eng A 1997;94:234. Smida T, Bosansky J. Kovove Mater 2002;40:146. Wei X, Moriarty JA. Phys Rev B 1996;54:6941. Guo YF, Wang CY, Zhao DL. Mater Sci Eng A 2003;349:29. Chen Q, Huang Y, Qiao L, Chu W. Sci China (Ser. E) 1999;42:1. Hai S, Tadmor EB. Acta Mater 2003;51:117.
407
[8] Bain EC, Dunkirk NY. The nature of martensite. Trans AIME 1924;70:25. [9] Wayman CM. The nature of martensite phase transformations nondiffusive. In: Cahn RW, Haasen P, editors. Physical metallurgy. Amsterdam: North-Holland; 1983. p. 1031 [chapter 15]. [10] Bancroft D, Peterson E, Minshall S. J Appl Phys 1956;27:291. [11] Rueff JP, Krisch M, Cai YQ, Kaprolat A, Hanfland M, Lorenzen M, et al. Phys Rev B 1999;60:14510. [12] Barker L, Hollenbach R. J Appl Phys 1974;45:4872. [13] Taylor R, Pasternak M, Jeanloz R. J Appl Phys 1991;69:6126. [14] Von Bargen N, Boehler R. High Pressure Res 1990;6:133. [15] Fujiwara H, Inomoto H, Sanada R, Ameyama K. Scripta Mater 2001;44:2039. [16] Zhang J, Guyot F. Phys Chem Miner 1999;26:419. [17] Kadau K, Germann TC, Lohmdahl PS, Holian BL. Science 2002;296:1681. [18] Nishimura K, Miyazaki N. Comput Model Eng Sci 2001;2:143. [19] Farkas D. Phil Mag 2005;85:387. [20] Latapie L, Farkas D. Modelling Simul Mater Sci Eng 2003;11: 745. [21] (a) Finnis MW, Sinclair JE. Phil Mag A 1984;50:45; (b) Finnis MW, Sinclair JE. Phil Mag A 1986;53:161. [22] Sih GC, Liebowitz H. Mathematical theories of brittle facture. In: Fracture: an advanced treatise. New York: Academic Press; 1968. [23] Riffkin J, Center for simulation, University of Connecticut, CT. Available from: http://www.ims.uconn.edu/centers/simul/. [24] Mao H, Bassett W, Takahashi T. J Appl Phys 1967;38:272. [25] Nagae Y. Mater Sci Eng A 2004;387–389:665. [26] Lee KY, Ray JR. Phys Rev B 1989;39:565. [27] Cheung KS, Harrison RJ, Yip S. J Appl Phys 1992;71:4009. [28] Chang J, Cai W, Bulatov VV, Yip S. Mater Sci Eng A 2001;309– 310:160. [29] Xu DS, Chang JP, Li J, Yang R, Li D, Yip S. Mater Sci Eng A 2004;387–389:840. [30] Machova´ A, Beltz GE, Chang M. Modelling Simul Mater Sci Eng 1999;7:949. [31] Lasˇek J, Chra´ska T, Krˇecˇek P, Bartusˇka P. Scripta Mater 1997;37:897.