Molecular dynamics study of crack processes associated with dislocation nucleated at the tip

Materials Science and Engineering, A 176 (1994) 263-269

263

Molecular dynamics study of crack processes associated with dislocation nucleated at the tip Hiroshi Kitagawa, Akihiro Nakatani and Yoji Shibutani Department of Mechanical Engineering, Osaka University, 2-1 Yamada-Oka, Suita, Osaka 565 (Japan)

Abstract Atomic scale changes of structure around a crack-tip in an f.c.c, crystal under in-plane shear (mode II) loading are analyzed by molecular dynamics simulation (MD simulation), and its counterpart in the continuum crystal plasticity model is discussed. A ductile fracture process involving dislocation nucleation from the stressed tip is always observed. The nucleated dislocation is driven away from the crack and a dislocation-free zone develops in the near-tip region. Time averages of local stress in the near-tip region, both before and after dislocation nucleation, coincide well with the linear elastic prediction. The critical stress intensity factor estimated by MD simulation agrees well with Rice's theoretical prediction derived from the unstable stacking energy concept. The temperature dependence of the critical factor can be explained primarily as a thermally activated process.

1. Introduction The fracture process around the crack tip is sensitive to the local atomic configuration and mechanical action on the crack due to the applied load. One of the elemental processes which appears in the atomic structure around a crack, which is followed by ductile fracture of crystalline material, is dislocation nucleation at the tip and emission of the crack out of this. Several computer simulations based on molecular dynamics simulation (MD simulation) have been carried out to study such an atomistic process of cracking, and details of the atomic structure and the dynamics appearing around the crack tip have gradually been explained (see for example ref. 1 ). In a previous work [2], it was shown that in mode I cracking in an f.c.c.-crystal, crack tip deformation accompanied by dislocation activation is not necessary and is rather rare. A typical rare case is a crack front with [101] direction and surface on the (010) plane [2]. A ductile crack opening with tip blunting is recognized in this case, which has in good agreement with results of the continuum model obtained by finite element method (FEM) analysis based on a crystal plasticity theory [3]. However, ductile deformation with emissions of screw dislocations from the tip is observed in every case of mode III (anti-plane shear) loading [4], and we found an interesting counterpart in a theoretical result of localized crack tip deformation in a continuum model given by Rice and Nikolic [4, 5]. In this paper, dynamic transition of the atomic structure around the tip in an f.c.c, crystal under mode II (in0921-5093/94/$7.00 ,%'DI 0921-5093(93)02520-D

plane shear) loading is studied by MD simulation. The "N-body" potential proposed by Finnis and Sinclair [6] is used to describe the interatomic interaction. The critical condition of dislocation activation and its temperature dependence are discussed, and the interesting correspondence concerning microdeformation mechanisms in the near-tip field predicted by FEM analysis is found in an atomic scale structural rearrangement involving dislocation activation at the crack tip.

2. Atomic model for molecular dynamics simulation A crack is assumed to be in an f.c.c, single crystal and is subjected to mode II loading as shown in Fig. 1. The atomic model for MD simulation has size 20~f2~a~ x 20/,[2ao in the x-y plane (a0 is the lattice constant) as seen in Fig. 2, which is regarded as two layers out of the periodic structure in the z-direction. The material is supposed to be copper single crystal, the interatomic interaction of which is described by the "N-body" potential proposed by Finnis and Sinclair [6], together with the material properties presented by parameters given by Ackland et al. [7]. Three models are introduced for analysis. In each model the initial crack has its front and its plane as follows: model 1, in [110] on (001); model 2, [ 110]( ] 11 ); model 3, [ 110]( ] 10). The crack is idealized as a vacancy sheet, the details of which are different for each model as shown in Fig. 3. The symbols D and F are used to indicate the type of boundary condition. © 1994 - Elsevier Sequoia. All rights reserved

H. Kitagawa et al.

264

/

Molecular dynamics study of crack processes

The external load is added on the boundary atoms (indicated by solid circles in Fig. 2) by their displacements (denoted by D) or the forces acting on them (denoted by F), which are estimated by an anisotropic elastic field with K. singularity [8] The magnitude of time step used in integration of the equation of motion is 10 -~5 s (fs) and the change in external field during the time step is assumed to be K, = 0.005 MPa m ~/2 p s - I.

~ Q

al Direction

ack Surface Direction of

~

~

i

rack gxtellsi°n

Fig. 1. Crack under mode II loading.

3. Continuum model for finite element method analysis To produce a phenomenological description, crack tip deformation in a plane strain field in an f.c.c, single crystal is analyzed by the FEM. The consitutive properties of the material are supposed to be described by Asaro's crystal plasticity theory [9]. All the slip systems which it is possible to activate, a total of 12, are taken into consideration. A slight strain rate sensitivity is introduced for convenience of numerical treatment. Details of the FEM analysis are the same as discussed in the previous paper [2]. The finite element mesh used in the analysis is shown in Fig. 4. The atomic model comprises a perfect crystal. This means that there is no dislocation multiplication mechanism in the model and it seems to be hard to generate a dislocation in any domain except the crack tip (and the crack surface). In order to introduce such structural completeness of the atomic model, FEM analysis is also performed using the model in which the possibility of crystal slipping is restricted to slip planes which pass the crack tip.

-X-X+X-?X-X.X-X+X.X-X-X-X.X-X'X-X-X-X':" °.-.~.~.%%:.L%¢.:.%%o.%Lo.o.%:.%%c?.~.%o.%=.~.%:.L%:.%L~.*.'.

•

:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: "&%%°.e.¢.~.%%°.%%%°.%%%e.~.%%%=.%%%e.%%%:fl.e.e.=j.*.%'2 "2.e.%=.*.:.=.%%%%~.%~.%~.%%:.~.%%~.L%%~?.~.%Le.%LLLL~.'. ".'p.%~.%:.=.*p.:.%%%=.%%%%L%%%%~.e.%%%~.%LLLLLL:.¢.*.

I(l] prescribed

o o

:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: • .-.o.~.%L%o.o.%e.%%%~.%L~.%:.%%:.~.o.%~.~.%%%¢?.~.~.~j.L%-.

•

================================================================================= =o=,oo~,oo:==:o:a°..

Nodes 1654 Elements 1592 (×4)

°°oo.oo*eo*o*~oo~.oo.~o.o.o.oooo~.o.o.o*e~.e°=¢o*~°°

iii!iiiiiigigiiii?iiiiii:iiiii#i!ii?:ili!ililiiiiiiii!iiiiiiiigliiiiiiii °o°°.,.°°°.°.,°°.°.°°°°.-.°.°.°°°.°.°.-°°°.°°°.°.°.,°°.-,°...°.-o,

Fig. 2. Atomic model. lOOx

• ,

o

•

,

a

~

o

,

o

,

•

,

,

o

o

,

o

........

o • . •

. o

.

.

. . . . o o , . . . .

I

•

.

o

,

o

•

o

•

. o

.

o

.

o :.

. o . o . o . Q .

o

,

.

.

o . o

. o . , ,.

(iii) D2

•

c

•

c

•

c

•

/

c

o . o . c

J /

~ . o

.

444 x

.* *. . o . * -

. * . o -

o - * . * . * . ~-

o.

J

.

. o . o . o.

- * - * .

o . ~ . o - * .

o

o . o . o . o . •

o . o . o . o . o , °

o

. o . o . a . a

a . * . * . o . o ,

:

.

o

(ii) D3, F3

o . , . o . Ù . o ,

n . o . *

,

,

o . o . o . ¢ o

: . o . o

~ . , . o . o . o .

. * . o . a . o

. o

(i) DI, F1 3

o

.

o

.

o

o . o . Q . c

o

~ ,

~

, . , . o )

o . o . ~

~.

\

°

,I

(iv) F2

Fig. 3. Details of atomic configuration in near-tip region.

Fig. 4. Finite element mesh.

H. Kimgawa et al.

/

,~lolecular dynami('s study qf crack processes

4. Results and discussion

4.1. A t o m i c structure around the crack tip Typical atomic structures appearing around the crack tip are shown in Fig. 5. The initial temperature is set at 40 K. To show the rearranged atomic feature more easily, lines are drawn through the atoms of the upper layer in the I 111 ] direction. In every case, partial dislocation is generated at the tip and moves away quickly on the / 111 } plane without any brittle fracture process. The moving-away directions are 0= + 125.6 ° (model 1), 0 ° (model 2) and 35.26°/144.74 ° (model 3) from the x-axis, which agree with the dense packing directions of the f.c.c, crystal. The temperature or the kinetic energy of the system rises a little just after dislocation nucleation in the displacement specified case (D). However, in the force Fixed Displa(em{~nt. B.C. (D)

"':"

:':

"c:c.' :: ....

"

. _..:, - :

:

specified case (F), vibration of the whole system appears and its amplitude increases with the external load, and the temperature increases jumpwise after dislocation nucleation. Thus. the fi)rce specified boundary condition gives dynamically too complex results to describe the essence of the fracture process at the tip. In the following discussion, we concentrate on the results obtained by the displacement specified boundary condition. 4.2. Estimation o['the critical value of K// The elastic stress distribution around the crack tip including the K n singularity is expressed by

r,,,(r.

A£11

F(0)

Fixed Fo~{e B.C. (F)

,:± • • :.:::: , " -

0 551XIPa~ n]

2MP:u

n

(S) .

i'i"i~!'i:i':i:£~i:i~:i:?'L~i:i'i'':'~:i:i'i:£'!:~7!!~';i£ !~:? ~

,'."P",'.~," , ' A % ' , ; , ? ; ' P U , P,]~V,:~ 2 P

;.i?-:i7 ;:i~'.i:7!:17!?.2'.i?:!ii,i'.177:.;!t.i: ....

•

. . . . . . . . . . . . . . . . .

..

[~

.

rl - I} 5 5 M P a ~

• '.',..:

'}i"1"' ,:;,,,.,,

<', m

.

'.,'P,; ,'.'.'V',,

l~ II

,',,,

0 '2qMPa~ m

(c)

Ixn = o 55MI)a,,m

265

[~ tl = 0 2 9 M P a ~ Ti]

Fig. 5. Atomic structure appearing in the deformed state: (a) model 1, (b)model 2, (c) model 3.

(])

H. Kitagawa et al. / Moleculardynamics study of crack processes

266

If the dislocation is generated at the tip of some plane at an angle 0 from the crack plane, under the condition that the shear stress at the point at distance r 0 from the tip reaches a certain critical value T0, we have a critical condition,

K~; -

l:0(2nr()) m F(O)

(2)

Values of KII cr and ~'0(2,7r0) 1/2 obtained from MD simulation are given in Table 1. The latter values which should be basically constant are distributed between 0.15 and 0.39. Variation of the details of atomic arrangement at the tip for each model seems to be the main reason for this dispersion. According to Rice [10], when the dislocation is emitted in the forefront direction of the crack, i.e. model 2, the critical condition of dislocation nucleation is estimated by (2~Fus ~/2 KllCr = k 1 - v

(3)

in which Fus is the unstable stacking energy. If the FS potential used in MD simulation is applied, ~, v and Fu~ are evaluated as ~ = 4 1 GPa, v=0.32, and Fo~ = 0.37 J m-2, and introduction of these values into eqn. (3) gives KHcr= 0.21 MPa m U2. This result agrees well with KIIcr = 0.20 MPa m U2 found in Table 1.

4.3. Comparison with finite element method analysis Figure 6 shows the deformed mesh and distribution of largely slipped crystal systems for models 1 and 3. The short line segments indicate the slip direction and its magnitude, the latter by the thickness of the line. If it is assumed that the full slip systems can be activated in the whole region, a so-called kink-mode slipping band [ 11 ] appears, the direction of which is perpendicular to the slip direction developed. This deformation does not appear in nature in the atomic model because of completeness of the lattice. If the potential slip system is restricted to a system whose slip plane passes the

crack tip, the kink-mode vanishes and the direction of the localized shear band agrees exactly with the dislocation moving-away direction.

4.4. Stress distribution Stress on the atomic scale can be evaluated through statistical average of the force on each atom from the surrounding. The stress component with respect to an atom i is evaluated by 1

i

Oafl = ~

ij

neighbors

E

ij

Veff,(rij) rar(~

j# i

(4)

rq

where if2 is atomic volume, Veff is the effective potential [6], rO is the distance between the i and j atoms and ra ij is the component of the vector pointing from the j to the i atom. Ordinarily, because of thermal fluctuation, the above equation gives a value which is unsteady with time and varies diversely with spatial location, but a simple mean with respect to some (short) time interval behaves smoothly. Moreover, before dislocation nucleation it agrees well with the linear elastic solution for the KII singularity field as a whole. The plots of Fig. 7 show the results of model 2. However, a nucleated dislocation experiences a repulsive force from the KII singularity stress field, and then moves away quickly and piles up at a certain location determined by the external mechanical constraint. Thus, as seen clearly in Fig. 7, a dislocation-free zone is produced in front of the crack tip. These plots are also explained well by the linear elastic solution of the interaction field produced by singularities of the crack tip and edge-dislocation with Burgers' vector b[12],

x-xl

rx>, (2zrx)l/2 t - 2 : ~ ( l _ v )

An example of the good correspondence between the above theoretical prediction and the results of MD simulation is shown in Fig. 8.

TABLE 1. Critical condition of dislocation nucleation Model

1 2 3

0 (deg)

- 125.26 125.26 0.00 35.26 144.74

1/F( O)

- 1.240 - 1.240 1.000 1.765 - 1.575

D (MPa m 1/2)

KItc'

r<)(2 gr0) I/2 (MPa m I/2)

K. ~' (MPa m 1/2)

To(2~ro) 1/2 (MPa m 1/2)

0.42 0.48 0.20 0.53 0.43

0.34 0.39 0.20 0.30 0.27

0.19 -0.28 0.28 --

0.15 -0.28 0.16 --

H. Kitagawa et al.

/

267

Molecular dynamics study of crack processes

/ ,,

,~4',,

7

.\

\\~

I D e f o r m e d M(,,~h( x 10)

She;Lr Slip

Fig. 6. Deformed mesh and distribution of the active slip system (finite element analysis): (a) model 1, (b) model 3.

I

I

I

I

I

k

I

a. ©

¢ I = 0.i75 ~ 0.200 M P a v ~

GQ

I

i

[

~e

¢

,0.

..a ct~

~_

i

. . . .

i

o o o o

K. = 0.4%~ o.45oMm,/~

~,

~9 GO

¢.O

0.075

~

O.lOO

0.025 ~ 0.050

~?"

0.275 ~ 0.300 --

-4

I

4

Linear Elastic Solution I

I

I

8

0.225 ~ 0.250 I

-4

12

x Into]

I

I

4

~o I

I

I

8

I

12

z [nm]

Fig. 7. Stress distribution in the forefront region of the crack tip (model 2): (a'l before dislocation nucleation, (b) after dislocation emission.

4.5. Discussion of the temperature dependence of the dislocation nucleation condition Roughly speaking, the dislocation nucleates at the instant when the resolved shear stress at the crack tip excited by the external force overcomes the potential barrier produced by the interatomic interaction as seen

in Section 4.2. Thermal fluctuation may promote the process, the tendency of which should be described by Arrhenius' equation. If the inherent activation energy to be overcome is independent of temperature and the exciting energy due to the external agency is proportional to the square

268

H. Kitagawa et aL

/

Molecular dynamics study of crack processes

of K., the probability of dislocation emission can be decided by a condition expressed by

observed, which might be the effect of characteristics caused by dynamic properties of the system.

(6)

UI,- A (K,,~) 2 = B T

where A and B are constants, U 0 is the inherent activation energy and T is the absolute temperature. T h e temperature dependence of Ku cr with the initial temperature obtained by M D simulation is shown in Fig. 9. T h e model used and method of calculation are same as model 2 (D) stated above, but the number of atoms included in the system is somewhat smaller. A clear temperature dependence of KIIcr is seen in Fig. 9, which is primarily explained as a thermally activated process. T h e dotted curve shows the least-square fitting based on eqn. (6). However, it should be noted that a stepwise relation between KII cr and T is

....

/f

° , ,

• ° " '° °.

%o.:.. o¢* ".

5. Concluding remarks Rearrangement of the atomic scale structure around the crack tip in an f.c.c, crystal under in-plane shear (mode II) loading has been discussed on the basis of MD simulation and its counterpart in the continuum crystal model obtained by F E M analysis. T h e main results obtained concerning the atomistic fracture process are summarized as follows. A ductile process involving dislocation nucleation from the crack tip is always observed. T h e incipient dislocation moves away quickly from the tip and a dislocation-free zone remains in the near-tip region. T i m e averages of local stress in the near-tip region coincide well with linear elastic prediction. T h e critical stress intensity factor estimated by MD simulation agrees well with Rice's theoretical prediction. T h e temperature dependence of the dislocation nucleation condition can be explained primarily as a thermally activated process.

References

O"J ,zZ 0

_,

, [ 0

5

15 m [nm] ( KII = 0.325 ~ 0.350MPav/-~ ) i0

Fig. 8. Stress distribution in the forefront region (model 2). comparison of theoretical prediction and MD simulation.

~=0.25m '

,

1 B. deCelis, A. S. Argon and S. Yip, Molecular dynamics simulation of crack tip processes in alpha-iron and copper, J. Appl. Phys., 54 (9)(1983) 4864-4878. M. Mullins, Atomic simulation of cracks under mixed mode loading, Int. J. Fract., 24 (1984) 189-196. S. Kohlhoff, S. Schmauder and P. Grumbsch, Coupled atomistic-continuum calculations of near interface cracking in metal/ceramic composites, metal-ceramic interface, ActaScr. MetalL Proc. Set. 4, Pergamon, 1990, pp. 63-70. 2 H. Kitagawa and A. Nakatani, Crystal plasticity model and molecular dynamics (numerical simulation of plane strain inhomogeneous deformation). In R. Wang and D. C. Drucker (eds.), Constitutive Relation for Finite Deformation of Polycrystalline Metals, Proc. IUTAM-Symp., Beijing, 1991,

Springer, Berlin, 1992, pp. 220-233. 3 H. Kitagawa and A. Nakatani, Microstructural aspects of crack extension in a crystalline material (a molecular dynamics study). In M. Jono and T. Inoue (eds.), Mechanical

,

Behavior of Materials-V1, Proc, 1CM-V1, Kyoto, 1991,

4 ._~ 0.2

S "~-. c~% -,.o?o~a

E

O..

5

•~ 0 . 1 5

6

0

o. ~6

. . . .

16o Temperature

. . . .

2oo

T IK]

Fig. 9. Temperature dependence of the critical value of K. (model 2).

7

Pergamon, Oxford, 1991, pp. 111-116. H. Kitagawa and A. Nakatani, Study of computational modelling for materials with crystalline structure, II molecular dynamics simulation of microscopic crack tip field under anti-plane loading. Proc. Jpn. Soc. Mech. Eng., Ser. A, 59 (557)(1993) 256-262. J. R. Rice and R. Nikolic, Anti-plane shear cracks in ideally plastic crystals, J. Mech. Phys. Solids, 33 (6)(1985) 595-622. M. W. Finnis and J. E. Sinclair, A simple empirical N-body potential for transition metals, Philos. Mag. A, 50 (1) (1984) 45-55 (erratum, 53 (1986) 161). G. J. Ackland, G. Tichy, V. Vitek and M. W. Finnis, Simple N-body potentials for the noble metals and nickels, Philos. Mag. A, 56 (6) (1987) 735-756.

tt. Kimgawu et al.

/

Molecuhlr dynamics study ~( cruck processes

8 G. C. Sih, R C. Paris and G. R. Irwin, On cracks on rectilinearly anisotropic bodies, Int. J. Fract. Mech., 1 (1965'3 189. 9 R. J. Asaro, Micromechanics of Co'stals and Polycr)'smls, Advances in Applied Mechanics, Vol. 23, Academicr Press, New York, 1983, p. 1. 10 J. R. Rice, Dislocation nucleation from a crack tip: an analy-

269

sis based on the peierls concept, J. Mech. /'los. Solids, 40 (2) i1992) 239-271. 11 R. Moham M. Ortiz and C. F. Shih, An analysis of cracks in ductile single crystals l, anti-plane shear. J. Mech. Phys. SoIMs, 40 (2',(1992) 1951-1976. 12 S. M. Ohr, Dislocation-crack interaction, J. Phys. ('hem. Solids, 48 ( 1 I ) (1987) 1007-1014.