Random fuse model in the presence of self-healing

Self-regeneration is a fundamental property of biological materials, leading to enhanced mechanical strength and toughness if subjected to stress and fatigue. Numerous efforts have been devoted to emulate this property and various self-healing materials have been designed with the aim of a practical adoption in construction and mechanical engineering. To achieve this, it is important to understand how damage evolution and fracture propagation are modified by self-healing and to evaluate how mechanical behaviour is affected before failure. In this paper, we implement for the first time a self-healing procedure in the random fuse model, whose characteristic scaling properties have been widely studied in the literature on damage evolution modelling. We identify some characteristic signatures of self-healing, showing that it can delay the failure of a material undergoing loading, but it also lead to a hard-to-predict, more catastrophic breakdown.


Introduction
An interesting feature of many biological materials is the capacity to repair their microscopic traumas occurring during fatigue cycles and stress events [1]. This property is called self-healing, which leads to an increase in the time to complete fracture. Moreover, self-healing combined with a multi-level hierarchical structure, which is a recurring feature in biological materials, allows to improve the global mechanical properties of the system [2,3].
These observations have triggered research on artificial self-healing materials, aiming to improve the materials adopted in civil and structural engineering, with a potentially significant impact on many practical applications. From the seminal work of White et al [4], by now many techniques have been introduced to design and manufacture artificial self-healing materials [5][6][7][8][9][10][11][12][13][14]. A comprehensive review can be found in references [15][16][17].
The interest for these systems is not only limited to the investigation of artificial self-healing materials, but also to the understanding of the behaviour near the final failure of structures undergoing continuous cycles of fatigue, damage and reparation. The unexpected collapse of large scale structures (e.g. [18]), demonstrates the importance of investigating damage evolution in structures of relevance for civil and structural engineering.
Since the first experimental realisation of a self-healing material in 2001 [4], relatively few theoretical/ numerical studies of self-healing systems have been performed in the literature. Self-healing processes have been included in the classical fibre bundle model (FBM) [19,20] and its extension to the hierarchical case [21,22]. The self-healing formulation of FBM has captured in general the synergy between self-healing and hierarchy that leads to an increase of strength and fracture toughness [23].
Other approaches have concentrated on specific aspects of experiments, e.g. the modelization of the fracturing of micro-capsules containing a healing agent, and subsequent flow of the latter. For example, Verberg et al used a hybrid approach with a coupled lattice Boltzmann model and a lattice spring model (LSM) to describe the motion of microcapsules on a substrate with an adhesive coating under the effect of an imposed flow [24]. Bluhm et al also proposed a modelling approach for self-healing polymer composites based on the theory of porous media, considering embedded catalysts, the liquid healing agents, the solid healed material and the gas phase in the cracks [25]. Li et al presented a model to derive effective mechanical properties for the microcapsules in cement-matrix self-healing materials [26], while Mauludin et al treated the problem using cohesive elements [27]. On the other hand, Salib et al developed a hybrid computational approach using LSM and the hierarchical Bell model to investigate the mechanical properties and self-healing behaviour of nanogel particles connected by stable and labile bonds [28], and the same approach was used to address the problem of designing strong and tough biomimetic polymer networks with the capability of reforming links in their chain [29]. Guo et al focused on simulating the crack healing mechanism in hydrogels with chemical and physical crosslinks [30], and Alsheghri et al proposed a cohesive zone damage-healing model to be implemented in a finite element scheme to describe self-healing polymers [31]. Other approaches have concentrated on modelling of self-healing at the molecular scale, e.g. Maiti et al [32] used coarse-grained molecular dynamics to study the behaviour of self-healing polymers and compute parameters such as the local elastic modulus, reaction rates and cure kinetics, for a continuum macroscopic scale model. A review of early attempts of applying numerical methods to self-healing materials is given in [33], and a more recent review on computational modelling of hierarchical materials, including self-healing, is provided in [34].
Our objective here is to adopt a different approach to most of the cited works, similar to that in [35], i.e. to address the modifications of the fracture mechanics behaviour of a self-healing material, including from a statistical standpoint, on a macroscopic scale, independently of the specific healing mechanism at play. This is important to understand the structure and the mechanical behaviour of biological materials subjected to fatigue loading and characterised by self-healing and hierarchy, e.g. human tendons [36,37], but modelling predictions can also help the design of artificial self-healing materials [38]. For these objectives, it is fundamental to focus on the improvement of the material mechanical properties, in particular the overall strength and toughness, the statistical properties across the hierarchical levels and the effectiveness of different self-healing mechanisms. Much work remains to be done in this field, in particular to understand how self-healing modifies the characteristic scaling laws of the system and its fracture processes. These features are fundamental in the practical application of self-healing materials.
We thus adopt an approach analogous to that presented in [23,37], except the that numerical evaluation tool here is the random fuse model (RFM). The RFM has been extensively used in the literature to model damage evolution [39,40], and provides scaling laws as a function of the system size that are well-known [41,42] in the absence of self-healing. A limited number of works, starting from Cowie [43,44], also considered repair effects in the RFM, mostly in connection with the post-slip healing of fracture faults in earthquakes or fatigue failure of asphalt [45], but not in self-healing materials as such. Our aim is thus to understand the modification of the RFM scaling laws by introducing a self-healing mechanism whose observation may help to understand the universal signature of self-healing processes. Since it has been already demonstrated that self-healing improves the global mechanical strength and toughness of the system, we will focus on other effects occurring before failure, showing that the final rupture of a self-healing system can occur in a more catastrophic manner, despite the longer lifetime.
The paper is organised as follows: in section 2, we introduce the RFM algorithm with the self-healing extension. In section 3, we present the results for various observables, in particular in section 3.1, the RFM characteristic current-voltage curve, which is analogous to the stress-strain curve in mechanical systems; in section 3.2, the results for maximum current and voltage, in 3.3 results of avalanche distributions and in 3.4 the crack roughness. Finally, in section 4, a final discussion of results is presented.

Random fuse model
We implement the RFM in a square lattice of L × L resistors with conductivity c=1 ( figure 1(a)). As shown in the figure A1, the choice of a square lattice instead of a diamond lattice, used elsewhere in the literature, does not modify the overall qualitative behaviour. A voltage difference V is applied along the horizontal direction between the two sides of the lattice, while periodic boundary conditions are set in the vertical direction. Each resistor has a breaking threshold x j on the current, randomly extracted at the beginning of the simulation from a Weibull probability distribution where the scale parameter is set to λ=1 and the shape parameter to k=2.
The algorithm is the following: starting from zero, the voltage is progressively increased by steps of ΔV. At each step, the Kirchoff system of equations is solved to find the local currents I j flowing in each resistor j and, consequently, the total current I flowing through the lattice. If some local currents I j exceed the thresholds, the resistor with maximum ratio I j /x j is eliminated and its conductivity is set to c=0. The system of equations is solved again with the same voltage to find how the currents are redistributed. This procedure is repeated until there are no more resistors whose threshold is exceeded. After this, the voltage can be increased by another step ΔV and the whole algorithm is repeated. The overall system fails when a fracture of broken resistors extends across the whole lattice and the total current drops to zero. This algorithm is implemented using an in-house developed C++ code. The Kirchoff system of equations is solved by inverting the matrix of the linear system using the C++ library Armadillo, directly embedded in the code, adopting a LU decomposition algorithm. In the following simulations, we set the voltage step ΔV=0.001. We have verified that this does not affect the results and the duration time of the simulations. If ruptures do not occur at some potential, the system of equations is not solved again; and since it is linear, currents in all fuses are scaled proportionally by the same factor of the potential, i.e. (V+ΔV )/V.
Thus, in a single simulation we obtain a characteristic function I(V ), as in figure 1(b), showing the total current as a function of the voltage (corresponding to stress as a function of strain). From this curve, many observables can be extracted: the peak value is the maximum current I max (strength), the corresponding voltage of the peak is V(I max ), while the final voltage is V max . The total dissipated energy E is the area under the curve. We denote with n tot the number of total fuses of the system, which can be calculated as = n L L 2 tot 2 , and n b is the number of broken fuses at a given voltage. Thus, we define the total damage of the lattice as d=n b /n tot .
To gain insight in the statistical properties of the system, one can study the avalanche distribution: the number s of broken resistors at a given voltage is the size of a single avalanche event, so that a histogram displaying the number of avalanches of size s can be obtained from a single simulation. A reliable estimate of all these observables can be obtained through repeated simulations, and we denote by á ñ ... this statistical average.

Self-healing
We now suppose that the broken resistors can be restored by some self-healing process. We are addressing the self-healing process from a generic theoretical point of view, without considering any specific mechanism. The focus is on the effect of the self-healing process, i.e. the repair of a local microcrack, modelled as the reattachment of a broken link of the lattice, and on the time delay of this repair. These two features are common to all self-healing processes, with different parameter values, and can be included in the framework of the RFM, independently of the specific features of the practical method adopted to achieve self-healing in applications. We define the self-healing rate η as the ratio between the number of restored resistors, namely n sh , and of broken ones n b , h = n n sh b . Thus, for η=0 the case without self-healing is recovered, and for η=1, all fuses are repaired. In our model, the self-healing rate η corresponds to the probability that a broken resistor will be repaired.
In [23], link restoration was assumed instantaneous, while in this paper we assume that there can be a delay in the regeneration of a broken resistor, i.e. broken links are restored after a certain amount of time. Since the RFM is a quasistatic model, there is no real time evolution, but the applied voltage can be considered equivalent to time: if a time unit is associated with a voltage step, simulations with a linearly increasing voltage correspond to a linearly increasing time (in the RFM, any rupture and current redistribution is considered instantaneous with respect to the voltage step). Thus, a delay time can be expressed in voltage units of the simulation. If a resistor is broken at voltage V 0 and its regeneration delay is V sh , it is restored when the applied voltage is V 0 +V sh . In a real life phenomenon, this is connected to a characteristic time of the self-healing process. For example, in some applications, microcracks are repaired by a healing agent contained in small capsules inside the material itself. In this case, the delay time would be the characteristic time scale for the healing agent to polymerise and close the microcraks. This time scale is present in principle in all self-healing processes, so it can be added phenomenologically in the model from a generic point of view.
We assign the delay of each resistor by extracting it from a Weibull distribution W(x, Δ sh , k sh ), where the scale parameter Δ sh is the characteristic delay and k sh is the shape parameter. When a resistor is restored, its conductivity is set again to c=1 and a new breaking threshold must be extracted from the Weibull distribution of the thresholds W(x, λ ′ , k), but the scale parameters must be rescaled so that the resistor is not instantaneously burnt. Thus, we set l l ¢ = + I j 0 , where I j 0 is the current flowing through it at the moment of the previous breakage. In this way, the global statistical properties of the resistors are unchanged by self-healing. A resistor can be broken and restored any number of times. The rest of the algorithm works as in the case without self-healing.
In this paper, we focus on the effects of the delay in the case η=1, i.e. all resistors can potentially be restored. If not explicitly stated, all simulations are therefore performed with η=1. The self-healing delay, expressed in voltage units, can be compared with the average value of V max of the system. For l á ñ V sh max  the case of 'fast' self-healing occurs, i.e. a broken resistor has a high probability of being restored before the failure of the system, and vice versa, for l á ñ  V sh max the case of 'slow' self-healing occurs. We will explore both of these limits. It is convenient to define the healing rapidity G º D 1 sh , so that the case without self-healing is obtained for Γ=0. In [23], an analytical calculation was performed to describe the evolution of damage in a FBM with selfhealing in the limit of the Daniel's theory [19,20], and to provide an approximate expression for its strength. The resulting expressions are approximated, since no account for damage or healing localisation is present, but they still provide a useful analytical benchmark. According to [23], the relation between stress σ (current I) and strain ò (voltage V ) can be written as: where ò 0 and m are the scale and the shape parameter of the Weibull threshold distribution, respectively, and E is the Young's modulus. In the presence of a delay, a modified version of equation (1) can be provided if an effective strain-dependent self-healing rate η eff is adopted in place of η, accounting for an average characteristic time delay Γ −1 in the healing process. In particular, results are fitted with a good approximation by the phenomenological law: which provides a simple expression with the correct asymptotic behaviour: the case Γ=0, i.e. an infinite time delay, is equivalent to the case without self-healing, and the case G  ¥ is equivalent to an instantaneous selfhealing as in [23]. For intermediate values of the delay, the maximum strength is reduced with respect to the instantaneous self-healing.

Characteristic I(V ) curve
The maximum voltage and the peak current both increase with self-healing. This is shown in figure 2 by the characteristic I(V ) curve for various self-healing rates. In particular, by increasing the value of Γ, the curves approach the ideally linear case I∝V, in which the last catastrophic fracture event, leading to the failure of the system, occurs immediately after the peak current and involves a large avalanche of ruptures. Thus, a signature of the self-healing process in the RFM is an increase in maximum strength and ultimate strain, but also a vanishing plastic phase and perfectly brittle behaviour with a catastrophic final event.
To be more quantitative, in figure 3(a) we report the behaviour of the mean total damage á ñ d of the system as a function of applied voltage. In standard RFM simulations, the damage approximately increases linearly up to the final larger event. For slow self-healing values, i.e. D á ñ  V sh max , the curve lies over the standard case. For fast self-healing the slope of the curve progressively decreases and the spike due to the last catastrophic event increases, and for the fastest self-healing case the spread of damage is inhibited up to failure. In figure 3(b), we report the statistics of the number of broken fuses and restored ones at the end of the test. By increasing Γ, the number of broken fuses decreases and most of them are fractured at the last rupture event. These results suggest that with a fast self-healing rate the system damage may appear limited despite a catastrophic rupture avalanche is imminent.
The abrupt failure of the system is not due to a weakness of the restored links, since on average only about ten percent of the final broken links have been previously restored. Links are repaired in an unstressed state, so that they have larger thresholds. Thus, we must conclude that with a self-healing delay Γ2, the system may display a limited damage, but be in a critical state in which a single rupture can trigger the total failure.
We observe also that an imminent failure cannot be easily predicted from the statistical distribution of the currents. We define the effective current as sum over all unbroken links of the local currents normalised by their threshold values, i.e with the notation of section 2 = å -I n n I j x j j eff tot b 1 ( ) ( ) ( ). This is an average measure of the imminence to the failure of the links. Another observable is also the fraction of links whose ratio I( j)/x( j) is above the percentage p, namely r p , so that, for example, r 0.95 indicates a link current of 95% of its threshold. This is an estimate of how many links are very close to failure. For both observables, results before the failure in the presence of self-healing display only a small percentage increase with respect to the standard case, as shown in figure 4, so that apparently these kinds of observables are not effective in indicating that the system is in a critical state close to failure. Average global quantities cannot capture this phenomenon because, in this case, the  fracture propagates due to local current spikes (i.e. stress concentrations) around links that are progressively broken by the final avalanche. We have verified that this behaviour is not an artifact due to the choice of the elementary step ΔV. Thus, we must conclude that this is a critical feature of the system for large self-healing values.

Peak current and maximum voltage
The scaling law of the peak current with the total size of the lattice in RFM is consistent with á ñ = + a I cL c max 1 2 , where c 1 and c 2 are constants and the exponent estimate is α ; 0.96 [46], obtained both with diamond and triangular lattices. Our data in the standard case, obtained using a square lattice in the range 12L160, leads to α=0.954(1), which is consistent with the value of the literature. The other constants are c 1 =0.308(2) and c 2 =0.487 (10).
A similar scaling law holds for the maximum voltage, which may be associated with the average lifetime of the system. For 30L160, we find a scaling law á ñ = + b V cL c max 3 4 , with c 3 =0. 39(1), c 4 =2.2(1) and exponent β=0.941 (8). Tables of datasets used for these fits are reported in the tables A1-A3. These results are stable, i.e. removing the larger lattices from the fitted range does not modify the results. This demonstrates that finite size effects are negligible.
In the presence of self-healing, these scaling laws are modified and cannot be fitted with a simple power-law behaviour. In figure 5, we report the maximum current increase as a function of Γ for different lattice sides L. For slow self-healing values, the increase of the maximum current with respect to the case without self-healing follows a power-law h D º á ñ -á ñ~G where the exponent value is a ; 2. This power-law behaviour in the presence of self-healing is found in the approximated theory [23,47], and the quadratic behaviour is consistent with equations (1), (2) for small Γ. The collapse on the same curve for all L is clear from Γ L2 . This can be derived by considering that V max ≈0.38L and the condition of slow self-healing is Γ V max 1, so that combining them we find the same condition.
For fast self-healing values, data depends on the lattice side L, which explains the modification of the scaling law of á ñ I max . However, the slope of the curve is decreasing with respect to the smallest values. This implies that by reducing the delay D sh , the self-healing process is progressively less effective.
In figure 6, we report for a comparison of the maximum current as a function of Γ for various η values and L=30. These results are qualitatively consistent with a fit derived from inserting η eff of equation (2) into equation (1) and calculating the maximum strength. A quantitative good agreement is found for the smaller selfhealing rate and Γλ. Moreover, results for various η collapse on the same curve, at least for slow self-healing values, by plotting ΔI/η, which is captured by equation (2).

Avalanche distributions
When a resistor is burned, the currents are redistributed inside the lattice and rupture avalanches can be triggered. The number s of resistors burned at the same voltage is the size of a single avalanche event. By iterating the simulations, the distribution p(s) of the avalanche sizes can be sampled, representing the probability of having an avalanche of size s. The avalanche distribution displays a decreasing power-law behaviour with a final peak due to the last catastrophic event involving many resistors ( figure 7(a)). If the final avalanche is removed, the distribution can be written in the form of a decreasing power-law with an upper cut-off: where τ is the power-law exponent, and s 0 is the cut-off, which can be written as a function of the lattice size as s 0 =L D , D being the fractal dimension of the avalanche. The exponent τ can be analytically calculated in the FBM, where it assumes exactly the value τ=5/2. In the RFM, it has been evaluated in [41], resulting in τ=2.75 for a diamond lattice and τ=3.05 for a triangular lattice, with fractal dimension of about D=1.18 for both. A robust method to extract these results from data is suggested in [41,48]: once the last catastrophic event is removed, the qth moment of the distribution º á ñ M s q q is evaluated as a function of the lattice size. These quantities have a power-law behaviour~s M L q q , from which an exponent σ q can be evaluated. If the distribution coincides with equation (3), then the following conditions should be verified: if q<τ−1, σ q =0, Figure 5. Rescaled plot of the peak current increase with respect to the standard case for various self-healing delay values. The data collapse on the same power-law for Γ L2 , i.e. for slow self-healing values. The dashed line is the power-law fit on these data, whose result for the exponent is a=1.96 (7). For fast self-healing values, the slope of all curves decreases. Figure 6. Results of the increase of the peak current as a function of Γ obtained for various self-healing rates η. For each series, the dashed line is the fit performed by inserting the expression (2) into equation (1) and calculating its maximum. In the inset, the same plot rescaled by η is reported. while for q>τ−1 σ q =qD+D(1−τ). Thus, from the plot of σ q as a function of q, both τ and D can be evaluated.
Our data are consistent with these results: we obtain D=1.07(4) and τ=3.00 (8). This is close to the value found for a triangular lattice, which is analogous to our case with a non-uniform starting current distribution on the lattice. With these results, the collapse of the avalanche distributions on the same scaling function can be obtained, as shown in the inset of figure 7(a).
In the presence of self-healing, the avalanche distribution is modified, as shown in figure 7(b) for various self-healing delay values. As observed in the previous sections, the peak due to the last catastrophic events increases with Γ, and for fast self-healing values, there are substantial modifications. If we repeat the procedure to extract τ and D, for fast self-healing values, the data rescaled do not collapse on the same curve. On the contrary, for slow self-healing there is a perfect collapse and, for example, for Γ=0.05 we obtain D=0.95 (2) and τ=2.86 (8). In the inset of figure 7(b) we report the comparison of the avalanche distribution once the last events are removed. The slope of the distributions decreases with fast self-healing values. An estimate can be obtained by fitting the first points of the avalanche distribution, where the cut-off function is less influential, and averaging it over the lattice sides. Results are shown in table 1.
We can conclude that in general the exponent of the avalanche distribution decreases with the addition of self-healing, i.e. there is a higher probability of larger rupture avalanches. For slow self-healing values, equation (3) is still valid, while for fast self-healing values there are substantial modifications to the cut-off expression.

Crack roughness
At the end of the simulation, the lattice is split into two regions by a winding crack that is roughly perpendicular to the flow of current, i.e. along the vertical direction. It is possible to study the 'fractality' of the crack in a quantitative way by means of different observables [49], which should give similar results provided that the crack  is indeed self-affine. For comparison between different self-healing rates, we choose the crack width whose scaling law has been investigated in [41]. The first step is to obtain a single-valued function f (y) to describe the crack, where yä[0, L] is the vertical coordinate of the lattice and f (y)ä[0, L] is the coordinate of the crack point along the horizontal direction. In our case, we identify the links belonging to the crack as those connecting points of the two split regions of the lattice. In order to obtain a single value for f (y), from this set we take only the vertical link that has minimum f (y) for each value of y. This is a conventional choice that does not affect the scaling laws.
Hence, we can define the crack width as where the average f is performed over an interval of size lL and á ñ ... is the statistical average over all the samples. The global width º W w L ( ) is expected to scale as a power-law W∼L ζ , where ζ is the fractal exponent of the crack.
In the RFM, an anomalous scaling has been observed [41], i.e. there is another local exponent ζ loc so that for regions of size l=L, the scaling law is~z z z w l l L loc loc ( ) , with ζ=0.80 and ζ loc =0.7 for a diamond lattice. In the case without self-healing, our data are in good agreement with this. In figure 8(a), the crack width is plotted as a function of l for various lattice sizes. From the fit of the global width W we obtain ζ=0.791 (1). The exponent ζ loc can be estimated from a fit on the region l=L of the case L=100 and we obtain ζ loc ; 0.70.
In table 2, results introducing self-healing are reported. The power-law behaviour of the crack width is confirmed except for the largest values of Γ, whose goodness of fit decays in the same range of Lä [12,100]. This can be explained by the larger rupture avalanche at the system failure, so that the crack does not propagate as in the standard case, rather it is formed almost simultaneously. In any case, the ζ exponent increases, since the lattice regions with restored links have smaller rupture probability, so that the final crack propagates around them with a winding path. Thus, on large lattices there is also a larger probability that the winding path amplitudes increase and, in other words, that the crack width increases with L at greater rate.
On the contrary, we do not observe significant variations of the local exponent ζ loc , since the previous effect is less influential on local scales. The invariance of the local scaling properties is confirmed by another  observable, i.e. the total crack length l c , reported in figure 8(b). This can be calculated as the number of links connecting points belonging to the two split parts of the lattice. In the standard case, we find a power lawbehaviour l c =l 1 L δ +l 0 , with δ=1. 16(1). This is another measure of crack tortuosity, since it indicates superlinear length increase of a one-dimensional curve with system size. In the self-healing case, the value of total crack length increases, but the scaling exponent δ does not exhibit substantial modifications. This implies that the tortuosity of the crack on local scales is not affected significantly by self-healing.

Conclusions
We have implemented for the first time a formulation of the RFM that includes self-healing processes and have investigated the modifications occurring for various statistical properties. In general, we must distinguish two regimes of self-healing, indicated for simplicity as 'fast' self-healing regime, if the healing of broken links occurs over a much smaller time scale than the life time of the system, and 'slow' self-healing regime, if the two are comparable.
For slow self-healing values, there is an improvement in the maximum strength and life time, with small modifications to the statistics of damaged links, avalanche distributions and magnitude of the final events leading to the failure of the system. A scaling law for the increase of the peak current with respect to the standard case is found, ΔI∼(Γ L) a with exponent a≈2. The avalanche distribution without final events before failure follows the same scaling law p(s)∼s − τ with a decreasing τ exponent.
For fast self-healing values, the damage of the system at any applied potential is reduced by increasing the value of the healing rapidity Γ, and the characteristic voltage-current curve tends to the ideal elastic case leading to a final catastrophic event, whose average size is increased. This is a crucial aspect, since the system may not appear critically damaged even up to a few steps before the catastrophic failure. The maximum current increases with a reduced slope, so that self-healing becomes progressively less effective. The avalanche distribution, in particular the cut-off expression, is modified and the scaling law for the standard case cannot be used to collapse the data set. Further studies are required to identify a new expression, however results herein show that the power-law exponent has a decreasing trend, i.e. that there is a higher probability of larger avalanches even before system failure.
The self-healing process does not substantially modify the exponent of the scaling laws for the local crack width and the total fracture length. An increase of the scaling exponent of the global crack width is found, which can be ascribed to the presence of regions of restored links not intersected by the crack. For this reason, we expect this to be a universal signature of a self-healing process.
An interesting extension of this study can be to determine reliable precursors to the failure of the system in the presence of fast self-healing, which is of obvious relevance for the practical applications of self-healing materials. Since the effectiveness of self-healing is progressively reduced, our results suggest that a maximum self-healing rate would not necessarily be an optimal solution in all cases.
A biological implication of this study is that self-healing occurring on smaller time scales than the damage rate increases strength and toughness but can lead to abrupt, brittle failures. This could be useful for lizards like salamanders, which use autotomy (i.e. self-amputation) as a defense mechanism to escape predators: in this case, a sharp, rapid fracture process is required when the tail is subjected to a high stress, such as a predator bite [50]. Figure A1. Comparison between square and diamond lattice results for the relative increase of the peak current ΔI/I max (left) and the final number of broken resistors n b /n tot (right) for various self-healing Γ values. Despite small relative variations, the qualitative behaviours are similar.