# Bounds for the time to failure of hierarchical systems of fracture

## Abstract

For years limited Monte Carlo simulations have led to the suspicion that the time to failure of hierarchically organized load-transfer models of fracture is non-zero for sets of infinite size. This fact could have a profound significance in engineering practice and also in geophysics. Here, we develop an exact algebraic iterative method to compute the successive time intervals for individual breaking in systems of height $`n`$ in terms of the information calculated in the previous height $`n1`$. As a byproduct of this method, rigorous lower and higher bounds for the time to failure of very large systems are easily obtained. The asymptotic behavior of the resulting lower bound leads to the evidence that the above mentioned suspicion is actually true.

preprint: HEP/123-qed  

There are few mechanical problems more complex and difficult to cast into a definite physical and theoretical treatment than the range of phenomena associated with fracture. Furthermore, there are few problems with a wider field of application: material science, engineering, rock mechanics, seismology and earthquake occurrence. Our understanding of fracture processes in heterogeneous materials has recently improved with the development of simple deterministic and stochastic algorithms to simulate the process of quasistatic loading and static fatigue . The load-transfer models used here are called fiber-bundle models, because they arose in close connection with the strength of bundles of textile fibers. Since Daniels’ and Coleman’s  seminal works, there has been a long tradition in the use of these models to analyze failure of heterogeneous materials . Fiber-bundle models differ by how the load carried by a failed element is distributed among the surviving elements in the set. In the simplest scheme, called ELS for equal load sharing, the load supported by failed elements is shared equally among all surviving elements. Important from the geophysical point of view is the hierarchical load-sharing (HLS) rule introduced by Turcotte and collaborators in the seismological literature . In this load-transfer modality the scale invariance of the fracture process is directly taken into account by means of a scheme of load transfer following the branches of a fractal (Cayley) tree with a fixed coordination number $`c`$. The HLS geometry nicely simulates the Green’s function associated with the elastic redistribution of stress adjacent to a rupture. In the static case, the strength of an HLS system tends to zero as the size of the system, $`N`$, approaches infinity, though very slowly, as $`1/log(logN)`$ . The dynamic HLS model was introduced in the geophysical literature in reference ; their specific aim was to find out if the chain of partial failure events preceding the total failure of the set resemble a log-periodic sequence . In the analysis of , it appeared that, contrary to the static model, the dynamic HLS model seemed to have a non-zero time-to-failure as the size of the system tends to infinity. In this paper we provide evidence that this conjecture is true by means of an exact algebraic iterative method where trees of $`n`$ levels, $`N=c^n`$, are solved using the information acquired in the previous calculation of trees of $`n1`$ levels. As a byproduct of this method, we obtain a rigorous lower bound for the time to failure of infinite size sets which results in being non-zero.

In  we showed how in fiber-bundle dynamical models, using a Weibull exponential shape function and a power-law breaking rule, one can devise a probabilistic approach which is equivalent to the usual approach  where the random lifetimes are fixed at the beginning and the process evolves deterministically. Without loss of generality, from the probabilistic perspective the set of elements with initial individual loads $`\sigma =\sigma _0=1`$ and suffering successive casualties is equivalent to an inhomogeneous sample of radioactive nuclei each with a decay width $`\sigma ^\rho `$; $`\rho `$ being the so called Weibull index, which in materials science adopts typical values between 2 and 10. As time passes and failures occur, loads are transferred, $`\sigma `$ changes and the decay width of the surviving elements grow. In the probabilistic approach, in each time step defined as:

$$\delta =\frac{1}{\underset{j}{}\sigma _j^\rho },$$
(1)   
one element of the sample decays. The index $`j`$ runs along all the surviving elements. The probability of one specific element, $`m`$, to fail is $`p_m=\sigma _m^\rho \delta `$. Eq. (1) is the ordinary link between the mean time interval for one element to decay in a radioactive sample and the total decay width of the sample as defined above. Due to this analogy, radioactivity terms will frequently be used in this communication. Note that we use loads and times without dimensions. The time to failure, $`T`$, of a set is the sum of the $`N`$ $`\delta `$s. For the ELS case, the value of $`\delta `$ as defined in Eq. (1) easily leads to $`T=1/\rho `$, which is the correct result for the time to failure of an infinite ELS set . In this communication we will apply these ideas to the HLS case obtaining the $`\delta `$s algebraically without having to use Monte Carlo simulations.

To give a perspective of what is going on in the rupture process of a hierarchical set we have drawn in Fig. 1 the three smallest cases for trees of coordination $`c=2`$. Denoting by $`n`$ the number of levels, or height of the tree, i.e. $`N=2^n`$, we have considered $`n=0,1`$ and $`2`$. The integers within parenthesis $`(r)`$ account for the number of failures existing in the tree. When there are several non-equivalent configurations corresponding to a given $`r`$, they are labeled as $`(r,s)`$, i.e, we add a new index $`s`$. The total load is conserved except at the end, when the tree collapses. Referring to the high symmetry of loaded fractal trees, note that each of the configurations explicitly drawn in Fig. 1 represents all those that can be brought to coincidence by the permutation of two legs joined at an apex, at any level in the height hierarchy. Hence we call them non-equivalent configurations or merely configurations. In general, each configuration $`(r,s)`$ is characterized by its probability $`p(r,s)`$, $`\underset{s}{}p(r,s)=1`$, and its decay width $`\mathrm{\Gamma }(r,s)`$. The time step for one-element breaking at the stage $`r`$, is given by

$$\delta _r=\underset{s}{}p(r,s)\frac{1}{\mathrm{\Gamma }(r,s)}$$
(2)   
This is the necessary generalization of Eq. (1) due to the appearance of non-equivalent configurations for the same $`r`$ during the decay process of the tree. In cases of branching, the probability that a configuration chooses a specific direction is equal to the ratio between the partial decay width in that direction and the total width of the parent configuration. We will compute at a glance the $`\delta `$s of Fig. 1 in order to later analyze the general case. To be specific, we will always use $`\rho =2`$. For $`n=0`$, we have $`\mathrm{\Gamma }(0)=1^2`$ and $`\delta _0=\frac{1}{1^2}=1=T`$. For $`n=1`$, $`\mathrm{\Gamma }(0)=1^2+1^2=2`$, $`\delta _0=\frac{1}{2}`$, $`\mathrm{\Gamma }(1)=2^2`$, $`\delta _1=\frac{1}{4}`$ and hence $`T=\frac{1}{2}+\frac{1}{4}=\frac{3}{4}`$. For $`n=2`$, $`\mathrm{\Gamma }(0)=1^2+1^2+1^2+1^2=4`$, $`\delta _0=\frac{1}{4}`$, $`\mathrm{\Gamma }(1)=2^2+1^2+1^2=6`$, $`\delta _1=\frac{1}{6}`$; now we face a branching, the probability of the transition $`(1)(2,1)`$ is $`\frac{4}{6}`$ and the probability of the transition $`(1)(2,2)`$ is $`\frac{2}{6}`$, on the other hand $`\mathrm{\Gamma }(2,1)=2^2+2^2=8=\mathrm{\Gamma }(2,2)`$, hence $`\delta _2=\frac{4}{6}\frac{1}{8}+\frac{2}{6}\frac{1}{8}=\frac{1}{8}`$. Finally $`\delta _3=\frac{1}{16}`$ and the addition of $`\delta `$s gives $`T=\frac{29}{48}`$.

Now we define the replica of a configuration belonging to a given $`n`$, as the same configuration but with the loads doubled (this is because we are using $`c=2`$). The replica of a given configuration will be recognized by a prime sign. In other words $`(r,s)^{}`$ is the replica of $`(r,s)`$. Any decay width, partial or total, related to $`(r,s)^{}`$ is automatically obtained by multiplying the corresponding value of $`(r,s)`$ by the common factor $`c^\rho =2^\rho =4`$. This is a consequence of the power-law breaking rule assumed. The need to define the replicas stems from the observation that any configuration appearing in a stage of breaking $`r`$ of a given $`n`$, can be built as the juxtaposition of two configurations of the level $`n1`$, including also the replicas of the level $`n1`$ as ingredients of the game. In Fig. 1, one can observe the explicit structure of the configurations of $`N=4`$ as a juxtaposition of those of $`N=2`$ and its replicas. From this perspective we notice that as the configurations for the height $`n`$ are formed by juxtaposing two configurations of the already solved height $`n1`$ and its replicas, with the restriction that the total load must be equal to $`2^n`$ because the total load is conserved, the single-element breaking transitions occurring can only be of three types. One case a) corresponds to the breaking of one element in a half of the tree while the other half remains as an unaffected spectator. Another case b) corresponds to the decay of the last surviving element in a half of the tree, which provokes its collapse and the corresponding doubling of the load borne by the other half. Finally, the third case c) corresponds to the scenario in which one half of the tree has already collapsed and in the other half one break occurs. In this third case the decay width is obtained from the information of the replicas. This holds for any height $`n`$, allows the computation of all the partial decay widths in a tree of height $`n`$ from those obtained in the height $`n1`$ and will be illustrated in Fig. 3.

Using henceforth a convenient symbolic notation, in Fig. 2a we have built what we call the primary width diagram for $`n=3`$ resulting from the juxtapositions of pairs of configurations of the previous diagram of $`n=2`$ and its replicas. It is implicitly understood that time flows downward along with the breaking process. In Fig. 2a the boxes represent $`n=3`$ configurations, and at their right, the two quantities in parentheses indicate the two $`n=2`$ configurations forming them. The value of the partial width connecting a parent and a daughter is written at the end of the corresponding arrow. The sum of all the partial widths of a parent configuration in a branching is always equal to the total decay width, $`\mathrm{\Gamma }`$, of the parent. From this primary width diagram one deduces the probability of any primary configuration at any stage $`r`$ of breaking, and consequently $`\delta _r`$ is obtained using Eq. (2). Finally, by adding all the $`\delta `$s we calculate $`T(n=3)`$. To illustrate the connection between the explicit configurations as those of Fig. 1 and the somewhat hermetic notation of Fig. 2, in Fig. 3, we have shown three explicit examples relating them. Fig. 3 is selfexplanatory. The three examples correspond to the cases a), b) and c) mentioned before.

By iterating this procedure, that is by forming the primary diagram of the $`n+1`$ height by juxtaposing the configurations of the primary diagram of the height $`n`$, we can, in principle, exactly obtain the total time to failure of trees of successively doubled size. Two examples are $`T(n=3)=\frac{63451}{123200}`$ and $`T(n=4)=\frac{21216889046182831}{46300977698976000}`$. The problem arises from the vast amount of configurations one has to handle in the successive primary diagrams. This fact eventually blocks the possibility of obtaining exact results for trees high enough as to be able to gauge the asymptotic behavior of $`T`$ in HLS sets. That is why, taking advantage of the general perspective gained with the exact algebraic method, henceforth our more modest goal will be to set rigorous bounds for $`T`$. This task is much simpler.

To set bounds, we will define effective diagrams in which for each $`r`$ there is only one configuration which results from fusing in some specific appropriate way all the configurations labeled by the different $`s`$ values; see Fig. 2b. These effective configurations are connected by effective decay widths denoted by $`a_r`$. Thus, the effective diagram for any $`n`$ is calculated from its primary diagram, and then is used to build a primary diagram of the next height $`n+1`$. For $`n=0`$, 1 and 2, $`a_r=\mathrm{\Gamma }(r)`$ because $`r`$, the $`\mathrm{\Gamma }(r,s)`$ are equal to a unique value $`\mathrm{\Gamma }(r)`$. The point is how to define $`a_r`$ for $`n3`$, so that all the $`\delta (r)(n4)`$ and hence $`T(n4)`$ are systematically lower (or higher) than its exact result. This goal is easily accomplished using,

$$a_r=\mathrm{\Gamma }_{max}(r)(or\mathrm{\Gamma }_{min}(r))$$
(3)   
i.e., by assuming that the configuration of maximum (minimum) width saturates the single element decay of the stage $`r`$. Better rigorous boundings are obtained using the geometric mean (or the harmonic mean), namely

$$a_r=\underset{s}{}\mathrm{\Gamma }(r,s)^{p(r,s)}(or\frac{1}{\underset{s}{}p(r,s)\frac{1}{\mathrm{\Gamma }(r,s)}})$$
(4)   
The bounds obtained from these formulae are plotted in Fig. 4. It is clear that those based on Eq. (3) are poor, in fact the lower bound goes quickly to zero. On the other hand, those based on Eq. (4) are stringent and useful. The value obtained for the lower (higher) bound to $`T`$, from Eq. (4), will be called $`T_l`$ $`(T_h)`$. The arithmetic mean, i.e., $`a_r=\underset{s}{}p(r,s)\mathrm{\Gamma }(r,s)`$ also leads to a lower bound but it is not as good as that coming from the geometric mean. The bounds emerging from these three means are rigorous because when configurations of different decay width are fused, the form of the generated function appearing in the calculus of the $`\delta `$s are concave (convex). This will be explained elsewhere.

We have fitted the data points of $`T_l`$ by an exponential function of the form $`T_l=T_{l,\mathrm{}}+ae^{b(nn_0)}`$, $`T_{l,\mathrm{}}`$, $`a`$, $`b`$ and $`n_0`$ being fitting parameters. The success of this fitting is crucial because this exponential decay to a non-zero limit is the hallmark of our claim. By choosing subsets formed only by points representing big systems we observe a clear saturation of the asymptotic time-to-failure towards $`T_{l,\mathrm{}}=0.32537\pm 0.00001`$. An analogous analysis of $`T_h`$ leads to $`T_{h,\mathrm{}}=0.33984\pm 0.00001`$. Similar fittings of the Monte Carlo data points are inconclusive, due to the intrinsic noisiness of the MC results and the limited size of the simulated sets ($`N<2^{16}`$ elements). What this result implies is that a system with a hierarchical scheme of load transfer and a power-law breaking rule $`(c=2,\rho =2)`$ has a time-to-failure for sets of infinite size, $`T_{\mathrm{}}`$, such that $`0.32537T_{\mathrm{}}0.33984`$. Thus, there is an associated zero-probability of failing for $`T<T_{\mathrm{}}`$ and a probability equal to one of failing for $`T>T_{\mathrm{}}`$. The critical point behavior is thus confirmed. Invoking conventional universality arguments one deduces that this property holds for a hierarchical structure of any coordination. The case of dynamical HLS sets using an exponential breaking rule will be reported shortly.

A.F.P is grateful to J. Asín, J.M. Carnicer and L. Moral for clarifying discussions. M.V-P. thanks Diego Vázquez-Prada for discussions. Y.M thanks the AECI for financial support. This work was supported in part by the Spanish DGICYT.