The dependence of the fractal dimension of crack on material for brittle fracture in two dimensions

Physics Letters A 308 (2003) 455–460 www.elsevier.com/locate/pla

The dependence of the fractal dimension of crack on material for brittle fracture in two dimensions S.G. Wang International Centre for Materials Physics, Institute of Metal Research, Chinese Academy of Sciences, Shenyang 110016, PR China Received 8 August 2002; received in revised form 6 December 2002; accepted 9 December 2002 Communicated by V.M. Agranovich

Abstract Two-dimension spring square network model involving the nearest and second neighbour bonds is studied by numerical simulation with Gauss–Seidel iteration for brittle fracture Mode I. The material-dependent parameter, G0 /G, is introduced in this simulation, G and G0 are the parameters related to the elastic moduli of the nearest and second neighbour bonds, respectively. Brittle fracture of materials with different G0 /G are simulated under the same quenched disorder. The fractal dimension of crack fluctuates with G0 /G. It seems that it is impossible for different materials to have a universal fractal dimension of crack in brittle fracture. However, it is possible to have different so-called universal fractal dimensions.  2003 Elsevier Science B.V. All rights reserved. PACS: 05.45.Df; 62.20.Mk; 81.40.Np Keywords: Fracture; Crack; Fractal dimension

1. Introduction In recent years, there have been extensive research and great development in quantitative fractography [1]. In experiments, some investigators have obtained the fractal dimensions of crack for different materials in two dimensions, such as 1.33 and 1.31 [2] for Cer.ging and Un.inc brittle fracture on micron scale, respectively; 1.1–1.2 [3] for intermetallic V3 Au containing different oxygen content brittle fracture on micron scale; from 1.05 to 1.33 [4] for six alumina materials and five glass-ceramics brittle fracture; from 1.179 to 1.447 [5] for 0.02 wt% carbon pure iron

E-mail address: [email protected] (S.G. Wang).

(at −196 ◦ C) and 11 wt% Cr-2 wt% Ni steel brittle fracture (at −80 ◦ C) on micron scale and from 1.23 to 1.86 [6] for single crystal silicon brittle fracture on nanometer scale etc. Although the fracture surfaces of different materials may have the same fractal dimension, many investigators inclined to note that fractal dimension depends on materials, microstructures and fracture modes etc. [7,8]. By computer simulations, in two dimensions the fractal dimensions are 1.58 ± 0.12 [9], 1.65 ± 0.05 [10], 1.7 [11] and 1.2 [12], the authors of Refs. [9–12] consider fractal dimension as a universal value for fracture surfaces. However, is it possible that there have different universal fractal dimensions for brittle fracture surfaces? and why can we obtain the same fractal dimension in experiments for different materials? In this Letter, the dependence of fractal di-

0375-9601/03/$ – see front matter  2003 Elsevier Science B.V. All rights reserved. doi:10.1016/S0375-9601(02)01774-7

456

S.G. Wang / Physics Letters A 308 (2003) 455–460

tively,

Fig. 1. The two-dimension elastic model for brittle fracture, the node in centre is node (i, j ).

Kn (i, j ) = G + W1 f1 (i, j ),

(1)

Gn (i, j ) = G0 + W2 f2 (i, j ).

(2)

G, W1 and G0 , W2 , are the parameters related to the elastic moduli of the nearest and second neighbour bonds, respectively. f1 (i, j ) and f2 (i, j ) are the random numbers distributed uniformly on [0 ∼ 1.0), W1 f1 (i, j ) and W2 f2 (i, j ) denote the quenched disorders of the elastic moduli of the nearest and second neighbour bonds, respectively. Twenty kinds of disorder configurations are taken into account in this simulation. ux (i, j ) and uy (i, j ) denote the coordinates of node (i, j ) in x-direction (i) and y-direction (j ), respectively. The force equilibrium equation of each node (i, j ) in Fig. 1 is m=4


  am Rm (i, j ) − 1 Km (i, j )em

m=1

+

n=4


√   cn rn (i, j ) − 2 Gn (i, j )sn = 0,

n=1 m=4


Fig. 2. The fracture Mode I.

  bm Rm (i, j ) − 1 Km (i, j )fm

m=1

mension of crack on materials is investigated by computer simulation for brittle fracture surface. The conclusion of this Letter can contribute to the two problems. The material-dependent parameter, G0 /G, is introduced in this simulation for brittle fracture (G and G0 are the parameters related to the elastic moduli of the nearest and second neighbour bonds, respectively). In this Letter, two-dimension square spring network [10,12] (L×L) (Fig. 1) model is employed to simulate crack propagation for brittle fracture Mode I (Fig. 2) for different G0 /G with the same quenched disorder. G0 /G of twelve bcc and fcc metals and eight hcp metals are calculated, it is a material-dependent parameter. Two morphologies of brittle fracture surfaces are shown.

+

n=4


√   dn rn (i, j ) − 2 Gn (i, j )tn = 0,

(3)

n=1

where am , bm and cn , dn are the director factors of the nearest and second neighbour bonds in x-direction (i) and y-direction (j ), respectively, they are 1 or −1. em , fm and sn , tn are the component factors for the nearest and second neighbour bonds in x-direction and y-direction for node (i, j ), respectively. Rm (i, j ) and rn (i, j ) are the nearest neighbour and second neighbour distances for the node (i, j ), respectively. They are Rm (i, j )m=1–4  2  2  = ux (i, j ) − umx (i, j ) + uy (i, j ) − umy (i, j ) ,

2. Simulation details

rn (i, j )n=1–4  2  2  = ux (i, j ) − unx (i, j ) + uy (i, j ) − uny (i, j ) ,

In Fig. 1, the elastic moduli of the nearest and second neighbour bonds for each node (i, j ) are Kn (i, j )(n=1,2,3,4) and Gn (i, j )(n=1,2,3,4), respec-

(4) umx (i, j )m=1–4 and umy (i, j )m=1–4 are the coordinates of four the first neighbour atoms for node (i, j ).

S.G. Wang / Physics Letters A 308 (2003) 455–460

unx (i, j )n=1–4 and uny (i, j )n=1–4 are the coordinates of four the second neighbour atoms for node (i, j ). At first, the strain σ is applied in the specimen showed in Fig. 2, form the first row, j = 1, the tip of crack, to 20th row j ∈ [1, 20], different strains, which are smaller than σ were applied, it decreases with j increasing. Simultaneously, ux (i, j ) of nodes from 1st to 20th row j ∈ [1, 20], if i ∈ [1, L2 − 2], increase by certain value, the larger i is, the smaller amplitude of increment is; if i ∈ [ L2 + 1, L], ux (i, j ) of these nodes decrease certain value, the larger i is, the larger amplitude of decrement is; if i is L2 − 1 or L2 , ux (i, j ) of these nodes increase the schematic value, it is negative or positive. uy (i, j ) of each node from 1st to 20th row j ∈ [1, 20] increases certain value, the amplitudes of increment decreases with j increasing. For each node (i, j ) from 1st to 20th row j ∈ [1, 20], Eq. (3) is solved, and then the nearest neighbour distance of each node in 1st row was calculated. If the distance of the nearest neighbour exceeds the critical value lc , the bond between two nearest neighbour nodes is broken irreversibly, and then the two nodes will displace certain distance along the bond in reverse directions. In fact, the strain field around crack tip is long ranged, Eq. (3) should be solved for the whole system with the boundary condition for the stress perpendicular the crack surface equal to zero. However, the farther distance away from form the crack tip, the less effect of node on crack propagation is, and it will be time consuming if Eq. (3) is solved for the whole system. In this Letter, let it be supposed that the atoms after (j + 21)th rows away form crack tip have little effect on crack propagation if the crack tip lie in j th row. When the crack is in j th row, from 1st to (j + 20)th row, different strains are applied in y-direction, uy (i, j ) of each node will increase like above; ux (i, j ) of each node will increase or decrease for each node like above also. The amplitude of increment of strain in y-direction (j ) is about the two order of magnitude larger than the amplitude of increment or decrement of strain in x-direction (i) in this Letter. For the same suppose mentioned above, Eq. (3) is solved for each node (i, j ) from j th to (j + 20)th row, and then the nearest neighbour distance of each node of j th row was calculated. If the distance exceeds the critical value lc , the bond between them is broken irreversibly, and then the two nodes will displace certain distance

457

along the bond in reverse directions. However, in solving Eq. (3), we only consider the effect of the nearest neighbour and second neighbour atoms on Eq. (3), and neglect the effect of other neighbour atoms on Eq. (3). It would be very troublesome to solve Eq. (3), if the effect of other neighbour atoms on Eq. (3) are considered. Because the direct interaction of other neighbour atoms is very weak, the effect of other neighbour atoms on Eq. (3) can be neglected in general, although error can be caused in solving Eq. (3). Eq. (3) is solved with Gauss–Seidel iteration method, which has been applied in the simulation of percolation process in hydraulic engineering and can applies to the simulation of fracture process [13]. The stopping condition is that residual error is less than 10−6 , it is sufficient for this simulation. Fractal features of these cracks are examined by the method in Ref. [9]. Nc is the total number of the bonds of network, it is proportional to the area of the network in the continuum limit; Nf is the total number of broken bonds, it is proportional to the length of perimeter of fractal crack. The following relation approximately holds [9]: Nf ∝ Nct .

(5)

The fractal dimension D of crack may be elevated as follows [9] (area)1/2 ∝ (length)1/D

(6)

and then D = 2t,

(7)

t can be obtained by the logarithmic relationship between Nf and Nc form Eq. (5), the logarithmic relationship between Nf and Nc for G0 /G = 0.3914 shows in Fig. 3, the slope of line is equal to t. Form Fig. 3, one can use Eqs. (5)–(7) to define the fractal dimension of crack. Two morphologies of brittle fracture surface (800 × 800) are shown in Figs. 4 and 5. The values of G0 /G are 0.00532 and 0.542 in Figs. 4 and 5, respectively; fractal dimensions of the cracks are 1.766 and 1.538 in Figs. 4 and 5, respectively, and these cracks of brittle fracture are river line with fractal structure. The curve of fractal dimension of crack versus G0 /G is shown in Fig. 6, which contains seventy five different G0 /G. In Fig. 6,

458

S.G. Wang / Physics Letters A 308 (2003) 455–460

Fig. 3. The logarithmic relationship between Nf and Nc for G0 /G = 0.3914.

Fig. 5. The morphology of brittle fracture surface (G0 /G = 0.542, D = 1.538), the arrow denotes the propagation direction of crack.

Fig. 4. The morphology of brittle fracture surface (G0 /G = 0.00532, D = 1.766), the arrow denotes the propagation direction of crack.

for each fractal dimension of G0 /G twenty kinds of disorder configurations are taken into account, so the data of fractal dimension in Fig. 6 are reliable. Fractal dimensions range from 1.536 to 1.925 when G0 /G varied systemically, we obtain them under the same disorder. Fig. 6 means that fractal dimension varies with G0 /G, however, there may have the different socalled universal fractal dimension for different G0 /G [14]. In Fig. 6, L is 800. From Fig. 6, the fluctuation

Fig. 6. The relationship between fractal dimension D and G0 /G.

of fractal dimension results from different G0 /G. The reason why we can obtain Fig. 6 may be that brittle fracture is random process and the sensitivity of different G0 /G to the same quenched disorder is different.

3. Calculation of G0 /G To our knowledge, G0 /G is first introduced in this simulation for brittle fracture. It is calculated with

S.G. Wang / Physics Letters A 308 (2003) 455–460

Morse potential function [15]   ϕ(i, j ) = A e−2α(rij −r0 ) − 2e−α(rij −r0 ) ,

(8)

where α and A are the constants with dimensions of reciprocal distance and energy, respectively, and r0 is the equilibrium distance of approach of the two atoms, rij is the distance between atom i and atom j   2  ∂ ϕ(i, j )   , G =  (9) ∂rij2 rij =r01   2  ∂ ϕ(i, j )    , G0 =  (10) ∂rij2 rij =r02 Table 1 The values of α, r0 , r02 , G0 /G for six bcc and six fcc metals Agd Nid Cud Ald Pbd Cad Bac Nac Moc Wc Crc α-Fee a b c d e

α a (Å−1 )

r0 a (Å)

r02 b (Å)

G0 /Gc

1.369 1.4199 1.3588 1.1646 1.1836 0.8054 0.6570 0.5899 1.5079 1.4116 1.5721 1.3885

3.115 2.780 2.866 3.253 3.733 4.569 5.373 5.336 2.976 3.032 2.754 2.845

4.09 3.52 3.61 4.05 4.95 5.58 5.02 4.225 3.15 3.16 2.88 2.87

0.0536 0.0355 0.0362 0.0245 0.0579 0.0134 0.333 0.466 0.133 0.181 0.174 0.238

From Ref. [15]. From Ref. [16]. Present calculation. fcc metal. bcc metal.

459

  1 G0  2e−2α(r02− 2 r0 ) − e−αr02  = 1 , G 2e−2α(r01− 2 r0 ) − e−αr01

(11)

where r01 and r02 are the distances of between two nearest neighbour atoms and between two second neighbour atoms, respectively. We can obtain G0 /G by Eq. (11), the value of α, r0 , r02 and G0 /G of six bcc and six fcc metals are tabulated in Table 1. For hcp metals, α is calculated with the following equation [17] α=

3Cv , 4γ E0 r0

(12)

where γ and Cv denote the thermal expansion coefficient and specific heat at constant volume, respectively. E0 denotes the latent heat of vaporization. The value of γ , Cv , r0 , r02 and G0 /G of eight hcp metals are tabulated in Table 2. In view of Eq. (11), Table 1 and Table 2, G0 /G is a material-dependent parameter, it is more than a few the ratio of bond energy between the nearest and second neighbour bonds [20]. It is difficult to calculate the value of G0 /G for alloy, however, it may be evaluated roughly with the ratio of bond energy between the nearest and second neighbour bonds. Only when G0 /G → 0 and elastic modulus is a constant, Eq. (3) can become Laplace equation [10,21] 'u = 0.

(13)

In Fig. 6 fractal dimension varies with G0 /G. However, it is possible to have different so-called universal constants [22].

Table 2 The value of α, E0 , Cv , r01 , r02 , γ , G0 /G for eight hcp metals Y Ti Re Mg Zn Er Ho Lu a b c d

α d (Å−1 )

E0 a (MJ/Kg)

Cv a (MJ/Kg K)

r01 a (Å)

r02 a (Å)

γ b ×106

G0 /Gd

2.202 0.792 2.504 4.538 6.046 1.913 1.939 1.370

4.777 9.83 3.417 5.30 1.782 1.896 1.824 2.444

0.499 0.300 0.257 1.025 0.382 0.1678 0.1646 0.1531

3.557 2.890c 2.746 3.196 2.659 3.470c 3.49c 3.43c

3.648 2.950c 2.760 3.209 2.665 3.56c 3.58c 3.50c

12.0 10.09 5.40 25.7 31.7 12.3 10.7 15.3

0.521 0.865 0.916 0.835 0.896 0.576 0.651 0.742

From Ref. [18]. From Ref. [19]. From Ref. [16]. Present calculation.

460

S.G. Wang / Physics Letters A 308 (2003) 455–460

4. Discussion and conclusion

References

It is not surprised that some researchers obtained different universal values in experiments and simulations, however, others concluded that fractal dimension depends on materials, microstructures and fracture modes etc. in experiments. The reason may be that the material parameter, G0 /G, did not varied systemically like Fig. 6 in experiments. One may obtain a universal fractal dimension in several specimens in experiments [23]. In computer simulations, the effect of the second neighbour bond on Eq. (3) is neglected [24] and Eq. (13) was used for brittle fracture [10,21]. However, this may draw the conclusion of universality [9,21]. No material-dependent parameters were introduced in simulation of brittle fracture for disorder media [9–11,22], only the very general parameters, such as the dimension of space rather than the microstructure were introduced, which may imply the conclusion of universality [8]. There may have different socalled universal fractal dimensions in some cases from Fig. 6. Form Tables 1, 2, Fig. 6 and Eq. (11), the fractal dimension of crack for brittle fracture is related to G0 /G, and G0 /G is a material-dependent parameter, it seems that fractal dimension cannot be a universal constant independent of materials. It seems that it is impossible to have a universal fractal dimension of fracture surface for different materials with the same quenched disorder. However, it is possible to have different so-called universal constants in some cases.

[1] B.B. Mandelbrot, D.E. Passoja, A.J. Paullay, Nature 308 (1984) 721. [2] G.R. Baran, C.R. Carmes, D. Wehbl, J. Am. Ceram. Soc. 75 (1992) 2687. [3] Y. Fahmy, J.C. Russ, C.C. Koch, J. Mater. Res. 6 (1991) 1856. [4] J.J. Mecholsky, D.E. Passoja, K.S. Feinberg-Ringel, J. Am. Ceram. Soc. 72 (1989) 60. [5] X.W. Li, J.F. Tian, H.H. Su, Z.G. Wang, J. Mater. Sci. Lett. 15 (1996) 2137. [6] M.W. Mitchell, D.A. Bonnell, J. Mater. Res. 5 (1990) 2244. [7] C.W. Lung, S.G. Wang, Q.Y. Long, Physica B 279 (2000) 139. [8] V.Y. Milman, N.A. Stelmashenko, R. Blumenfeld, Prog. Mater. Sci. 38 (1994) 425. [9] H. Takayasu, Phys. Rev. Lett. 54 (1985) 1099. [10] H. Takayasu, Prog. Theor. Phys. 74 (1985) 1343. [11] H.J. Herrmann, A. Hansen, S. Roux, Phys. Rev. B 39 (1989) 1678. [12] Y. Termonia, P. Meakin, Nature 320 (1986) 429. [13] H.F. Wang, M.P. Anderson, Introduction to Groundwater Modeling, Freeman, San Francisco, 1982. [14] B.B. Mandelbrot, The Fractal Geometry of Nature, Freeman, San Francisco, 1982. [15] L.A. Girifalco, V.G. Weizer, Phys. Rev. 114 (1959) 687. [16] C. Kittel, Introduction to Solid State Physics, 7th Edition, Wiley, 1996. [17] S.S. Mitra, S.K. Joshi, J. Chem. Phys. 34 (1961) 1462. [18] Metals Handbook, Vol. 2, 9th Edition, American Society of Metals, 1988. [19] D.B. Sirdeshmukh, J. Appl. Phys. 37 (1966) 3910. [20] S.G. Wang, E.K. Tian, C.W. Lung, J. Chem. Phys. Solids 61 (2000) 1295. [21] H.J. Herrmann, S. Roux, Statistical Models for the Fracture of Disordered Media, North-Holland, Amsterdam, 1990. [22] A. Hansen, E.L. Hinrichsen, K.J. M˘aløy, S. Roux, Phys. Rev. Lett. 71 (1993) 205. [23] K.J. M˘aløy, A. Hansen, E.L. Hinrichsen, S. Roux, Phys. Rev. Lett. 68 (1992) 213. [24] M. Sahimi, S. Arbabi, Phys. Rev. B 47 (1993) 713.

Acknowledgements The author is grateful for the support of the National Natural Science Foundation of China under Grant No. 19874064.