Exact Solution to an Interacting Extreme-Value Problem: The Pure-Flaw Model

Simple models play a key role in the microstructural analysis of mechanical failure in composites and other materials with complex and often disordered microstructures. Although equal load-sharing-models are amenable to rigorous statistical analysis, problems with local load enhancements near failed regions of the material have so far resisted exact analysis. Here we show for the first time, that some of the simpler of these local-load-sharing models can be solved exactly using a sub-stochastic matrix method. For diluted fiber bundles with local load sharing, it is possible to find a compact expression for the characteristic equation of the sub-stochastic matrix, and from it obtain an asymptotic expansion for the largest eigenvalue of the matrix. This in turn gives the asymptotic behavior of the size effect and statistics of the fiber-bundle models. We summarize these results, and show that the important features of the exact result can be obtained from a single scaling analysis we had developed previously. We also compare the statistics of fracture with the appropriate limiting extreme-value survival distribution (a Gumbel distribution), and, as previously indicated by Harlow and Phoenix, note that the Gumbel distribution performs quite poorly in this problem. We comment on the physical origin of this discrepancy.


Introduction
It has been known since the pioneering work of Chaplin [1] and well known since the classic work of Griffith [2] that randomly occurring flaws or weak links effectively determine the observed tensile strength of materials.Early on it was realized that the dependence of failure upon the weakest part of a material structure gives rise to non-Gaussian statistics for fracture stress and strain.These developments lead to the classical period of the development of the statistics of extremes by mathematicians such as Dodd [3], Frechet [4], Fisher and Tippett [5], von Mises [6], Gnedenko [7], and Gumbel [8], Following the work of Duxbury et al. [9][10][11], there have been many attempts to use random network models to determine the statistics and size dependence of material breakdown [12][13][14][15][16].These calculations have in many cases elucidated the general behavior and size dependence of breakdown, but few exact results have been produced.
Perhaps the simplest model that shows the statistics of brittle failure has been the pure-flaw, chainof-bundles model of Harlow and Phoenbc [17] which has been studied by Harlow [18] and more recently by Hartow and Phoenix [19].In this model there is a series or chain of m structurally and statistically independent bundles of n elements each as shown in Fig. 1, where the vertically applied uniaxial stress is shared by the surviving vertical fibers (bonds).Each element or fiber is independently present with probability p and absent with probability/= 1-p.The survival probability of the chain of bundles is then the survival probability of a single bundle raised to the power m.The main difficulty in this analysis is calculating the survival probability of a single bundle.The extension of this theory to the survival of two-dimensional networks is straightforward and amounts to the approximation that cracks or flaws only exist and break along the direction transverse to the direction of the applied stress.Following Harlow and Phoenix, we assume the local-load-sharing mode! for a flaw (crack) of length n to be which is to say that the entire force applied to the cluster is concentrated at the tip (on the fibers adjacent to each end of the flaw or vacancy cluster).Failure of any surviving bond (all of which have the same strength) leads to a rip which causes failure of the entire bundle.Solution of this mode!requires finding the bond (weakest link) which experiences the largest stress enhancement and that stress which would break this most stressed bond.Duxbury, Leath, and Beale [11] showed how a one-dimensional model such as this could be used as a simple model for fracture or breakdown of a two-dimensional network.If one considers that cracks or flaws only exist and break horizontally, then the two-dimensional model becomes that illustrated in Fig. 2 (i.e., no horizontal bonds break).Then if we imj)ose spiral boundary conditions, where the last site in a row is connected to the first site of the previous row on the other side of the sample, then the NxN network problem is reduced to a one-dimensional chain of N^ fibers (or bonds) in parallel like that in Fig. 1. ■ Fig. 2. A two dimensiunal lattice with only horizontal Cracks.Spiral boundary conditions identify each site on the right edge with the site on the left edge of the previous row.

Single-Cluster Calculation: The Sub-Stochastic Transition-Matrix Method
In this calculation one assumes that the weakest link is the pure flaw or cluster of vacancies of the largest size that exists in the sample.The survival probability is then closely related to the probability that in a sample of length L{=N^) that there are no clusters of vacant bonds (flaws) of size greater than some prescribed value n.Using a generatingfunction technique Duxbury, L-eath, and Beale [11] calculated exactly the asymptotic form of the probability to be c,in)^[]-pr^y in the limit of large L, It is now possible to rederive this result while introducing the sub-stochastic transition-matrix method.Following Harlow [18] we define all possible endings of a fiber bundle of length L-t-1 and the way in which those endings may be generated from a bundle of length L and the probabilities of those endings.Since there are no allowed clusters larger than n, the allowed bundle endings or distinctive endings at a particular site are an occupied site (1) followed by 0^r^« vacancies (0) so that these distinct endings are spanned by the basis vector tj>L^ = {4>{ih'hm><i>om <^ia..(i))t where the last element contains n zeros.
Then the probability increment for going from a cluster of size r to a cluster of size r' on the next site is included in the matrix product Af<fe = which is the same as the matrix M operating; times on the probability vector 4\> for the starting site, or The probability CL(n) that there are only clusters up to size n in the entire bundle (or network) of size L is thus where the sum is over all the elements of t^L-The simplest and most natural boundary conditions are periodic ones where Ci.{n) becomes the trace since the 1st and L th sites must be the same, where AA are the eigenvalues of M. We find the eigenvalues of M via its characteristic equation where a =plk and b =f/X = (\-p)/X.A cofactor expansion of the determinant £>" about its last row, gives immediately the recursion relation where D"-i is the n X « determinant for clusters up to size (n -1), With Du = (a -1), the solution, upon iteration of Eq. (7a), is the characteristic equation Summing this geometric series we obtain This equation is the characteristic equation of M times (A -/), so there is an additional spurious root at/(since we are interested in the largest root, this does not affect the analysts).Since M is primitive and non-negative, its largest eigenvalue is real and distinct and it can easily be seen that all the eigenvalues are less than 1 and Amu approaches 1 for large n.Thus we set Amw = 1 -c and expand Eq. ( 8) to lowest order in e and /", which gives us A™«l-p/""'+C>Cr^).
In order to find the failure probability as a function of applied stress, we use the load-sharing rule Eq. ( 1), coupled with the fact that the failure of the bond carrying the largest local stress nucleates catastrophic failure, and thereby use the relation ah/a-l+-j^, (11) where a\, is the breaking strength of a single fiber.Note that we could have used a variety of other load-sharing rules here, and for example the same expression with n raised to an arbitrary power is also of physical interest.This result combined with Eq. ( 10) yields the probability S{a) that a fiber bundle will survive at stress a S(iT) = (l-pf'^ r (12) For large n and L, this becomes the modified Gumbel form, introduced previously [10,11] in the analysis of the random fuse network.Although Ci(n) in Eq. ( 10) becomes a Gumbel distribution in n, the substitution of «(o-) from Eq. ( 11) produces a modified Gumbel distribution that is significantly different from a Gumbel form in a in the high-reliability tail of the distribution.This modification is di.scussed further in Sec. 4,

Double-Cluster Calculation
Several authors [12,16,19,20] have suggested that the most critical defect is not a single cluster of n vacancies, but rather a double cluster (double colinear crack) of n vacancies separated by a single occupied site located anywhere within the H -I-1 adjacent sites.Such a double crack is shown in Fig. 3a.This candidate for the most critical crack is appealing because the stress enhancement at the interior occupied site grows as n in network models rather than as n"^ (as for the edges of a single crack in a two-dimensional network) and because the increased entropy of the (n +1) locations of the occupied site makes it more probable.Thus, following Harlow and Phoenix [19], we consider the probability of bundles of length L with repeated double cracks (and single cracks when the occupied site is at either end) not exceeding n vacant sites in any two adjacent cracks separated by a single site.These repeated double cracks are shown in Fig. 3b.Harlow [18] showed that this problem is amenable to analysis by the sub-stochastic transition-matrix method.The load-sharing rule is still given by Eq. (1) as before but now n is the sum of the number of vacant sites immediately on the left and right of any isolated intact bond or fiber.Thus it is necessary to keep track of not only the number of vacancies in the cluster being considered but also those in the previous vacant cluster.There are now (fi +l)(n -(-2)/2 distinct endings that must be considered at a site (or bundle ending); these are given by the basis vector ^ = <^i),<^io)...<^iu.,.o>; <^ioi).<^ioin)...^ioio..,o);</*(K»i)...<^iooio...n>;...;0(io...oi)' where there are no more than n total vacancies in any element.With this ordering of states the n =4 substochastic matrix for this problem, for example, is given by Eq. (13).
A great simplification in the characteristic equation W"<fr = A<^, (14) where ^ are the eigenvectors of M", is possible since most of the rows of M" contain only a single non-zero element /.This gives, for example /<^i) = A<^iQ, or &<^i) = (/"/A)^i) = <^(i(j)-By such relations, we can eliminate all the rows of M" except the rows with/>'s corresponding to the reduced basis vector <^ = <^i),<^iui),<^i«ii),...,i^ioo...oi).
where a =p/A, and b =//A.This A/' matrix can be considered as a new transition-matrix which adds a cluster at a time rather than a bond or fiber at a time.Thus we have the characteristic determinant equation.
forn ^3 and odd.While for n even S4, we find.

2* --z-z+'^-z '"^
1 -ab""^ ■ In Eqs.(20a) and (20b), The key quantity s is given by Eq. ( 19) above.Equations (20a,b) are the exact expressions for the characteristic equation of the original M in Eqs. ( 13) and ( 14).Again we find that the largest eigenvalue of Af is near 1 for large «, so we set A = 1 -e and expand Eqs.(20a,b) and find that in both cases, to lowest order in tand/".e = [(n+2)p = -/j]f-*'+(9(P"^), Comparing this double-cluster result to the singlecluster result Eq. ( 9) we find the expected (n +2) from the possible locations of the single bond in a double cluster of size (n +1).The (-p) is a correction to properly handle the single-cluster cases as well as the double-cluster case, since these are included whenever the isolated bond is located at one end of the double cluster.It is only important for small n.In order to check and better understand the asymptotic form Eq. ( 21) of ABU, and the importance of the other eigenvalues A, we have made several numerical evaluations of the various equations.First, we have numerically found the largest eigenvalue Amu of the original, full, substochastic transition-matrix M as given by Eq. ( 13).Using the unit vector as a starting vector we repeatedly apply the matrix M to it.Since the largest eigenvalue is unique, this process converges exponentially to the largest eigenvalue.We found in general that convergence occurred to six significant figures with at most 50 matrix products (even for matrices Af of dimension (n + ])(n +2)/2 = 10,000.The sparsity of M with this iterative procedure eliminated matrix-storage problems.The results of this iterative procedure for CLin) versus n are shown in Fig. 4 as solid lines, and the dots give Am« with Amax as given by the asymptotic form Eq, ( 21).Good agreement is seen for allp, with a small deviation in the p =0.2 data.However, a more stringent test is needed in the high-reliability (large n) tail.Thus, in Fig. 5 we plot the quantity ^j^=j^M(n+2y-pl The solid straight line versus (n +1) is the asymptotic result, which is linear in n and this is compared with the iterated numerical values (circular symbols), for/7=0.5 and 0.2.For large n, in all cases the two calculations agree.But for small values of p <0.5 there appears a minimum in (1 Amm)//^*' versus n which corresponds to higher order terms inf.In particular the next order term in A^" is 0(f^^) which would appear as a Oif'^') correction in Fig. 5. Finally, we test for the accuracy of neglecting all but the maximum eigenvalue A4«.It is possible to directly, numerically evaluate the trace of M'-by iteration since it only requires storage of the matrix and one vector at any time.Using this method, we have evaluated the quantity il-itr{M'-)y"-)/r\ (23) which should converge to (1 -Amw)//"*' when the largest eigenvalue is dominant.A numerical test of this convergence is shown in Fig. 6, for L = 1000, and shows that for large lattice sizes, the most important corrections to Ci(n) are the higher order contributions to Am^i, which are of 0(f^^), rather than the neglected smaller eigenvalues of M, which are relatively unimportant here.Finally, we obtain the asymptotic form Can)={1-[(n +2)p'-pT*' + 0(f'"'')}'-.(24) Thus, following the same arguments as at the end of Sec. 2, we find that to leading order the survival probability of the entire network or chain-of-bundles is which is again of the form of a modified Gumbel distribution with slightly different coefficients from Eq. ( 12), there are additive (log(logJL)) corrections in the double-cluster case.The logarithmic scaling law in turn implies that the failure statistics is of the double-exponential form given in Eq. ( 28).The WeibuH and Frechet distributions always give power-law size scaling.These qualitative arguments are very powerful and are confirmed by the rather elaborate, exact calculation described here.

Fig.
Fig. t.A one dimensional array of intact bonds (fibers) and flaws (vacancies).The tensile stress u, is applied vertically.

Fig. 3 .
Fig. 3. a) A double duster of €n toul vacancies plus one isolated occupied bond, b) Repeated, overlapping double clusters; each pair of clusters an indicated by the brackets contains € n total vacancies plus one isolated bond.