Cracks in random brittle solids: From fiber bundles to continuum mechanics Sylvain Patinet, Damien Vandembroucq, Alex Hansen, Stéphane Roux To cite this version: Sylvain Patinet, Damien Vandembroucq, Alex Hansen, Stéphane Roux. Cracks in random brittle solids: From fiber bundles to continuum mechanics. The European Physical Journal. Special Topics, 2014, ￿10.1140/epjst/e2014-02268-9￿. ￿hal-01226393￿ HAL Id: hal-01226393 https://hal.science/hal-01226393 Submitted on 9 Nov 2015 HAL is a multi-disciplinary open access L’archive ouverte pluridisciplinaire HAL, est archive for the deposit and dissemination of sci- destinée au dépôt et à la diffusion de documents entific research documents, whether they are pub- scientifiques de niveau recherche, publiés ou non, lished or not. The documents may come from émanant des établissements d’enseignement et de teaching and research institutions in France or recherche français ou étrangers, des laboratoires abroad, or from public or private research centers. publics ou privés. EPJ manuscript No. (will be inserted by the editor) Cracks in random brittle solids: From fiber bundles to continuum mechanics Sylvain Patinet1,a, Damien Vandembroucq1,b, Alex Hansen,2,c, and Stéphane Roux3,d 1 Laboratoire PMMH, ESPCI/CNRS-UMR 7636/Univ. Paris 6 UPMC/Univ. Paris 7 Diderot, 10 rue Vauquelin,F-75231 Paris cedex 05, France 2 Institutt for fysikk, NTNU,N-7491 Trondheim, Norway 3 LMT-Cachan, ENS-Cachan/CNRS/PRES UniverSud Paris, 61 Av. Président Wilson, F- 94235 Cachan cedex,France Abstract. Statisticalmodelsareessentialtogetabetterunderstanding oftheroleof disorderinbrittledisorderedsolids. Fiberbundlemodels playaspecialroleasaparadigm,withaverygoodbalanceofsimplicity andnon-trivialeffects.Weintroducehereavariantofthefiberbundle model where the load is transferred among the fibers through a very compliant membrane. This Soft Membrane fiber bundle mode reduces totheclassicalLocalLoadSharingfiberbundlemodelin1D.Highlight- ing thecontinuum limit of themodel allows to compute an equivalent toughnessforthefiberbundleandhencediscussnucleationofacritical defect. The computation of the toughness allows for drawing a simple connection with crack front propagation (depinning) models. 1 Introduction Itisneedlessto emphasizethe importanceofbrittle fractureforsafetyandreliability issues. In fact, the condition for crack propagation is well understood since Griffith pioneeringwork[1].Criticalenergyreleaseratesorequivalentlystressintensityfactors are known to be material characteristics defining the onset of crack propagation for an existing crack. However, a resisting issue concerns crack initiation. Indeed, for a longperiod,thecommonviewwasthatcrackswouldinitiatefrompreexistingdefects. Brittlenesswasinvokedto considerthatassoonas initiateda crackwouldpropagate inanunstablefashiondowntocompletefailure.Thus,amazingly,solidmechanicsfor the most part disappeared from strength predictions. Because the initiation defect was local, the failure load remained a statistical issue, with probability distribution whichsimply resultedfromthat of the defect distributionafter an elementaryelastic computation (together with some assumptions for the defect geometry). This was formalized within “weakest links” type of approaches, that can be traced back to a e-mail: [email protected] b e-mail: [email protected] c e-mail: [email protected] d e-mail: [email protected] 2 Will beinserted by theeditor the pioneering work of Weibull.[2] These approaches had also the merit of shedding somelightonsalientfeaturesunexpectedinthemechanicsofhomogeneousmaterials, such as the occurrence of systematic size effects. A recent review on this topic and extensions thereof can be found in Ref. [3]. Alternative approaches based on a purely deterministic picture was proposed by Francfort and Marigo[4], as an extension of Griffith’s theory. This theory, although appealing, raisessubtle mathematicalissues — which may in fine resume to physical questions — related to how smallscales may or may not disappear atthe continuum limitforsuchsingularproblems.Anotherdirectionofattackconsistsofintroducinga morecomplicatedpictureforthecracktip,througha“fractureprocesszone” whereby, atmicroscopicscales,adisplacementdiscontinuitycanstillallowforthetransmission of stress. Such models following the initial work of Barenblatt [5] illustrate the fact thatinitiationandpropagationcanbeincorporatedinthesamepictureifinaddition totoughness,acriticalstressisalsoconsidered.Thisisyetanotherwaytostressthat some microscopic length scale will survive in this macroscopic picture. Leguillon [6] for instance proposes such a dual criterion for accounting for crack initiation in the spirit of these fracture process zone models. Although this latter class of theories emphasizes the importance of small scales, materials are always treated as if they were homogeneous. Thus, from a continuum mechanics point of view, the choice is left between ap- proaches that either emphasize the geometry of defects (without mechanical mani- festation other than linear elasticity) or focus on fine details of local variations of mechanical free energy (not to say damage) but ignores material variability. It is therefore interesting to consider discrete models (where no ambiguity is left about microscopic limits) and where disorder is explicitly considered. On the one hand, such models may be considered as too simplistic for portraying any specific materials, however, on the other hand, they may shed some light on statistical size effects that may be expected. This is very important as scaling is not neutral. If interpolation based on Weibull distribution is fair when the same sizes and level of probabilities as determined in identification, it becomes very fragile when extrapola- tion is needed. Hans Herrmann[7–9]has been one of the pioneers in this matter with firstpapers on electrical analog, “fuse networks”, to brittle material published in 1985. However, a direct numerical simulation reveals to be a difficult pathway to reaching definitive conclusions,asresultsarepronetoslowcorrectionsto scaling,andbecomeexpensive in computational power as one moves to two or three dimensions. [10] Nevertheless, those simulations could demonstrate quite clearly that the weakest link approach could not be applied at the smallest scales. Therefore, it was already clear from the start that neither approach based on prior defect statistics only nor on a more sophisticated interaction model but within a homogeneous medium description were applicable,andhencethemeritofthestatisticallawsandtheattachedcorresponding size effects could be questioned. To progress along these lines, fiber bundle models [11,12] play a key role. Indeed despite their simple definition, and their introduction about 90 years ago, an ex- haustive characterizationof their behaviorhas only beenachievedinsome particular cases[13–15]. These models consist of a collection of elastic brittle fibers with random strength loaded in parallel. One very important feature of these models leading to different scaling regimes is the way load is distributed amongst the unbroken fibers. In the “Equal Load Sharing” (or ELS), [12] load is equally distributed over all surviving fibers. This corresponds to clamping all fibers onto a rigid substrate onto which the external load is applied. Because of this even distribution, the dimensionality of the Will beinserted by theeditor 3 substrate plays no role. This limit can thus be seen as a mean-field model which is analytically tractable. The opposite limit of Local Load Sharing, (or LLS), is based on a local redis- tribution of the load carried by a fiber onto its unbroken neighboring fibers [15,16]. HerethebreakdownscenarioisverydifferentfromELS,inthesensethatafterafirst diffuse damage regime, (where fibers breakage is spread over the entire system), a localized region hosts a crack which — after its nucleation — grows without limits. Macroscopically,thisfailureappearsasbrittlebutafterasignificantnumberoffibers have been broken. However, the LLS model suffers from some ambiguity on the way theloadisredistributedafterfiberbreakage.Thisambiguitybecomesevenmorepro- nounced in two dimensions and more. In one dimension, some aspects of this model areanalyticallytractable[15,17].However,mostly,numericalmethods mustbe used. There are also other more realistic variants of the fiber bundle model. The Soft Clamp Model [18–21]has the fibers attached to clamps that respond elastically. The clamps are infinite half spaces and the force distribution among the fibers are due to the elastic response of the clamps. The transition between ELS and LLS has also been studied by varying contin- uously the sharing rule trough a simple power law stress redistribution function σ r−γ, where r is the distance from the broken fiber [22–24]. It was shown that ∼ the stress-transferfunction decay exponent dramatically impacts breakingprocesses, ultimatestrengthandcreepregime.InthelimitingcaseofELS,thebreakingoffibers is a completely random nucleation process while, for LLS, it becomes a nucleation problem. In the present study, we modify the Soft Clamp Model by replacing the infinite halfspaceclampsbythinelasticmembranes.ThiscorrespondstotheLLSlimit,that is to say an infinite exponent for the stress-transferfunction (γ ) in [22–24].We →∞ notethatbothinthe SoftClampmodelandthe Membranemodel, oneofthe clamps canbereplacedbyaninfinitelyrigidone.Forsymmetryreasons,thismodelwillhave the same behavior as were the two clamps equal. Otherphysicallyrealisticmodelsorlimits couldhavebeenconsidered.Inparticu- lar,anelastic plate describedthroughLove-Kirchhofftheoryhas naturallyabending stiffness. Indeed, there is no objection to considering such a model, (and it has al- ready been considered for instance in [25]) however, in 1D where the plate reduced to an elastic beam, it does not coincide with the LLS model precisely because of the flexuraleffectandcontinuityofthe rotation.This introducedadampened oscillatory behavior close to a crack. Because our aim is here to introduce a natural extension of LLS model in 2D that suffers no ambiguity in its definition, we will stick to the definition of the Soft Membrane limit. Let us also mention a related model, namely that of shear delamination of a thin film from a rigid substrate as introduced by Zaiser et al.[26–28]. A modeling proposed in [26] shows that the in-plane displacement obeys an anisotropic Laplace equation, that can be somehow compared to the Poisson equation ruling the out- of-plane displacement for the membrane under distributed load. The source term that is present in the Poisson equation is the result of the external loading and the balance equation leading to damage and crack propagation. A similar effect is naturally present in other modelings [18,19,25,29]. Since it controls the instability of crack propagation which is the hallmark of the LLS model, it is the analogy with models that have no distributed loadings that should be considered carefully. Inthepresentstudy,weshowthatitispossibletodefineandcomputeatoughness allowing for understanding the post-peak behavior, i.e. after the maximum load has beenreached.Moreover,thislimitisalsoshowntomatchthatofacrackfrontpinning model. The simplicity of the model finally allows us to study its behavior in the 4 Will beinserted by theeditor continuumlimitviaaperturbativeapproach,drawingpromisingpathfromnucleation to crack propagation in brittle heterogeneous materials. z u = −1 z f x f x a) broken fibers b) f Lap[u] = −1 ϕ=∆u.n z y n x c) Fig. 1. Load sharing rules in a fiber bundle between clamps that are elastic membranes. The lower membrane is infinitely stiff, whereas the upper membrane responds elastically. A macroscopic load F is applied to the system through a homogeneous pressure onto the membrane and is shared between unbroken fibers subjected to local forces f. a) When the uppermembrane is stiff, the model behavesas theELS model. b) In one dimension, a very soft upper membrane (or rather, string) will cause the model to behave as the LLS model. c) In two dimensions, a very soft upper membrane will transfer an additional load only to fibersadjacent toalready brokenones. 2 Soft Membrane Fiber Bundle in One Dimension ThefiberbundleinitiallyconsistsofN fibersdistributedhomogeneouslyoveralinear domain of size L = N, at positions labeled by integer x coordinates. To reduce the importance of edge effects, periodic boundary conditions are implemented. Each fiber is assumed to be elastic-brittle. Initially, all fibers are present. Then, through the application of load, fibers will be broken one at a time. A sequence index, t, is introduced which simply counts the number of brokenfibers. For simplicity, t will be referred to as “time”. Let us emphasize that the model itself does not depend on the physical time. In the elastic regime, a macroscopic load F is applied in the system, and is shared between unbroken fibers as illustrated in Fig. 1. A fiber at position x is subjected to a local force f(x,t). If the fiber is broken, then f(x,t) = 0. If not, because of elasticity, f(x,t) is proportional to the external load, f(x,t) = Fϕ(x,t), where ϕ(x,t) characterizes the load redistribution on the fiber, and that depends on the distribution of broken fibers and hence on t. Each fiber is characterized by its strength, f (x), independent of time. At each instant, the macroscopicbreaking load c isthevalueofF suchthatoneandonlyonefiberreachesitsthreshold,f(x,t)=f (x), c hence f (x) c F (t)=min , (1) c x ϕ(x,t) Will beinserted by theeditor 5 which corresponds to an extremal dynamics. Theloadsharingruleis whatallowsus todefine ϕ(x,t). Forinstance,ELSmeans that all surviving fiber support the same force as shown in Fig. 1a. Overall balance requires that ϕ(x,t) = 1, and hence ϕ(x,t) = 1/(N t) if the fiber is unbroken, − and 0 else. The membPrane model is equivalent to the LLS rule in 1D when the membrane (string)isverysoft.WenowdiscusstheLLSruleinsomedetail.Asabovementioned that balance imposes ϕ(x,t) = 1. In 1D, it is rather straightforward to arrive at the fact that for any interval of m consecutive fibers broken, an additional force of m/2N is transfered toPthe two surviving fibers nearest neighbor to the end of this interval. This comes in addition to the 1/N force that it always present. Note that a fiber may be neighborof two intervalsof m andm failed bonds respectively for left l r and right, and hence its load partition coefficient will be 2+m +m l r ϕ(x,t)= (2) 2N The Soft Membrane model is now considered. In this case, the deflection u(x) of this string under the applied loading has to be solved. If we consider that the fibers are much stiffer than the string, if the fiber is present at position x, then u(x) = 0. Otherwise, Lap[u](x)= 1 (3) − whereLapis the discretelaplacianLap[u](x)=u(x+1) 2u(x)+u(x 1).The load transferredby a membrane to its boundary writes fT =− u n where−n is anouter −∇ · unit vector, normal to the membrane boundary (see Fig. 1). In 1D, this means that a fiber locatedat xreceivesanadditionaltransferfT(x)=[u(x+1) u(x)]+[u(x − − 1) u(x)] where we separated the contributions of the left and right neighbors. In − this case, the load applied to a surviving fiber x (which is such that u(x)=0) is ϕ(x,t)=(1/N) 1+ u(y) (4) ! y n.n.x X where a summation over nearest neighbors (n.n.) of x is involved. The occurrence of the second order differential operator Lap comes from a conservation (balance) equation so that the total load applied onto the system is finally transferred to the unbroken fibers. The deflection u(x) of a membrane of extremities x and x writes i f u(x)=(x x )(x x)/2 and the load transfersfT(x )=fT(x )=(x x 1)/2. i f i f f i − − − − For a fiber surrounded by two intervals of m and m bonds, one thus recovers the l r verysame loadpartitioncoefficientthan inthe 1D LLSfiber modelgivenabove(Eq. (2). Figure 2a shows the classical 1D model load versus time t. We considered here a uniformdistributionoffiberstrengthin [0;1].Thecriticalloadisreachedatt=t∗ = 1256.Thus,theweakestlinkassumptionisindeedfarfrombeingapplicable.Figure2b showsthe redistributionfunctionatthe peakforce.We seethatthe largestfiberload ϕ comes from an interval of 10 consecutive broken fibers. 3 Soft Membrane Fiber Bundle in Two and More Dimensions Thegeneralizationofthe1Dsoftmembranemodeltotwodimensions(2D)isstraight forward.ThesystemisaL Lsquare,sothatthenumberoffibersisN =L2.Afiber × locationisdesignatedeitherbya2Dvectorxorbyitsintegercoordinates(x,y).The question is now to define the redistribution function. 6 Will beinserted by theeditor 0.12 6 0.1 5 0.08 4 ) F(t)c0.06 ∗ϕ(x,t23 0.04 1 0.02 0 0 0 1000 2000 3000 0 2000 4000 6000 8000 10000 t x Fig. 2. (a) Load versus time in 1D (N = 10000 and periodic boundary conditions). The critical time is the one where the maximum force is reached, t∗ = 1256. (b) Local force distribution ϕ(x,t∗) at the critical load. The membrane deflection, u(x), obeys u(x) =0 intact if fiber at x is (5) Lap[u](x)= 1 broken (cid:26) − (cid:26) where Lap[u] is again a finite difference laplacian operator Lap[u]= (u(y) u(x)) (6) − y n.n.x X The redistribution function on a surviving fiber in x is simply ϕ(x,t)=(1/N) 1+ u(y) (7) ! y n.n.x X i.e., precisely the same expression as in the 1D model. Figures3,4and5illustratethe loadversustime,thedisplacementandlocalforce per fiber at the critical state, and the failure time map at a late stage.It is apparent that, as expected, a critical defect is nucleated and finally grows with no possible arrest. This growth takes the form of a circular crack as seen in Fig. 4. Let us stress againthatthefractionofbrokenfibersatthepeakstressisaratherlargeproportion, t∗/N 23%. Again the weakest link approach is inapplicable. Simultaneously, when ≈ considering the overallmacroscopicforce versus time, Fig.3a, there are no fore-signs of nucleation of the critical defect. This observation can be rationalized, in the spirit of homogenization, by consid- ering a uniform equivalent threshold, feq. Considering a circular cluster of broken c bonds, the load supported by the broken fibers πR2F is distributed over the edge, leading to a force per fiber of R(t)F/2 and hence F (t) = 2feq/R(t). The radius is c c related to time as πR(t)2 =t t . Thus 0 − 2√πfeq F (t)= c (8) c √t t 0 − Figure 3b shows that this analysis holds after nucleation of the critical defect, as 1/F (t)2 has indeed a lower envelope obeying an affine law. The straight line fit c allows to define an effective threshold equal to feq 0.34. c ≈ Will beinserted by theeditor 7 0.2 2500 0.15 2000 F(t)c 0.1 21/F(t)c11050000 0.05 500 0 0 9000 9500 10000 0 2000 4000 6000 8000 10000 t t Fig. 3. (a) Load versus time for the Soft Membrane fiber bundle model in 2D. (b) Post- criticalanalysisofcriticalloadversustime.Theredlinebelongstothecircularcrackanalysis. 50 50 100 100 150 150 200 200 50 100 150 200 50 100 150 200 Fig. 4. Membranedisplacementu(x)inthe2DSoftMembranefiberbundleat thecritical load(left),andatalaterstage(right).Thecenterofthecriticaldefectcanbeguessed from thelate stage. Theperiodic boundary conditions are also visible. 5 50 4 100 3 2 150 1 200 0 50 100 150 200 Fig. 5. Local redistribution function ϕ(x,t∗) at thecritical load. 4 Measurement of Toughness Havingunderstoodthelatescenarioofacrackpropagation,onemaystudyitdirectly: a specific crack pattern already present in the initial state can be designed. For in- stance, half of the domain (say 0 y < L/2) consists of broken fibers in the initial ≤ state. The problem is thus initially similar to the 1D case. The membrane deflexion 8 Will beinserted by theeditor field is u(x,y)=(1/2)y(L/2 y) (9) − and hence the force acting along the edges is the total load (1/2)L2F divided by the front length (on both edges) 2L or FL/4. It can also be computed as Lap[u]. In our case, for L = 200, the force acting on the front edge is 50 times higher than in the bulk. Thus most (but not all) of the activity is taking place at the front. Note, however,that as the front propagatesas a whole, the overloadat the front increases, andhencethediffusedamageregimeisonlytakingplaceinaninitialtransientregime, which vanishes as L . →∞ Periodic boundary conditions being implemented, one may stop one of the two fronts(y =Linourcase)byinsertingalineoffiberswithveryhighstrength.Figure6 shows the displacement field and force acting on the fibers after some progression starting from the previous initial conditions. It is observed that a slightly roughened front progress from y = L/2 onward. Figure 7 shows the map of the time at which the fibers were broken. One can easily observe indeed the above mentioned front propagation, together with some “dust” which corresponds to the early damage of bondsbeingbrokenaheadofthefront(inthepresentcase,allbondshavingastrength less than 2/L 0.01 are broken in a first stage). Note that the dynamics of the ∼ preexistingcrackdoesnotproceedbyavalanchessincewe performextremaldynamic simulations, i.e. each “time” step only one fiber is broken. 20 20 40 40 60 60 80 80 y100 y100 120 120 140 140 160 160 180 180 200 200 50 100 150 200 50 100 150 200 x x Fig. 6. Displacement u(x,y) (left) and local force f(x,y) (right) after some progression. Figure 8 shows the load versus time. One observes again clearly the first initial (diffuse) damage regime where fibers are not located at the front. The interesting regime is the front progression stage. In the continuum, if the fiber strength feff c were uniform, one wouldexpect a uniformprogressiony =y +t/L (where y =L/2 0 0 in ourcase),so that the forceexertedatthe frontwouldbe (y +t/L)F/2and hence 0 the critical load would be 2feff F (t)= c (10) c (y +t/L) 0 Henceonecandefinean(instantaneous)effectivebondstrength,feff,fromthecriti- c calforcethroughF (t)(y +t/L)/2.ThisquantityisshowninFigure8.Themaximum c 0 value of the latter strength appears as a well-defined quantity which is the macro- scopic effective strength. The histogram shown in Fig. 8 is another way of revealing the threshold. It allows us to estimate feff 0.46 (11) c ≈ Will beinserted by theeditor 9 100 150 t 200 50 100 150 200 x Fig. 7. Map of failure time where thedownwards front progression is easily observed. as comparedto the elementary (mean field) approximationof f =0.5 that ignores c h i the effect of load redistribution. 0.5 250 0.4 200 efff(t)00..23 n(f)110500 50 0.1 0 0 0 0.1 0.2 0.3 0.4 0.5 100 110 120 130 140 150 f t Fig. 8. (a)Effectivefiberstrengthfeff,deducedfrom thecritical loadbyamultiplication c of half the mean crack length. (b) histogram of the instantaneous critical load, where the asymptotic toughness is theupperbound of this distribution. It is noteworthy that this estimate is also slightly larger than the one estimated from the random nucleation case. This may be due to the fact that a somewhat larger amount of damage has been experienced in the plain 2D lattice at the onset of nucleation. Estimating this damage at this point as roughly D = 0.23 (D is the fraction of missing bonds), leads to an effective decrease of the effective strength (when averagedover the front) down to (1 D)feff 0.34 very consistent with the − c ≈ previousestimate.To be complete, oneshouldalsocorrectthe aboveestimate by the damage level of the front type geometry, but because of the large system size, the damage is roughly D =0.01, and does not affect our estimate.