A Simplified Probabilistic Model for Nanocrack Propagation and Its Implications for Tail Distribution of Structural Strength

This paper presents a simplified probabilistic model for thermally activated nanocrack propagation. In the continuum limit, the probabilistic motion of the nanocrack tip is mathematically described by the Fokker-Planck equation. In the model, the drift velocity is explicitly related to the energy release rate at the crack tip through the transition rate theory. The model is applied to analyze the propagation of an edge crack in a nanoscale element. The element is considered to reach failure when the nanocrack propagates to a critical length. The solution of the Fokker-Planck equation indicates that both the strength and lifetime distributions of the nanoscale element exhibit a power-law tail behavior but with different exponents. Meanwhile, the model also yields a mean stress-life curve of the nanoscale element. When the applied stress is sufficiently large, the mean stress-life curve resembles the nasquin law for fatigue failure. nased on a recently developed finite weakest-link model as well as level excursion analysis of the failure statistics of quasi-brittle structures, it is argued that the simulated power-law tail of strength distribution of the nanoscale element has important implications for the tail behavior of the strength distribution of macroscopic structures. It provides a physical justification for the two-parameter Weibull distribution for strength statistics of large-scale quasi-brittle structures.


INTRODUCTION
One of the primary interests in reliability analysis of engineering structures is to determine the structural load capacity at low failure probabilities. Obviously neither brute force numerical simulations nor direct experiments could reveal the tail portion of the failure statistics. Therefore, analytical modeling becomes an essential tool for studying the strength distribution including its tail behavior.
In recent years, considerable efforts have been directed towards developing a physics-based framework for modeling the failure statistics of engineering structures that are made of brittle heterogeneous (a.k.a. quasi-brittle) materials, such as concrete, ceramics, composites, rock, etc. [17]. Previous studies suggested that a logical way to derive the functional form of the strength distribution of macroscopic quasi-brittle struc-tures is through the failure statistics of nanoscale structural elements [2,3,7]. The motivations for such a multiscale approach are twofold: (i) the macroscopic damage and failure of materials originate from crack growth at the nanoscale, and (ii) under most loading scenarios the motion of the nanocrack tip can be considered as quasi-stationary and consequently the corresponding probability of individual nanocrack jump can be calculated directly from the jump frequency.
Over the past two decades, significant advances have been achieved in atomistic modeling of nanoscale fracture [813]. Earlier attempts of atomistic modeling of material fracture relied largely on molecular dynamics simulations using some empirical interatomic potentials. In later studies, efforts have been directed to developing hybrid quantum mechanics/molecular dynamics models with an aim of improving both the accuracy and efficiency of the simulation [10,11,14 16]. Nevertheless, as mentioned earlier, direct atomistic computation would not be able to probe the tail behavior of the failure statistics of nanoscale structures even for typical loading durations of laboratory experiments.
Meanwhile, it has been noted that the nanocrack jump can be regarded as a thermally activated transition between two adjacent metastable states [2,3,7,17]. This transition has been studied theoretically for many decades. One of the most well-known results is the transition rate theory [18,19], which describes the frequency of transition between two metastable states for the case of a large energy barrier. The quasi-stationarity of the nanocrack jump implies that the jump frequency is proportional to the probability of the jump [3,7]. To calculate the overall failure probability of a nanoscale structure, such as a system of disordered nanoparticles or an atomic lattice, we would also need a mechanistic model. Recent studies [3,7] showed that the motion of the nanocrack tip can be modeled as a stress-driven random walk process described by the .okkerPlanck equation. However, the previous qualitative analysis focused only on the mean behavior of the system. Meanwhile, the analysis assumed a constant drift velocity, which is an oversimplification. .or a nanocrack, the drift velocity is directly related to the rate of the crack jump, which generally varies with the crack length [2,3,7].
In this study, we analyze the probabilistic behavior of nanocrack propagation by using the .okkerPlanck equation, in which the drift velocity is explicitly related to the crack geometry. The model is then used to calculate the failure probability of a nanoscale element as a function of both the applied stress and loading duration. The result of this analysis has direct conse-quences for reliability analysis of macroscopic quasibrittle structures.

MODEL .ORMULATION
Consider a crack in some nanoscale element, such as a system of disordered nanoparticles (.ig. 1a). We analyze the motion of the nanocrack tip in a one-dimensional setting (.ig. 1b): The movement of crack tip in the positive x-direction (i.e. forward jump) represents crack propagation while the movement in the negative x-direction (i.e. backward jump) represents crack healing.
The nanocrack experiences discrete jumps and the distance of each jump is equal to one nanoparticle spacing . α δ It has been suggested that the nanocrack jump can be treated as a random process, and these jumps are independent of each other [17]. Therefore, the probability, p(x, t + ∆t), that the crack tip occupies location x at time t + ∆t, which is often referred to as the occupational probability, can be calculated as ( , x − δ is the probability of forward jump from a x − δ to x, ( ) a n x + δ is the probability of backward jump from a x + δ to x, and r(x) is the probability that crack tip remains at location x. It is noted that all the transitional probabilities, m(x), n(x) and r(x), are considered to be time independent.
By expressing all the terms in Eq. (1) by using the Taylor expansion for their corresponding quantities at (x, t), we have .ig. 1. Propagation of a nanocrack in a system of disordered nanoparticles (a) and mathematical representation of the motion of the nanocrack tip (b) (color online).

Crack healing
Crack propagation , x p x t x m x n x t φ = − ∆ is the difference in frequencies of crack propagation and healing at the current tip location, and 2 ( ) [ ( ) ( )] x m x n x t φ = + ∆ is the frequency of leaving the current tip location.
Equation (3) represents the governing equation for the time evolution of the occupational probability. To incorporate the physics of nanocrack propagation into Eq. (3), we need to further relate the transitional probabilities (i.e. m(x) and n(x)) to the mechanics of nanoscale fracture. To this end, we first note that the nanocrack propagates through thermally activated breakages of nanoparticle connections, which are governed by the corresponding interparticle potential function. Depending on the potential function, there could exist a zone attached to the tip of the traction-free crack, in which the interparticle resistance decreases with an increasing crack opening. This zone is commonly called the fracture process zone (or the cohesive zone), which is an essential concept in nonlinear fracture mechanics pioneered by Barenblatt in 1959 [20, 21]. If the structure size is considerably larger than the cohesive zone size, it is legitimate to replace the original traction-free crack and the cohesive zone altogether by an equivalent elastic crack with its tip located approximately at the middle of the cohesive zone [22]. This approach is usually referred to as the equivalent linear elastic fracture mechanics. We adopt this concept in the present analysis, and therefore the nanocrack considered in the model represents an equivalent elastic crack.
The purpose of the present study is to investigate the general probabilistic behavior of nanocrack propagation, which does not rely on any particular potential function. We note that each thermally activated breakage of nanoparticle connections can be described by a jump between two metastable states, which represent the states of the structure before and after crack propagation by one nanoparticle spacing (.ig. 2). The activation energy barrier Q to be overcome is much larger than the thermal energy b (k T = 0.025 eV at room temperature, where b k Boltzmann constant and Tabsolute temperature). Based on the transition rate theory, the frequency of the jump over the energy barrier can be described by the Kramers equation [18] =6.626 × 10 34 J s, h is Planck constant. When the structure is subjected to an externally applied stress, the two adjacent metastable states would exhibit an energy bias (.ig. 2). Clearly, this energy bias would favor crack propagation.
Recent studies have shown that the energy bias between two adjacent metastable states due to the applied stress is significantly smaller than thermal energy b k T [3,7]. Therefore, the Kramers equation is still applicable. Based on .ig. 2, we have .ig. 2. Thermally activated jumps over energy barriers (color online).
Since we use the equivalent linear elastic fracture mechanics to model the nanocrack tip, ∆Q(x) can be related to the energy release rate at the crack tip [3,7], i.e.
where 0 l is the characteristic size of the nanoscale element, E is the elastic modulus, b is the thickness of the element in the transverse direction, G is the energy release rate function, σ is the applied nominal stress, and g(z) is the dimensionless energy release rate function. Substitution of Eqs.
In the present analysis, we consider that the nanoscale element would attain the critical state of stability once the crack propagates to a length of . a l At that point, the crack would propagate dynamically. Therefore, we may limit our analysis to domain x ∈ [0, ].

Ek T cbl
Note that the dimensionless stress σ is similar to the Péclet number used in the convectiondiffusion equation, which determines the relative dominance of the convection and diffusion behaviors. To solve Eq. (9), we need to supplement it with the initial and boundary conditions. Here we consider that, at τ = 0, the location of the nanocrack tip is uniformly distributed over a normalized length α (α < 1), i.e. 1 where H(x) is the Heaviside step function. When the crack tip reaches the left boundary (ξ = 0), the tip must either stay at the boundary or move to the right. Therefore, we can prescribe a reflecting boundary condition at ξ = 0. To express the reflecting boundary condition, we first rewrite Eq. (9) as , , gc p p ξ ξ τ = σ ξ ξ τ − ξ τ is the flux of the occupational probability. The reflecting boundary condition implies zero flux of the occupational probability, i.e. 2 , (0) (0, ) (0, ) g p p ξ σ τ − τ = 0. Since the energy release rate vanishes at zero crack length, i.e. g(0) = 0, the reflecting boundary condition at ξ = 0 reads , (0, ) 0.
When the crack tip arrives at the right boundary (ξ = 1), the nanoscale element reaches failure. Therefore, we can model the right boundary as an absorbing boundary, i.e. p(1, τ) = 0. (13) Equations (10), (12), and (13) provide a complete set of initial and boundary conditions for Eq. (9), and we can solve the time evolution of the occupational probability p(ξ, τ) for a given applied stress .
σ The main interest of the analysis is the failure probability of the nanoscale element at time τ, which can be calculated as It is noted that Eq. (14) contains both the probability distribution functions of strength and lifetime. .or a given duration of applied stress, the failure probability is equivalent to the strength distribution while for a given applied stress the failure probability represents the lifetime distribution.
In this study, we apply the model to analyze the propagation of an edge crack in a large nanoscale element. The corresponding dimensionless energy release rate function can be written as g(cξ) = 1.12 2 πξ, and Eq. (9) becomes ) Q x k T ∆ ≤0.1 is satisfied for a case where the critical crack length is about 3000 atomic spacing. .or the initial condition, we choose α = 0.95 to represent a nearly uniform distribution of the initial crack size (note that the occupational probability must vanish at the absorbing boundary at τ = 0). It should be mentioned that, based on Eq. (9), the choice of the energy release rate function (or crack geometry) would not affect the qualitative dependence of failure probability on the applied stress and loading duration, which is the main interest of this study.
Equations (19) and (20) Coefficients n C can be determined from the initial condition (Eq. (10) Equation (24) gives the failure probability at the stress-free condition. .igure 3 shows that the present analytical solution and the numerical solution of Eq. (16) agree well with each other. This analysis shows that, through stress-free random walk, the nanocrack could reach its critical state in a finite duration. However, the typical diffusion time scale d t is much longer than the usual service lifetime of the structure (for example, the diffusion time scale d t is about 12 million years for 0 b Q k = 1.4 × 10 4 K and a critical crack length of 3000 atomic spacings), and as a consequence, for structures under either short-time strength test (510 min) or long-term service loading (e.g. 100 years), the probability of stress-free failure of a nanoscale element is negligibly small.

SIMULATION RESULTS AND DISCUSSION
By solving Eq. (15), we obtain the failure probability of the nanoscale element as a function of loading duration for a given applied stress (.ig. 4a), which represents the lifetime distribution. .igure 4b shows the lifetime distributions in the low failure probability regime. It is seen that, for a wide range of low failure probabilities (10 13 ≤ f P ≤ 10 6 ), the tail distribution of the lifetime approximately follows a power law, i.e. f ( ) | . P σ τ ∝ τ (25) Meanwhile, for a very short loading duration (e.g. τ = 10 10 ), the failure probabilities for the cases of σ = 1 and 2 are very close to each other. The present initial condition (Eq. (10)) implies that, when the loading duration is very short, the diffusion process contributes significantly to the failure probability. Consequently, the applied stress would have a diminishing effect on the failure probability. Note that the shortest loading duration used in the current calculation is τ = 10 10 . If we consider a even shorter loading duration, we may expect that the failure probability for the case of σ = 4 would also approach the probabilities for the cases of σ = 1 and 2.
The solution of Eq. (15) can also be presented as the failure probability of the nanoscale element as a function of applied stress for a given loading duration, i.e. the strength distribution. .igure 5 shows the calculated strength distribution for some values of τ that are relevant to engineering applications. .or example, τ = 10 5 represents approximately 100-year loading and τ = 10 10 represents a laboratory static fatigue test (e.g. 20 hours). It is seen from .ig. 5 that, for a wide range of low probabilities (10 13 ≤ f P ≤ 10 6 ), the strength distribution exhibits a power-law dependence on σ with an exponent of 2, i.e. 2 f ( )| . P τ σ ∝ σ (26) Note that for extremely low failure probabilities, the strength distribution must deviate from the power law because the structure has a finite failure probability even at zero stress, which corresponds to the diffusion failure probability fd ( ).
P τ An essential result of the foregoing discussion is that, in the low failure probability regime (10 13 ≤ f P ≤ 10 6 ), both the strength and lifetime distributions exhibit a power-law behavior but with different exponents. By combining these two tail distributions, we obtain a power-law tail distribution of the overall failure statistics as 2 f ( , ) P c σ τ ≈ σ τ for small values of 2 , σ τ (27) where c = const. .igure 6 shows that Eq. (27) agrees well with the simulated dependence of f P on 2 σ τ for a range of low failure probabilities that is relevant to engineering practice. It should be noted that, for typical loading durations (e.g. 10 10 ≤ τ ≤ 10 5 ), this applicable range of failure probabilities implies that the applied stress is not too small. The deviation from the power-law behavior at low streses can be seen from .igs. 4b and 5. Therefore, it is clear that Eq. (27) does not represent the extreme left-tail behavior of the failure statistics as σ → 0, which is beyond the concern of most engineering designs. The present model indicates that there does not exist any stress threshold below which the failure probability is zero. Instead, it predicts a finite failure probability (i.e. the diffusion failure probability fd ), P albeit negligibly small, even when the stress is absent.
In addition to the tail behavior of the failure statistics, we also investigate the mean failure behavior predicted by the model. .igure 7 presents the relationship between the mean failure time τ and the applied stress . σ It is seen that, when the applied stress is small, the mean failure time asymptotically approaches a constant, which represents the mean diffusion failure time. As the applied stress increases, the mean failure time decreases considerably. When the stress is sufciently .ig. 5. Tail failure probability of the nanoscale element in terms of the applied stress for τ = 10 10  , C σ τ = (28) which represents a power-law stresslife relationship. Equation (28) is analogous to the well-known Basquin law for fatigue failure [24]. Experiments showed that the stresslife curve of macroscopic failure of engineering ceramics also follows such a power-law form except that the power-law exponent of the applied stress is much higher [25]. The difference in the exponents for nanoscale and macroscale failures can be explained by a multiscale energetic model of fracture kinetics [2,3,26].
To further explore the mean stresslife behavior (Eq. (28)), we replot the failure probability in terms of a single random variable 2 , σ τ as shown in .ig. 8. It is found that, when σ ≥ 3, all the failure probability curves can be collapsed onto a single curve except for the extreme left tail if we use 2 σ τ as the abscissa. In other words, if we denote 2 σ τ using a random variable ζ, then almost the entire failure probability curve can be described by is the cumulative distribution function of random variable ζ. Though Eq. (29) does not apply to the extreme far left tail, it would give an accurate mean value of ζ, and clearly we would have 2 . C σ τ =

IMPLICATIONS .OR STRUCTURAL .AILURE PROBABILITY
The foregoing analysis (Eq. (27)) indicates that, for engineering practice, the tail of strength distribution of a nanoscale element can be considered to follow a power-law function, i.e. . P ∝ σ This power-law tail distribution has profound implications for reliability analysis of quasi-brittle structures. Recent focus has been placed on structures failing under controlled loads at macrocrack initiation, which represents one of the most dangerous failure types [13]. .rom a statistical viewpoint, this type of structures can be modeled as a chain of material representative volume elements (RVEs), and consequently the strength distribution of the structure can be calculated by using a finite weakest-link model [1,2,7]. This model represents statistically the damage localization mechanism, which governs the failure of quasi-brittle structures.
It is convenient to express the structural strength in terms of the nominal strength, which is defined as n m ( ) P bd σ = m (P is the maximum load capacity of the structure, d is the characteristic structure size, and b is the width of the structure in the transverse direction). The probability distribution of the nominal strength can be expressed by where 1 ( ) P x is the strength distribution of one representative volume element, and i c is the constant such that i c σ is the maximum elastic principal stress in the ith representative volume element.
Recent studies showed that the strength distribution of the macroscopic representative volume element can be related to the strength distribution of its nanoscale elements through a statistical multiscale transition model [13,7]. The model consists of a combination of bundles and chains arranged in a hierarchical manner (.ig. 9). The chain and bundle models capture, in a statistical sense, the damage localization and load redistribution mechanisms, respectively. It can be proven that, for the chain model, the tail of strength distribution would follow a power law if the strength of each element has a power-law tail distribution, and the power-law exponent remains the same. .or the bundle model, the power-law tail of strength distribution would also be retained however the exponent would increase in proportion to the number of elements in the bundle [3,7,27]. Therefore, based on the hierarchical model, the power-law tail distribution of the strength of nanoscale elements implies that the tail distribution of the representative volume element strength would also exhibit a power law but with a larger exponent, i.e. 1 0 ( ) ( ) m P s σ = σ (m > 2). Meanwhile, the hierarchical model predicts that the core of the strength distribution of marcoscopic representative volume element could be approximated by a Gaussian distribution. The Gaussian distribution is a result of the asymptotic behavior of the bundle model, which does not depend on the strength distribution of its element [3,7]. Therefore, the only feature of the strength distribution of nanoscale element that is transmitted through the scales is the power-law tail behavior, while the bulk part of the strength distribution of nanoscale element is inconsequential for the strength distribution of the macroscopic representative volume element.
Based on Eq. (30), the power-law tail distribution of the strength of representative volume element directly implies that the tail distribution of structural strength follows a power law with zero threshold. .igure 10 shows the measured strength distributions of several quasi-brittle materials in the Weibull scale. It is seen that the tail distribution follows a straight line,  , , m s σ are constants, and 〈x〉 = max (x, 0). In this model, the strength threshold 0 σ is introduced empirically in order to improve the fitting of the measured strength histogram over the two-parameter Weibull distribution. In recent studies [13,7], it has been demonstrated that the finite weakest-link model with zero strength threshold can provide optimum fits of the strength histograms of many quasi-brittle materials. The deviation of the measured strength histograms from the two-parameter Weibull model is caused by the finiteness of the number of representative volume elements in the structure, but not by the zero strength threshold. Therefore, introducing a finite strength threshold is unnecessary. .urthermore, it was shown that, at the large size limit, the three-parameter Weibull model yields an incorrect size effect curve of the mean structural strength leading to an unsafe design [7,31].
Since the overall strength distribution of the structure is modeled by a weakest-link model, the strength distribution of large-size structures is governed by the left tail distribution of the representative volume element strength, which is also commonly referred to as the domain of attraction [32]. According to the theory of extreme value statistics, the power-law tail distribution of 1 ( ) P σ indicates that the strength distribution of large-size structures must follow the two-parameter Weibull distribution [1,3,7].
The foregoing analysis is anchored by the finite weakest-link model, which assumes that the strength of each representative volume element is statistically independent. A more general approach for modeling the strength statistics of quasi-brittle structures is through the level excursion analysis [6,33]. The analysis considers a random strength field t ( ) f x as well as a random stress field σ(x) caused by the random microstructure of the material. By employing a nonlocal failure criterion, one may consider that the structure attains its peak load capacity once, at any material point, the nonlocal stress reaches the local material strength, i.e.
where Z(x) is the dimensionless random stress field such that σZ(x) represents the random field ( ) x σ of the nonlocal stress that is used for the failure criterion, and Ω is the structure domain. The use of nonlocality in Eq. (32) reflects the finite size of the localized damage zone. Eq. (32) indicates that the failure probability of the structure can be regarded as the probability at which the random field ( ) x σ exceeds t ( ) f x at least once, which is also called the first passage probability. The first passage probability can be calculated by modeling the random level excursion as a modified Poisson process [6,32,33].
The main advantage of the level excursion analysis is that it can account for both random strength and stress fields with their individual autocorrelation features. In recent studies [6,33], Eq. (32) was used to investigate the strength distribution of quasi-brittle structures. One essential result is that, if the underlying probability distribution function of random strength field t ( ) f x has a power-law tail, the strength distribution of the structure would also have a power-law tail [33]. Meanwhile, it was also shown that Eq. (32) predicts a twoparameter Weibull distribution of nominal strength for large-size structures, which agrees with the theory of extreme value statistics.
Based on the foregoing discussion, we may conclude that the power-law tail distribution of macroscopic structural strength can be directly related to the tail distribution of nanoscale elements. The present probabilistic modeling of nanocrack propagation further reveals that this power-law tail is closely tied with the behavior of thermally activated crack jumps, which can be described by the transition rate theory (Eq. (4)). This analysis provides a physical justification for the power-law tail distribution of quasi-brittle structures, which in turn validates the Weibullian strength distribution for large-size quasi-brittle structures (or, equivalently, brittle structures).

CONCLUSIONS
The thermally activated jumps of the nanocrack tip can be modeled, in a continuous limit, by the .okker Planck equation. The drift velocity is calculated by the net frequency of the crack jump. Based on the transition rate theory, the jump frequency can be directly related to the energy release rate function, which varies with the crack length.
In the absence of applied stress, the nanoscale element has a finite failure probability. However, for typical loading conditions, this diffusion-governed failure probability is far smaller than the risk level that is of practical interest. .or a wide range of low failure probabilities that are relevant to engineering practice, the failure statistics of the nanoscale element exhibits a power-law behavior in terms of both the applied stress and loading duration, i.e. 2 f . P ∝ σ τ When the applied stress is sufciently large, the mean failure behavior of the nanoscale element is featured by an inverse power-law relation between the applied stress and mean lifetime 2 . − τ ∝ σ This stresslife relationship is analogous to the Basquin law for fatigue failure.
The recently developed multiscale statistical model of quasi-brittle structures indicates that the power-law tail distribution of the strength of nanoscale elements directly determines the tail behavior of the strength distribution of macroscopic structures. The present nanocrack propagation model shows that the powerlaw tail distribution of structural strength, which is the key assumption of the Weibull strength statistics, can be explained by the thermally activated propagation of nanocracks.