Stick-Slip Dynamics in Fiber Bundle Models with Variable Stiffness and Slip Number

We present an extension of fiber bundle models to describe the mechanical response of systems which undergo a sequence of stick-slip cycles taking into account the changing stiffness and the fluctuating number of slip events of local material elements. After completing all stick-slip cycles allowed, fibers can either ultimately break or can keep their final stiffness leading to softening or hardening of the bundle, respectively. Under the assumption of global load sharing we derive analytic expressions for the constitutive response of the bundle with both quenched and annealed disorder of the failure thresholds where consecutive slips occur. Our calculations revealed that on the macro-scale the bundle exhibits a plastic behavior, which gets more pronounced when fibers undergo a higher number of stick-slip cycles with a gradually degrading stiffness. Releasing the load a permanent deformation remains, which increases monotonically for hardening bundles with the maximum deformation reached before unloading starts, however, in the softening case a non-monotonous behavior is obtained. We found that the macroscopic response of hardening bundles is more sensitive to fluctuations of the number of stick-slip cycles allowed than of the softening ones. The quenched and annealed disorder of failure thresholds gives rise to the same qualitative macro-scale behavior, however, the plastic response is found to be stronger in the annealed case.


INTRODUCTION
Fibre bundle models (FBM) are one of the most important theoretical approaches to the damage and fracture of disordered materials [1]. In the framework of FBMs, the specimen is discretized as a bundle of parallel fibers which are subject to an external load along the fibers' direction [2,3]. The Young modulus of fibers is typically assumed to be constant so that materials' heterogeneity is entirely represented by the randomness of the strength of fibers. Even in their simplest form, FBMs provided a deep insight into the process of damaging of heterogeneous materials [2,4,5] making also possible to embed fracture processes into the general framework of statistical physics [1,6,7] and to clarify its analogy to phase transitions and critical phenomena [8][9][10][11][12][13].
Soon after the introduction of the basic concept of FBMs by Peires in 1927 [14], the model had been extended to capture time dependence and fatigue effects [15]. During the past decades subsequent developments of the model have demonstrated that varying the mechanical response [16] (brittle, plastic) and rheological (visco-elastic) behavior [17][18][19][20] of individual fibers, furthermore, the degree of strength disorder [21][22][23], range of load sharing (local, global) [11,24,25] following breaking events, and the way of loading [19,20,26] (quasi-static, creep, fatigue) the model is able to capture a broad spectrum of materials' behavior. Due to this flexibility, the model has gained a wide variety of applications from the fracture of fiber reinforced composites [25,27], through granular materials, where force chains were treated as load bearing fibers [28,29], to the rupture of biological materials [30].
Recently, we have proposed an extension of FBMs [31,32] to describe the mechanical response of systems with a complex micro-structure which respond to external loading by local rearrangements, for instance, of particles like in granular materials [28,33] and in agglomerates of dipolar particles [34], or by an activation of internal stored length such as spider silk [35]. A special area where stick-slip dynamics can be exploited for technological applications, is the development of nanocomposites, in particular carbon nanotube (CNT) composites, where CNTs are embedded in a polymeric matrix. One of the interesting properties of CNT nanocomposites is the ability to absorb vibrational energy which, together with high strength and fatigue tolerance, makes them perfect candidates for new multifunctional composite structures [36]. Recent experimental studies have revealed that the damping enhancement can be attributed to the CNT-matrix stick-slip caused by the severe elastic mismatch which leads to shear stress localization at the interfaces [37]. To represent the microscale complexity underlying stick-slip dynamics, fibers of our model were assumed to undergo a sequence of stick-slip events: when the local load reaches a failure threshold the fiber does not break, instead its load drops down to zero but the fiber has the ability of support load again. The model has been successfully applied to analyze the mechanical response of sheared granular materials [33] and the effect of root reinforcement on the stability of soils [29]. Complementing the stick-slip FBM with a healing mechanism it proved to be capable to describe the failure process of snow [38,39].
In the basic setup of the model it is assumed that fibers keep their original stiffness during their entire damage history [31,32], which is a crude simplification. It is reasonable to assume that restructuring events give rise to a degradation of the local stiffness or to hardening, e.g., due to damage or activating stored length inside the material, respectively. Additionally, fibers of the bundle are allowed to undergo the same number of slip events although in an extended sample the number of possible restructurings may have a spatial variation. To make the stick-slip FBM more realistic, in the present paper we resolve these limitations by allowing for the change of the stiffness of fibers as a result of slip events, and we capture the fluctuations of the number of stick-slip cycles allowed for fibers. Under the assumption of global load sharing, we derive analytic expressions for the macroscopic constitutive response of the bundle both for quenched and annealed disorder of failure thresholds where slip is activated. The stiffness change is treated in a multiplicative way which allows for a unified framework of stiffness degradation and stiffening. We demonstrate that the stick-slip dynamics results in a plastic behavior on the macro-scale and explore consequences of the new degrees of freedom of the model.

STICK-SLIP DYNAMICS WITH VARYING STIFFNESS
The model consists of N parallel fibers which are characterized by the same initial stiffness value E 1. Under an increasing external load σ fibers have a linearly elastic behavior up to a threshold load σ th . When the load on a fiber exceeds the failure threshold we assume that the fiber does not break, instead it slips and gets extended by increasing its equilibrium length until the load drops down to zero. Heterogeneity of the material is represented by the randomness of the slip thresholds, which are sampled from a probability distribution p(σ th ). The slip event is instantaneous in the sense that it does not take time, however, after slipping has been completed the fiber can support load again. As an important step of generalization of the model we let the fibers' stiffness change after slipping in a multiplicative way, i.e. the stiffness is updated as where a ≥ 0 is the stiffness parameter of the model. Note that the special case a 0 captures the immediate irreversible failure of the fiber right at the first slip which essentially results into the same dynamics as the classical fiber bundle model [2,4,40]. The parameter choice a 1 recovers the original stick-slip model where stiffness does not change during the loading history of fibers [31,32]. In our present study we focus on the parameter ranges 0 < a < 1 and a > 1, which represent stiffness degradation and stiffening of fibers following slipping, respectively. For practical applications of the model, stiffness degradation (a < 1) is typically caused by internal damage of material units represented by fibers. Stiffness increase (a < 1) occurs, for instance, in granular materials under compression where the restructuring of force chains may be accompanied by stiffening [28], and in biological materials like spider silk which respond to an increasing load by the activation of stored length [41].
After the slip has been completed, the fiber gets sticked so that it can support load again described by the constitutive law where ε denotes the strain of the fiber. Eq. 2 takes into account that the relaxed length of the fiber is extended with the strain threshold ε th σ th /E of slipping. This dynamics has the consequence that fibers may fulfill again the slipping condition and can eventually undergo a stick-slip sequence representing the gradual restructuring of the material. In order to describe such sequences we set the number of allowed slip events K ≥ 1, which is first assumed to have a fixed value for all the fibers of the bundle. It follows from Eq. 1 that after k stick-slip cycles the stiffness E ′ of the fiber has the value which can be greater or lower than the initial value E, for a > 1 and a < 1, respectively. It is a crucial question at which threshold loads the subsequent slip events occur. In the simplest case we can assume that the threshold value σ th is fixed for the entire history of a fiber which provides a representation of the quenched disorder of the materials' micro-structure. It can also occur that after restructuring events the local physical properties of the system change which can be captured by assigning a new threshold value to the fiber from the same probability distribution p(σ th ) each time a slipping occurred (annealed disorder).
For the interaction of fibers we assume global load sharing in such a way that fibers are stretched between two hard loading plates which ensure a global response of the entire system following local slipping events. However, global load sharing does not imply equal load of fibers in our stick-slip system since fibers which suffered a different number of slips have different stiffnesses and relaxed lengths, hence, they keep different loads. After fixing the type of disorder the threshold strains of slips can be determined from the corresponding stress thresholds. A fiber which has slipped k times at the consecutive strain thresholds ε (1) th , ε (2) th , . . . , ε (k) th up to the externally imposed strain ε, keeps the load where the sum of failure thresholds ε k th determines the relaxed length ε k 0 of the fiber. In the following we derive the macroscopic constitutive relation of the bundle for the cases of quenched and annealed disorder of slipping thresholds with global load sharing. Our main goal is to explore the consequences of the changing stiffness and of the fluctuations of the number of slip events fibers can experience.

Quenched Disorder of Failure Thresholds
Quenched disorder means that slips of a fiber always occur at the same stress threshold σ th , assigned to it in the initial state of the system. However, the corresponding threshold strains are not constant, which is illustrated in Figure 1, where the damage history of a single fiber is presented with a stiffness parameter a < 1. It can be observed that in spite of the constant threshold load σ th , the strain values ε (1) th , ε (2) th , ε (3) th , . . . where slips occur, gradually increase due to the degrading stiffness. It follows from Eq. 4 that the threshold strains ε (k) th (k 1, . . . , K) of consecutive slips of a fiber are determined by its initial strain threshold ε (1) th and by the stiffness parameter a of the model as The relaxed length ε (k) 0 of the fiber after the kth slip is the sum of all previous threshold strains which yields For the sum of the geometric series inside the brackets we introduce the shorthand notation S(a, k) so that Eq. 6 simplifies to ε (k) 0 ε (1) th S(a, k). Here the value of S(a, k) can be cast into the closed form S(a, k) (a − k − 1)/(a −1 − 1) for a ≠ 1. The above expressions are valid both for gradual degradation a < 1 and stiffening a > 1, resulting in an increasing and decreasing sequence of strain thresholds of slip events, respectively. Note that in the particular case of the original stick-slip model with a 1, the sum S(a, k) takes the value S(a, k) k.

Derivation of the Constitutive Equation With Varying Stiffness
To obtain a closed analytic form for the constitutive equation, we assume strain controlled loading of the bundle between two hard plates. At a given strain ε during the loading process, the bundle is a mixture of subsets of fibers which are either intact (no slip), or have suffered different number of slip events k, where 1 ≤ k ≤ K holds. Based on Eqs 4-6 the failure index k of fibers can be expressed in terms of their initial failure thresholds ε (1) th and the externally imposed strain ε as The value of a is set to a 0.8 so that gradual stiffness degradation occurs through the subsequent stick-slip periods. Since the stiffness is changing, the threshold strains of slipping are increasing ε (1) th , in spite of the fixed stress threshold σ th of the fiber.
where the failure index k 0 stands for intact fibers.
The macroscopic constitutive response of the system can be obtained by summing up the load σ k (ε) kept by the subsets of fibers which suffered exactly k slip events The partial loads σ k (ε) can be expressed in terms of the disorder distribution as where the first, second, and third terms provide the contribution of intact fibers, of the fibers which have suffered exactly 1 ≤ k < K slip events, and of the fibers which have completed all the K stickslip cycles, respectively. Note that the integration limits capture the separation of the subsets of fibers given by Eq. 7. In the limiting case of small deformation ε → 0 only the first term Eq. 9 has a finite contribution recovering the expected linear behavior σ(ε) ≈ Eε with the initial stiffness. In the opposite limit ε → + ∞ only the last term Eq. 11 survives which expresses that after the number K of allowed slip events the fibers still keep load and a linear behavior emerges With the asymptotic stiffness E a a K E. On the right hand side 〈ε (1) th 〉 denotes the average slip threshold in the initial state of the bundle. Between the two limits, at intermediate strains the second term Eq. 10 controls the macroscopic response as fibers gradually undergo more and more slip events.
The constitutive behavior of the stick-slip bundle is illustrated in Figure 2 for several values of the maximum number K of slip events at a fixed value of the stiffness parameter a 0.8. For the explicit calculations we considered exponentially distributed slip thresholds with the probability density function where the parameter λ is set to λ 1. This disorder distribution has the advantage that all the expressions of Eqs 9-11 can be obtained analytically. When presenting the results we rescale the stress σ and strain ε with the fracture strength σ 0 c and ε 0 c of a simple equal load sharing FBM of the same threshold distribution where fibers break when the load exceeds their strength. Note that the constitutive equation of this classical FBM coincides with σ 0 (ε) of Eq. 9. It can be seen that at intermediate strains a plastic behavior emerges which gets more pronounced for higher values of K, i.e. the hardening curves tend to asymptotic straight lines according to Eq. 12 which are preceded by longer and longer plateau regimes in spite of the degrading stiffness. The model can also account for the ultimate breaking of fibers after completing  Frontiers in Physics | www.frontiersin.org March 2021 | Volume 9 | Article 613493 the allowed number K of stick-slip cycles by omitting the term of σ K (ε) in Eq. 8. Fiber breaking implies that in the limit ε → ∞ the load bearing capacity of the bundle gradually diminishes so that σ(ε) → 0 follows, however, increasing the number of slip events K a broader plastic plateau emerges similar to the hardening bundles (see Figure 2).

Effect of the Varying Fiber Stiffness on the Macroscopic Response of Stick-Slip Bundles
The value of the stiffness parameter a has an important effect on the overall behaviour of the system. It can be inferred from Eqs 5, 6 that for the stiffening case a > 1 consecutive slip events rapidly follow each other so that the asymptotic regime Eq. 12 is reached at a relatively low strain. Decreasing a at a constant value of K, the plastic regime preceding hardening gets more-and-more extended (see Figure 3 for illustration). For the lowest values of the stiffness parameter a the hardening regime is approached only at very high strains which hinders the precise structure of the σ(ε) curves in Figure 3. That's why we further analyze the emergence of the plateau regime at low a values in Figure 4 using logarithmic scale on the horizontal axis presenting also the partial loads σ k (ε) of the subsets of fibers of different failure indices. It can be observed that lowering the stiffness parameter, the plastic plateau develops as a steady regime decorated with some oscillations. It can be seen that for higher a values, where the stiffness slowly degrades, the curves of σ k (ε) strongly overlap each other and their peak load rapidly decreases with k. As a decreases the steady state emerges because the σ k (ε) curves get more-and-more separated while their peak load increases approaching σ c of the classical FBM, which coincides with σ k 1 (ε). In the limit of a → 0 the steady stress regime disappears because the σ k (ε) functions become almost completely separated (see Figure 4D). The separation of the σ k (ε) curves in the low a limit also reveals that their functional form is essentially the same, determined by the constitutive equation of the classical fiber bundle model σ k 0 (ε). The overall shape of the constitutive curve of FBMs has a high degree of robustness for a broad class of disorder distributions [2], which implies the robustness of the macroscopic behaviour of the stick-slip FBM presented above, against the distribution of slip thresholds p(ε (1) th ). It is a very important consequence of the stick-slip dynamics that upon unloading σ → 0 the bundle a permanent deformation ε r remains which depends on the maximum deformation ε m reached before unloading sets on. If fibers keep their last stiffness value after K stick-slip cycles, the permanent deformation Frontiers in Physics | www.frontiersin.org March 2021 | Volume 9 | Article 613493 5 monotonically increases and it has an upper limit ε max r , which can be realized if the unloading starts along the asymptotic linear regime of a hardening bundle. The value of ε max r can be obtained from Eq. 12 using the condition σ(ε max r ) 0, which yields ε max r S(a, K)〈ε (1) th 〉.
For intermediate starting points ε m , the remaining deformation ε r can be determined by applying the condition σ(ε r ) 0 in the constitutive equation Eq. 8 taking also into account that along the unloading curve no slip events occur. This leads to the final form ε r K−1 k 1 a k S(a, k) εm / S(a,k) εm / S(a,k+1) ε 1 p(ε 1 )dε 1 + a K S(a, K) which converges to ε max r in the limit ε m → ∞. Note that the starting point of unloading ε m appears in the upper bound of the integrals. The permanent deformation ε r of hardening bundles is presented in Figure 5 for several values of the maximum slip number K. It follows from Eqs 14, 15 that for increasing the number of slip events K and decreasing the stiffness parameter a, hardening bundles store a higher plastic deformation, which is also supported by the numerical results of the figure.
Examples of unloading curves are presented in Figure 2 for the case of softening bundles, where fibers break after K slips. Since no slip can occur under a decreasing load, the unloading curves are always straight lines and their slope, i.e. the unloading modulus decreases with increasing ε m . The analytic expression of the unloading modulus E u coincides with the denominator of Eq. 15 multiplied by the initial stiffness E of fibers The expression shows that the unloading modulus E u at a given ε m is the weighted average of the moduli Ea k of the subsets of fibers with slip numbers k 0, 1, . . . , K, where the weights are determined by the distribution p(ε 1 ) of the slip thresholds. In the limit ε m → ∞, the unloading modulus E u converges to the asymptotic value E u → E a Ea K .
For the unloading modulus of softening bundles the last term inside the brackets has to be skipped to take into account the ultimate breaking of fibers. It can also be inferred from Figure 2 that the permanent deformation of softening bundles is not monotonous, i.e. the ε r (ε m ) function has a maximum and decreases when unloading along the tail of the softening regime of the constitutive curve σ(ε). This behavior can be realized in Figure 2 by the changing order of the end points of the unloading curves. Figure 5 demonstrates that even the maximum value of the remaining deformation of softening bundles is significantly smaller than the corresponding maximum permanent deformation ε max r of the hardening ones. It follows from Figure 4 that once a plastique plateau can be identified, its extension is practically the same for hardening (no breaking) and softening (breaking after K slips) bundles. At a given parameter set a, K, the extension of the plateau can be characterized by the asymptotic value ε max r of the remaining deformation Eq. 14 of hardening bundles. It is interesting to note that in the case of a > 1 the value of S(a, K) in Eq. 14 converges to S ≈ 1/(1 − a −1 ) for K → ∞, which implies a finite limit of the permanent deformation and hence, of the extension of the plateau, when fibers get stiffer after slip events. However, for the case of stiffness degradation a < 1, the sum S(a, K) does not have a finite limit for large K values giving rise to a monotonically broadening plateau as the maximum slip number K increases.

FLUCTUATING MAXIMUM NUMBER OF SLIP EVENTS
Fiber bundle models with stick-slip dynamics can be applied to understand damage accumulation and fracture in a large variety of systems where micro-scale restructuring plays a dominating role. However, in the present form of the model the maximum number of allowed stick-slip cycles is fixed for all the fibers, which is a very strong constraint and limits the applicability of the model under realistic conditions. In this section we extend the model to capture the effect of the fluctuating slip number K. We assume that K is a random variable which is sampled from a probability density function g(K) defined over the range 0 ≤ K < ∞ with the normalization condition Note that the possibility of fibers with K 0 is included, which implies that some fibers remain always intact during the loading process.
The constitutive equation σ(ε) of the stick-slip fiber bundle with fluctuating maximum number of slip events K can be obtained by averaging the contributions of subsets of fibers of constant K Eqs 9-11 with the distribution g(K) The first term of the right hand side of Eq. 3 represents the load kept by intact fibers taking into account that a fiber can be intact either because it can never slip K 0, or because it is damageable K > 0 but at the current strain it has not experienced any slip k 0. The second term stands for those fibers which have slipped exactly k times and still can undergo further stick-slip cycles k < K with the same maximum value K, while the last one is the contribution of the fibers which have completed all allowed stick-slip cycles and still support load with their final stiffness.
which we use for explicit calculations to demonstrate the outcomes of the generic derivations. Here the parameter 〈K〉 denotes the average of the slip number K inside the bundle. Assuming that flaws responsible for slip events occur in an uncorrelated manner along fibers, the Poissonian distribution controlled by the average 〈K〉 gives an adequate description of the statistics of the maximum slip number K [42]. Figure 6 presents a comparison of the constitutive curves of the stick-slip bundle with constant and Poisson distributed maximum number of slip events in such a way that the constant Ks are set to be equal to the average 〈K〉 of the Poissonian. It can be seen that the functional form of the two sets of constitutive curves is practically the same. For low values of K 〈K〉 the corresponding curves fall relatively close to each other, however, deviations increase with increasing 〈K〉. The reason is that the standard deviation of the Poissonian distribution growth as the square root of the average 〈K〉 √ so that at higher 〈K〉 the distribution Eq. 20 gets broader and the fluctuations of K become more relevant. The fluctuating K affects also the asymptotic form of σ(ε) of hardening bundles which yields and for the asymptotic Young modulus E a and for the maximum value of the permanent deformation ε max r , respectively. Note that the average value 〈a K 〉 is calculated with the distribution g(K) of the maximum slip number K, while the average threshold 〈ε (1) th 〉 is determined by the disordered strength of fibers p(ε th ).
The constitutive response of softening bundles with fluctuating maximum slip number can also be obtained from Eq. 3 by skipping the last term, which represents the contribution of fibers with failure index k K. Figure 7 compares the behavior of softening stick-slip bundles with constant and fluctuating maximum slip numbers where the average 〈K〉 was set to be equal to the fixed K values of the corresponding bundles. It is interesting to note that contrary to the hardening case, deviations of the corresponding constitutive curves are more pronounced for small deformations. The position of the maxima of the corresponding σ(ε) curves nearly coincide, however, bundles of fluctuating K keep a lower load than their constant K counterpart. As the deformation ε increases the curves cross each other indicating the higher load bearing capacity of bundles with fluctuating K. For large deformation ε → ∞ the stress must converge to zero with a constant maximum slip number, however, the final stiffness of bundles with fluctuating K is not zero. This is the effect of those fibers which are not allowed to slide K 0, resulting in a finite asymptotic stiffness E a g(0)E e −〈K〉 E. Increasing the average slip number 〈K〉 the fraction of unbreakable fibers exponentially goes to zero making the difference of constant and fluctuating K bundles larger for low values of K (see Figure 7).

ANNEALED DISORDER OF FAILURE THRESHOLDS
In materials with a complex micro-structure slip events may be followed by a change of the local material properties. In our fiber bundle model this behavior can be captured up to some extent by assigning a new threshold value to the fiber after each slip event from the same probability distribution p(σ th ). This type of annealed disorder leads to a constitutive behavior qualitatively similar to the case of the quenched one, however, with a more complicated dynamics. After k slip events the load σ kept by a fiber is given by Eq. 4 of the general model construction, however, now the stress threshold σ th is not fixed for the fiber, instead after each consecutive slip event a new threshold is generated σ (i) th from the same distribution p(σ th ). It has the consequence that the threshold strains ε (1) th , ε (2) th , . . . , ε k th of consecutive slips events can be obtained as where ε k σ (k) th /E are strain values which have the same distribution p(ε th ) for the entire history of the fiber. Eq. 24 shows that although the stress thresholds σ (k) th are generated with the same distribution, the corresponding strain values ε k still have to be transformed to obtain the strain thresholds, where fiber slips occur. The constitutive behavior of a single fiber and the relation of the variables ε (k) th and σ (k) th are illustrated in Figure 8. Based on the above expressions the constitutive equation σ(ε) of the entire bundle with a fixed number of allowed slip events K can be cast into the form where again the load bearing contributions of the subsets of fibers of different slip numbers k 0, 1, . . . , K are summed up: the first term represents the load kept by fibers which are intact k 0 at the deformation ε; the second one is the sum of the contributions of fibers which have undergone exactly 1 ≤ k < K slip events, and the last term captures the load bearing capacity of fibers after completing all the allowed K slip events. Note that the products of probability density functions occur due to the independence of consecutive slip events, additionally, the upper bounds of the integrals have also complex dependencies.
The constitutive responses of stick-slip bundles with quenched and annealed disorder are compared in Figure 9 for the case of hardening with fixed values of the number of slip events K. It can be observed that the two sets of curves have a qualitatively similar behavior. For K 1 there is no difference between the two types of disorder so that the corresponding curves must coincide, however, larger quantitative differences of quenched and annealed responses are observed for higher K values. The sequence of independent identically distributed failure thresholds of annealed disorder gives rise to broader plateau regimes with a higher value of the mean stress along the plateau. For large deformation both sets of curves converge to asymptotic straight lines whose slope E a a K E does not depend on the type of disorder. Due to the qualitative similarities of the macroscopic responses of quenched and annealed disorder stick-slip bundles, we skip the details of further comparisons.

DISCUSSION
We presented an extension of fiber bundle models of stick-slip dynamics incorporating the effect of stiffness change of the fibers after slip events, and the fluctuations of the number of stick-slip cycles fibers can experience under an increasing external load. Stick-slip dynamics implies that when the load of a fiber exceeds its local strength the fiber does not break, instead it slips which increases its relaxed length. As a consequence, the load of the fiber drops down to zero, however, the fiber retains its load bearing capacity. The stick-slip FBM has a high degree of complexity making it flexible to describe various materials' behaviors. In order to enhance the applicability of the model, we introduced a parameter which controls the change of fibers' stiffness after slip events in a multiplicative way allowing for both gradual degradation and stiffening. During their loading history fibers can undergo a sequence of stick-slip cycles. In our model construction the failure thresholds where slip is activated can be either fixed for the entire failure process of the fiber (quenched disorder), or it can be sampled from a probability distribution (annealed disorder) representing fixed structural disorder of materials and the effect of the local change of the material's behavior after slip events, respectively. As another novel element of our study, the maximum number of slip events is treated as a random variable inside the bundle sampled from a probability distribution. The total number of slip events allowed and the threshold loads where slip is activated are independent random variables. After completing all the stick-slip cycles a fiber can either keep its final stiffness or it can suffer ultimate breaking, which result in global hardening and softening of the bundle in the limit of high deformation, respectively. We analyzed the model in the mean field limit, i.e. global load sharing (GLS) was assumed. However, in stick-slip FBMs GLS does not imply equal load on fibers, because at a given strain fibers of the bundle can have different relaxed lengths and stiffness values. We derived closed analytic forms for the macroscopic constitutive response of the bundle both for quenched and annealed disorder of the slip thresholds. These results showed that on the macro-scale the bundle exhibits a plastic behavior, i.e. the σ(ε) curves develop a plateau regime which becomes broader with increasing number of slips. We analyzed in details the role of the changing stiffness of fibers in the emergence of the plastic plateau. Our calculations revealed that the stiffness parameter controls the degree of overlap of the contributions of fiber subsets of different failure index. In the limit of low stiffness parameter a ≪ 1, the plastic plateau is decorated with well separated maxima, while for slowly degrading stiffness a(1 and high values of the slip number K a smooth horizontal plateau is obtained.
Releasing the load on the bundle, a permanent deformations remains, which increases monotonically with the maximum deformation reached before unloading started for hardening bundles, while in the softening case a non-monotonous behavior is obtained. The permanent deformation of softening bundles proved to be smaller than that of their hardening counterpart. The asymptotic value of the permanent deformation of hardening bundles can be used to characterize the extension of the plateau regime at different stiffness parameters a. Our calculations revealed that increasing the number of allowed slips K the extension of the plateau has a finite limit for stiffening bundles a > 1, while it diverges for stiffness degradation a < 1.
We showed that fluctuations of the number of stick-slip cycles allowed for fibers affect the behavior of both hardening and softening bundles. Comparing hardening bundles of constant and fluctuating slip numbers with the same average value revealed that larger fluctuations result in a narrower plastic regime affecting also the asymptotic stiffness and the permanent deformation of the bundle. For softening bundles the FIGURE 8 | Constitutive behavior of a single fiber with annealed disorder of the slip thresholds. The value of the stiffness parameter is a 0.8 so that the fiber stiffness gradually decreases (compare to Figure 1). The consecutive thresholds σ (k) th are drawn from the same probability distribution p(σ th ), then the corresponding strain thresholds ε (k) th are obtained from Eq. 24. fluctuations of the slip number result in lower and higher load bearing capacities at small and large deformations, respectively, compared to their constant slip number counterparts. Annealed disorder of slip thresholds results in a qualitatively similar macroscopic response to quenched disorder, however, its description involves a higher mathematical complexity.
Here we focused on the effect of the changing stiffness and fluctuating slip number of fibers on the macroscopic behavior of a stick-slip FBM, and demonstrated that varying its parameters the model is capable to capture several aspects of the macro-scale consequences of stick-slip dynamics. The constitutive curves with the softening regimes and oscillations along the plateau can only be realized in strain controlled experiments. Under stress controlled loading slip events are followed by a load redistribution inside the bundle which can induce further slips and eventually can trigger an entire avalanche of slipping fibers until the bundle gets stabilized. As to the next, we are going to explore the effect of the varying stiffness and fluctuating slip number on dynamics and statistics of slip avalanches, which can have a relevance for the understanding of restructuring avalanches of granular materials and of earthquakes.
Although, the model is complex, still it could be further extended to fit to specific applications. For instance, when consecutive slip events result in the accumulation of internal damage of fibers, the degradation of fibers' stiffness may be accompanied by the reduction of fibers' strength. This correlation of local strength and stiffness can be captured by the model in such a way that at slip events fibers get a new failure threshold (annealed disorder) from a distribution which has the same functional form as the original one, however, its average is gradually reduced. Due to its flexibility, the model can serve as a starting point to develop more realistic micro-mechanical models of carbon nanotube reinforced polimeric composites where the high damping ability of the material originates from the stick-slip occurring at the CNT-matrix interface. Our approach can complement recently developed micro-mechanical models of CNT nanocomposites [43][44][45] providing an efficient framework for the representation of two sources of disorder (local strength and slip number), degradation and stiffening after slip events, and softening or hardening after the maximum slip number is reached. Work in this direction is in progress.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.