Current address : ]Arnold-Sommerfeld-Center for Theoretical Physics and Center
for NanoScience, Ludwig-Maximilians-Universität München, D-80333 München,
Germany

# How criticality meets bifurcation in compressive failure of disordered
solids

Ashwij Mayya  Estelle Berthier  [  Laurent Ponson<EMAIL_ADDRESS>Institut Jean Le Rond D’Alembert UMR 7190, Sorbonne Université, CNRS, Paris,
France

###### Abstract

Continuum mechanics describes compressive failure as a standard bifurcation in
the response of a material to an increasing load: damage, which initially
grows uniformly in the material, localizes within a thin band at failure. Yet,
experiments recording the acoustic activity preceding localization evidence
power-law distributed failure precursors of increasing size, suggesting that
compressive failure is a critical phenomenon. We examine here this apparent
contradiction by probing the spatial organization of the damage activity and
its evolution until localization during compression experiments of 2D cellular
solids. The intermittent damage evolution measured in our experiments is
adequately described by a non-stationary depinning equation derived from
damage mechanics and reminiscent of critical phenomena. In this description,
precursors are damage cascades emerging from the interplay between the
material’s disorder and the long-range stress redistributions following
individual damage events. Yet, the divergence of their characteristic size
close to failure, which we observe in our experiments, is not the signature of
a transition towards criticality. Instead, the system remains at a fixed
distance to the critical point at all stages of the damage evolution. The
divergence results from the progressive loss of stability of the material as
it is driven towards localization. Thus, our study shows that compressive
failure is a standard bifurcation for which the material disorder plays a
marginal role. It also shows that precursory activity constitute by-products
of the evolution towards localization and can therefore be used to assess the
residual lifetime of structures.

compressive failure $|$ localization instability $|$ criticality $|$ failure
precursors $|$ failure forecast

††preprint: APS/123-QED

##  I Introduction

Damage localization is the standard mode of failure of materials under
compression. Decoding this degradation process is therefore the cornerstone of
the design of reliable and safe structures such as buildings, bridges, tunnels
and a countless number of mechanical parts under compressive loading
conditions. During their life in service, these structures may progressively
lose their mechanical integrity. Comprehending damage evolution to predict
their remaining lifetime is an essential component of modern tools of
structural design and predictive maintenance. Yet, the appropriate theoretical
concepts for describing damage spreading and ultimately localization are still
vigorously debated and constitute an active topic of research [1, 2, 3, 4].

Continuum damage mechanics is a powerful approach for describing the
compressive failure of materials such as rocks, ceramics or mortar [5, 6, 7,
8]. In this framework, discrete damage mechanisms like microcrack growth are
described at a continuum scale through the degradation of the local elastic
stiffness of the material [9, 10, 11, 12, 13, 14, 15]. Beyond some critical
load level, this softening leads to a bifurcation from the homogeneous damage
field to a localized damage that only grows within a thin band and leading to
material failure [16, 17].

In parallel, and almost independently to the development of damage mechanics,
the intermittent dynamics of damage growth preceding compressive failure has
attracted a lot of attention. Acoustic emissions have been used as a
preferential means of experimental investigation. Experimental measurements
reveal that damage grows through bursts that display robust scale-free
statistics [6, 18, 19, 11, 20, 21, 22, 23, 24, 25, 26]. Accounting for
material disorder, various theoretical works [27, 28, 29, 30, 31] have
proposed to describe failure as a discontinuous (first-order) phase transition
where the precursory damage events emerge from the sweeping of an instability.
These ideas were primarily discussed in the context of toy models of failure,
using e.g. random fuse models. As a result, a direct comparison with the
statistical properties of precursors measured experimentally is not possible,
leaving unresolved the applicability of these concepts to real materials.
Motivated by the observation of an increase of the precursors’ size close to
localization [32, 33, 23, 34, 35, 26], an alternative scenario in which
compressive failure is described as a continuous (second-order) phase
transition was also proposed [36, 37, 38, 39, 33, 40, 41, 26]: similarly to a
large range of driven disordered elastic systems [42, 43, 44], the bursts of
activity characterizing the response of damaging materials are interpreted as
critical fluctuations, or avalanches, that are reminiscent of the so-called
depinning transition, a critical phenomenon emerging from the competition
between disorder and elastic interactions. Above a depinning threshold, damage
is thus expected to grow at some finite speed, eventually leading to failure.
When approaching this critical point, the material should then display scale-
free fluctuations with diverging length and time scales, a feature that has
been observed in some compression experiments [41, 26].

Despite the appeal of such a scenario, it comes in direct contradiction with
continuum damage models which describe compressive failure as a standard
bifurcation taking place in homogeneous solids, for which material disorder
and hence precursors play a minor role. The objective of this work is to
reconcile these two seemingly incompatible scenarios: Is compressive failure a
standard bifurcation or a depinning transition ? Are its precursors the result
of a transition towards a critical point ? These interrogations have crucial
engineering implications:  Do precursory damage events foretell impending
failure? Recently, we examined this issue in a 1D toy model [45]. However, the
generalization of our results to real materials was limited by the short-range
interactions and the system dimension considered in our model. Here, we follow
a different approach: we start from the in-depth characterization of the
damage precursors in a model experimental system. Subsequently, experimental
observations are confronted with the competing scenarios, yielding our two
main results :

  * •

The scale-free statistics of the precursory activity is reminiscent of the
critical avalanche dynamics during depinning of a driven disordered elastic
interface. Yet, in stark contrast with the critical transition scenario, the
divergence of the length and the time scales of precursors close to failure
results from the on-coming localization, a standard bifurcation taking place
in homogeneous solids.

  * •

The approach towards the bifurcation does not drive the damaging solid towards
a critical point. Instead, it remains at some finite distance from criticality
during the whole process of damage accumulation.

Our article is structured as follows. First, we present our experimental setup
(section II). We carry compression test of cohesionless soft cellular solids.
We show that this material behave like an elasto-damageable medium, ensuring
the applicability of our findings to other brittle materials. The
characteristic features of the precursory damage activity are then analyzed at
both the global scale (using the force-displacement response of the specimen)
and the local scale (using full-field measurements of the mechanical
quantities). This multi-scale characterization provides the complete (non-
stationary) statistical structure of precursors that follows robust scaling
laws. To rationalize these experimental observations, we use damage mechanics
that we extend to disordered solids (section III). We derive the evolution
equation of the damage field that reveals the complex connection with driven
disordered elastic interfaces and depinning transition. Our approach also
captures damage localization that is described as a standard instability. In
section IV, we come back to our experimental data and validate several non-
standard aspects of the proposed depinning model. The compatibility of both
competing scenarios with our experimental observations is discussed in section
V. We argue that damage spreading is reminiscent of the non-stationary
depinning dynamics of a driven disordered elastic interface that culminates in
a standard bifurcation at localization. Section VI provides a direct
engineering application of our work. We bring the experimental proof of
concept that precursors can be harnessed for predicting the residual lifetime
of structures. The analysis of the acoustic emissions recorded during our
experiments is presented in section VII. Their statistical similarity with the
acoustic precursors recorded during the failure of standard brittle solids
ensures the generality of our results. The numerical resolution of the
proposed damage model is carried in section VIII. It provides a comprehensive
interpretation of the statistics of failure precursors observed in our
experiments. Finally, the broader implications of the proposed non-stationary
depinning scenario o other phenomena including amorphous plasticity are
discussed in section IX. The methods employed in our experiments and our model
are briefly presented in the Annexes A and B. A thorough description is
provided in the Supplementary Information [46].

##  II Experimental investigation of failure precursors

###  II.1 Damage localization

Taking inspiration from Poirier et al.’s experiments [47], we perform
compression tests of 2D cohesionless soft cellular solids as shown in Fig.
1(a). The results presented below are based on ten different experiments
during which the evolution of the material microstructure and the associated
damage is tracked using full-field measurements. (see Appendix A1 for a
detailed description of our experimental set-up). Following Karimi et al.
[48], the effect of friction between cells is described at a mesoscopic
continuum scale by introducing an equivalent elasto-damageable medium. As
shows in the following sections and justified in detail in Appendix A2, such a
cohesionless cellular solid is a model system that mimics the jerky
dissipative response of brittle disordered solids under compression.
Dissipation taking place at the scale of the individual cells is tracked in
space and time using a high speed camera (see Appendix A4-5 and SI Sec. 1 and
2 for details). We thus circumvent the drawbacks inherent to X-ray tomography
that provides the detailed spatial structure of damage events in 3D materials
but at the cost of temporal resolution [49, 50]. We also overcome the
limitations of acoustic emissions that provide highly resolved time series but
with a rather poor spatial resolution [51, 52].

First, we focus on the average mechanical response of the specimen. A typical
force-displacement curve recorded during a compression test under displacement
control conditions is shown in Fig. 1(b). As evidenced from the snapshots of
the experiments taken at different load levels (top inset of Fig. 1(b) and
Supplementary Video S1), the specimen initially deforms rather uniformly, even
beyond the linear elastic regime. We consider the deviation to cell
circularity beyond the elastic limit as a measure of the damage level,
$d_{\circ}$. As shown in the inset of Fig. 1(c), the larger $d_{\circ}$, the
lower the Young’s modulus $E(d_{\circ})$ of the material, as expected for
elasto-damageable media. Together with its Poisson’s ratio
$\nu(d_{\circ})\simeq\nu_{\circ}$ that remains nearly constant, the function
$E(d_{\circ})$ describes the impact of damage on the mechanical response of
our 2D elasto-damageable solid (see Appendix A2). Damage is observed to grow
homogeneously from the elastic limit $F=F_{\mathrm{el}}$ until the peak load
$F=F_{\mathrm{c}}$ (refer to Table. S1 for a list of notations), except for
the cells close to the boundary where friction with the lateral wall prevails.
Tracking the collapsing of cells from their deviation to circularity
$d_{\circ}>0.5$, we observe that at peak load, no cell has collapsed yet,
except at the boundaries. After peak load, however, a band consisting of
progressively collapsing cells appears, initiating from the top corners of the
structure. The localization band is clearly visible further away from peak
load, as shown on the snapshot in the upper right corner of Fig. 1(b). The
onset loading $\Delta_{\mathrm{c}}$ of localization is inferred from the
evolution of the vertical strain
$\langle\epsilon_{\mathrm{yy}}\rangle_{\mathrm{\mathcal{R}_{2}}}$ averaged
over the bottom region $\mathcal{R}_{2}$ of the specimen (indicated in Fig.
1A): as displayed in the lower inset of Fig. 1(b), for
$\Delta>\Delta_{\mathrm{c}}$,
$\langle\epsilon_{\mathrm{yy}}\rangle_{\mathrm{\mathcal{R}_{2}}}$ saturates
and departs from the strain $\epsilon_{\mathrm{yy}}^{\mathrm{ext}}$ imposed by
the loading machine. On the contrary, the strain
$\langle\epsilon_{\mathrm{yy}}\rangle_{\mathrm{\mathcal{R}_{1}}}$ measured in
the upper region $\mathcal{R}_{1}$ follows the imposed strain. This confirms
that damage localization starts at peak load. To confirm this important
result, we investigate the spatial distribution of damage growth rate near
peak load (see Appendix A5). Before peak load, it is rather homogeneously
distributed, while righter after peak load, its value is several times larger
in the thin band where cells start to progressively collapse.

This behavior is consistent with the observations made in a wide range of
brittle solids [53, 54, 51, 55, 56]. It is also captured by the non-local
damage models recently proposed in Berthier et al. [57] and Dansereau et al.
[58] (see SI Sec. S3F [46] for the analytical prediction). If the experiments
were under force control conditions, as in most real-life structural
applications, a sudden collapse of the cells resulting in the catastrophic
failure of the specimen would also occur at peak load.

###  II.2 Precursors as cascades of damage events

Figure 1: (a) Schematic of the compression experiment depicting the front view
of the hexagonal packing of soft cellular solids. The box dimensions are
$205\mathrm{mm}\times 170\mathrm{mm}\times 30\mathrm{mm}$. Owing to the
displacement $\Delta$ applied to the specimen through a piston moving at a
rate $v_{\mathrm{ext}}$, the cells undergo a deformation that is recorded
using a high speed camera. The following figures describe the typical
mechanical response of the specimens, as observed in one of our experiments.
(b) Force-displacement response of the specimen. Top insets: The emergence of
a localization band of collapsed cells, corresponding to highly localized
deformations, is visible on the snapshots of the specimen taken at different
load levels. Bottom inset : The strain level averaged over the bottom region
$\mathcal{R}_{2}$ indicated in panel (a) saturates after localization for
$\Delta>\Delta_{\mathrm{c}}$ while the one measured in the upper region
$\mathcal{R}_{1}$ follows the linear trend
$\epsilon_{\mathrm{yy}}^{\mathrm{ext}}=\Delta/H=v_{\mathrm{ext}}\,t/L$ imposed
by the loading machine. (c) Construction of an equivalent force control
experiment from the mechanical response of the specimen measured under
displacement control. The start and the end of a damage precursor taking place
at a constant force $F_{0}$ are denoted by $\Delta_{\mathrm{ini}}$ and
$\Delta_{\mathrm{end}}$. The precursor size $S_{\mathrm{global}}=\Delta
E_{\mathrm{d}}$ defined as the dissipated energy during the event corresponds
to half the work of the external force $\Delta
W=F_{0}\,(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})$. Inset : The stiffness
degradation of the effective medium with increasing damage. (d) Variations of
the precursor size $S$, in terms of dissipated energy, with the distance to
failure $\delta$. Inset : Comparison of precursor sizes computed from the
global analysis and the local analysis. Figure 2: (a) Energy density maps,
$\rho(\vec{x})$, and their corresponding binary thresholded formats (inset) of
energy dissipation at various distances to failure in a typical experiment.
The common color bar depicting the scale for the dissipation energy density
$\rho$ ($\mathrm{J/m^{2}}$) is presented under the first panel. The spatial
extent of the cascades extracted using the auto-correlation of the thresholded
maps are denoted by the red bracket in the inset. (b) Evolution of clusters
within a cascade (second panel in (a), at $\delta\simeq 0.025$) with colors
representing unique cluster IDs. The center of mass of each cluster is shown
in the background with a marker of same color as the cluster. Here, the third
axis is time. The largest clusters appear in multiple slices of time. Scaling
of the size (in terms of dissipated energy) (c) and duration (d) of cascades,
normalized by their maximum values, as a function of their characteristic
spatial extent $\xi$, normalized by the system size $L$. Inset in (c) : Pair
correlation function of the centers of mass of the clusters within a cascade
permitting to extract the fractal dimension $d_{\mathrm{f}}\simeq 1.15$. (e)
Distribution of the cascade sizes obtained during the whole duration of
experiments (diamonds) and in the vicinity of final failure (circles). (f)
Distribution of waiting time between cascades. The waiting time is defined as
the difference in time-stamps of the arrival of two successive avalanches
$\tau_{w}\mathrm{(s)}=t_{i}(\Delta_{ini})-t_{i-1}(\Delta_{ini})$. Inset:
distribution of waiting times considered as the difference in value of force,
$\tau_{\mathrm{F}}\mathrm{(N)}=F_{i}-F_{i-1}$. (g) Fractal analysis of the
spatial distribution of the seeds of the precursors in a typical experiment
(inset), showing a behavior reminiscent of a spatial Poisson process
($d_{\mathrm{f}}\simeq 2$). The statistical analyses is based on data recorded
during ten experiments.

We now analyze the precursory damage activity taking place before peak load. A
closer examination of the force-displacement curve in Fig. 1(c) reveals sudden
force drops of amplitude much larger than the precision $\pm
0.05~{}\mathrm{N}$ of our load cell. These drops are followed by a linear
increase of the force, recovering the force drop albeit with a degraded
macroscopic stiffness (lower slope). The investigation of the specimen
response under different loading rates shows similar behaviors. This
alternating sequence of damage growth and elastic re-loading is reminiscent of
the avalanche dynamics observed in driven disordered elastic systems [42, 43,
44]. We construct the mechanical response of the specimen (in red) in an
equivalent force control experiment where displacement jumps (from
$\Delta_{\mathrm{ini}}$ to $\Delta_{\mathrm{end}}$) at constant force
correspond to a cascade of damage growth, also called an avalanche (Fig.
1(c)). The precursors defined this way can be shown to be statistically
similar to those that would be measured during an actual force control
experiment (see SI, Sec. S4 for a numerical validation). As under force
control conditions damage cascades take place at constant force, the work of
the loading machine $\Delta W$ during the event can be shown to contribute
equally to the increase $\Delta E_{\mathrm{el}}$ of the elastic energy and to
the dissipation $\Delta E_{\mathrm{d}}$ by damage (see Appendix A3). Hence,
$\Delta E_{\mathrm{d}}=\Delta
W/2=(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})\,F_{0}/2$ and we define this
quantity as the precursor size $S_{\mathrm{global}}$. The evolution of
$S_{\mathrm{global}}$ with the distance to failure
$\delta=(F_{\mathrm{c}}-F)/(F_{\mathrm{c}}-F_{\mathrm{el}})$ is shown in Fig.
1(d). The introduction of $F_{\mathrm{el}}$ ensures that $\delta=1$
corresponds to the beginning of the damage accumulation regime. However,
choosing another definition, as e.g., $(F_{\mathrm{c}}-F)/F_{\mathrm{c}}$,
does not modify our conclusions. We note that the cascading dynamics and its
amplification on approaching failure $(\delta\rightarrow 0)$ observed in our
2D cellular material under compression are reminiscent of the intermittent
damage activity evidenced by acoustic emissions in standard brittle materials
(see Fig. 1(d)).

Alternatively, we can also identify and characterize the precursors at the
local scale using our (time-resolved) full-field measurement of the
displacement and the damage field (see SI Sec. 1C and Figs. S1(c)-(f) [46] for
details on the local analysis). From these quantities, we compute the field of
stored elastic energy in the effective elasto-damageable medium. Considering
energy balance at the local scale, we can thus determine the dissipation
energy density $\rho(\vec{x},\Delta)$ that we integrate over an avalanche
$\rho(\vec{x})=\int_{\Delta_{\mathrm{ini}}}^{\Delta_{\mathrm{end}}}\rho(\vec{x},\Delta)\,d\Delta$.

Maps of dissipation energy density $\rho(\vec{x})$ depicting the complex
spatial structure of precursors are presented in Fig. 2(a). We observe a
diffuse pattern, yet containing locally well defined regions of varying
intensity and size. These clusters are reminiscent of the time and space
correlated structure of incremental damage events. In practice, highly
correlated individual damage events can be grouped together by implementing a
spatio-temporal clustering algorithm on the fields $\rho(\vec{x},\Delta)$
recorded during $\Delta_{\mathrm{ini}}<\Delta<\Delta_{\mathrm{end}}$. This
reveals the cluster-like structure of precursors, illustrated in Fig. 2(b)
(see also Video S2) for the precursor shown in the second panel of Fig. 2(a).

The energy $S_{\mathrm{local}}=\int\rho(\vec{x})d\vec{x}$ dissipated during
the cascade compares well with the precursor size $S_{\mathrm{global}}$
inferred from the force-displacement response, see Fig. 1(d). It is also in
good agreement with the precursor size computed using the field of dissipated
energy inferred from the incremental damage field, thus validating the
assumption of local energy balance (see SI Sec. S2, Figs. S2(d) and S2(e)
[46]). Finally, this agreement comes in support to the description of our
cohesionless cellular solid as an elasto-damageable medium.

###  II.3 Statistical characterization of precursors

We now explore the properties of the damage cascades observed in our
experiments. First, their spatial extent is determined from (thresholded) maps
of dissipation energy density (insets of Fig. 2(a)). The employed threshold
value $\rho^{\ast}$ is inferred from the distribution of local dissipation
densities that follows an exponential decay
$P(\rho)\propto\mathrm{e}^{-\rho/\rho^{\ast}}$ (see SI Sec. S1(c) and Fig.
S1(g) [46]). We extract the characteristic length $\xi$ from the 2D auto-
correlation of thresholded dissipation map (SI Sec. S1(d) and Fig. S1(h)
[46]). The length grows with the avalanche size as
$S\propto\xi^{d_{\mathrm{f}}}$, where $d_{\mathrm{f}}\simeq 1.07$ is the
fractal dimension, see Fig. 2(c). While cascades spread over the whole
specimen (see Fig. 2(a)), $\xi$ represents the spatial extent of the largest
clusters constituting the cascade. Interestingly, $\xi$ reaches the specimen
size $L$ on approaching failure, which implies an upper limit on the size of
the precursors. An independent estimate of the fractal dimension of the
precursors is obtained from the spatial distributions of the clusters. We
identify the location of their center of mass, as illustrated on the right-end
of Fig. 2(b), and then compute the correlation function $C(r)\propto
r^{d_{\mathrm{f}}}$ defined as the fraction of pairs of points whose
separation is less than $r$ [59]. This provides a fractal dimension
$d_{\mathrm{f}}\simeq 1.15$ (see inset of Fig. 2(c)) compatible with the one
obtained from the spatial distribution of the individual damage events. Note
that the spatial extent $\xi$ of a damage cascade is different from its size
$S$ that corresponds to the energy dissipated by damage during the cascade.

We now seek to determine the characteristic duration of damage cascades. We
come back to the force-displacement response and explore the sequence of load
drops observed within an avalanche. The precursor duration $T$ is defined as
the number of load drops (also see SI Sec. S1B and Figs. S1(a) and S1(b)
[46]). It scales with the characteristic length $\xi$ of the precursor as
$T\propto\xi^{z}$, see Fig. 2(d), with dynamic exponent $z\simeq 0.53$.

Thus, a damage cascade is characterized by its size, its spatial extent and
its duration. All three quantities are related to each other by scaling laws.
The probability distribution of these quantities is studied in Fig. 2(e),
where we focus on the distribution $P(S)$ of precursor sizes, the other
distributions $P(\xi)$ and $P(T)$ being inferred from the previous scaling
laws. Considering all the precursors ($\delta\in[0,1]$) or only the ones close
to localization $(\delta\in[0,0.01])$, both distributions follow a power-law
statistics but with two different exponents $\beta_{\mathrm{tot}}\simeq 2.34$
and $\beta\simeq 1.30$, respectively. This difference results from the
increase of the size of the largest precursors on approaching failure (see SI
Sec. S4A [46, 60]). The exponent $\beta$ will be connected later with the
marginal stability of the material elements [61, 62] that will be investigated
in Sec. IVA.

Finally, we characterize the correlations in the sequence of damage cascades.
The distribution $P(\tau_{\mathrm{w}})$ of waiting times separating two
successive damage cascades is shown in Fig. 2(f). It follows an exponential
law
$P(\tau_{\mathrm{w}})\propto\mathrm{e}^{-\tau_{\mathrm{w}}/\tau_{\mathrm{w}}^{\star}}$,
defining a characteristic waiting time $\tau_{\mathrm{w}}^{\star}$. This
result is at odd with the power-law distribution of waiting times separating
acoustic events in compression experiments [21, 22, 23]. To confirm our
observation, we perform additional experiments with a loading rate
$v_{\mathrm{ext}}$ twice smaller. Interestingly, we also measure an
exponential distribution, but with a characteristic waiting time about twice
larger (see Fig. 2(f)). We thus replace the waiting time $\tau_{\mathrm{w}}$
by the force increment $\tau_{\mathrm{F}}$ separating two successive
precursors, so that distributions corresponding to different loading rates
collapse on a single curve, as shown on the inset of Fig. 2(f). The
exponential distribution of waiting time, characteristic of uncorrelated
events described by a Poisson process, suggests that precursors are triggered
independently from each others. This is further confirmed by the spatial
distribution of the seeds (first damage event) of precursors that we define as
the center of mass of the cluster appearing at $t(\Delta_{ini})$. The fractal
analysis of the precursor seeds provides $C(r)\propto r^{2}$ a behavior
reminiscent of spatially uncorrelated events (see Fig. 2(g)).

##  III Theoretical modeling of compressive failure

The statistical features of the precursors measured in our experiments
strikingly remind the avalanche dynamics of elastic interfaces driven in
disordered media. In these models, an elastic interface responds to a
continuously increasing drive and exhibits scale free avalanches or crackling
noise [42, 43, 44]. The size, spatial extent and duration of the avalanches
are related by scaling laws with universal exponents that depend on the
interface elasticity and its dimension. For interfaces with long-range
elasticity, avalanches are formed by a set of correlated clusters that are
spatially disconnected similar to damage clusters within a cascade observed in
our experiments (see Fig. 2(b)). Taking inspiration from Weiss et al. [40] and
using the non-local theory proposed in Dansereau et al. [58], we derive below
an evolution equation of the damage field in the specimen for a compression
test under force control conditions. This theoretical formulation sheds light
on the connection with models of driven disordered elastic interfaces. We
provide here the main ingredients of the derivation of the damage evolution
law, the detailed calculations being presented in SI Sec. S3A-D [46]

First, we assume that the material behaves as an elasto-damageable solid. We
thus introduce a damage field $d(\vec{x},t)$ that describes the level of
damage accumulated in the specimen at the location $\vec{x}$ and time $t$.
Damage growth is inferred from a balance of energy, by comparing two
quantities: the local driving force $Y[d(\vec{x},t),t]$ which provides the
rate of elastic energy released for an incremental growth of damage, and the
damage resistance $Y_{\mathrm{c}}[d(\vec{x},t)]$ which provides the material
resistance to damage and corresponds to the rate of energy dissipated for an
incremental growth of damage. The first quantity is similar to the elastic
energy release rate introduced in fracture mechanics, which drives crack
propagation and is a quadratic function $Y\propto\sigma_{0}^{2}$ of the
nominal compressive stress $\sigma_{\circ}=F/(L\,W)$ applied by the test
machine [63, 64]. The damage resistance is equivalent to the fracture energy
introduced in fracture mechanics. Note that in our model, $Y_{\mathrm{c}}$
depends not only on $\vec{x}$ as precursors emerge from the material
inhomogeneities, but also on $d$ as a damage event in $\vec{x}$ may change the
subsequent failure resistance in the same material element. The damage then
increases in the material element $\vec{x}$ if the local value of the driving
force $Y(\vec{x})$ reaches the material resistance $Y_{\mathrm{c}}(\vec{x})$.

To describe the damage field fluctuations resulting from the material
heterogeneities, we introduce a reference damage level $d_{\circ}=\langle
d(\vec{x},t_{\circ})\rangle$ and the damage field perturbations $\Delta
d(\vec{x},t)=d(\vec{x},t)-d_{\circ}$ over the time $\delta t=t-t_{\circ}\ll
t_{\circ}$ to ensure that $\langle\Delta
d(\vec{x},t)\rangle_{\mathrm{\vec{x}}}\ll d_{\circ}$. The driving force and
the damage resistance can then be linearized as
$Y[d(\vec{x},t),\sigma_{\circ}]=Y_{\circ}(d_{\circ},\sigma_{\circ})+\Delta
Y[\Delta d(\vec{x},t),\sigma_{\circ}]$ and
$Y_{\mathrm{c}}[d(\vec{x},t)]=Y_{\mathrm{c\circ}}(d_{\circ})+\Delta
Y_{\mathrm{c}}[\Delta d(\vec{x},t)]$. The zero-order equation
$Y_{\circ}(d_{\circ},\sigma_{\circ})=Y_{\mathrm{c\circ}}(d_{\circ})$ provides
the relationship between the reference damage level $d_{\circ}$ and the
reference applied load $\sigma_{\circ}=\sigma(t_{\circ})$. In the following,
we investigate how the damage field perturbations
$\Delta\dot{d}(\vec{x},t)\propto\Delta Y[\Delta
d(\vec{x},t),\sigma_{\circ}]-\Delta Y_{\mathrm{c}}[\Delta d(\vec{x},t)]$
evolves over time. We write the total driving force as the sum of three terms
[40, 58],

$\begin{split}\Delta\dot{d}(\vec{x},t)\propto&\,\mathcal{K}(\sigma_{\circ})\left[v_{\mathrm{m}}(\sigma_{\circ})\,t-\Delta
d(\vec{x},t)\right]+\\\ &\psi(\sigma_{\circ})\ast[\Delta
d(\vec{x},t)-\langle\Delta
d\rangle_{\mathrm{\vec{x}}}]-y_{\mathrm{c}}[\vec{x},d(\vec{x},t)].\end{split}$
(1)

The first (local) term comprises the effect of the driving by the test
machine, where the driving speed $v_{\mathrm{m}}\propto v_{\mathrm{ext}}$ sets
the damage growth rate. Considering a pseudo-interface of position $\Delta
d(\vec{x},t)$, this term acts as a rigid plate moving at a speed
$v_{\mathrm{m}}$ and pulling on the interface with springs of stiffness
$\mathcal{K}$ (see SI Fig. S3 [46] for a schematic representation). The second
term is non-local. It describes the interactions within the specimen, its
values in $\vec{x}$ depends on the damage level $\Delta d(\vec{x},t)$
everywhere in the specimen. In practice, the kernel $\psi(\sigma_{\circ})$
(provided in Eq. (7)) describes the spatial structure of the redistribution of
driving force taking place in the aftermath of an individual damage event. It
decays as $\psi\propto 1/r^{2}$ [40, 58]. Its angular dependence is shown in
Fig. 3(a) for the particular case of uni-axial compression. It exhibits a
quadrupolar symmetry with non-positive regions (in blue). Hence only a
fraction of the neighboring elements are reloaded in the aftermath of a damage
event while the others (located above and below the damaged element) are
actually unloaded. The third term represents the effect of material disorder.
Its spatial average is close to zero as the contribution of the hardening
$\langle\Delta Y_{\mathrm{c}}\rangle_{\mathrm{\vec{x}}}=\eta\,\langle\Delta
d\rangle_{\mathrm{\vec{x}}}$ (where $\eta$ is a hardening parameter and
observed in our experiments) is taken into account in the first term (see SI
Sec. S2 and Figs. S2(b) and S2(c) [46]). The presence of $d(\vec{x},t)$ as an
argument of the disorder term implies that the evolution equation (1) is
strongly non-linear, leading to the rich phenomenology that we now discuss.

Equation (1) provides a clear connection between damage evolution and
disordered elastic interfaces: the accumulated damage field is analogous to a
2D elastic interface $\Delta d(\vec{x},t)$ driven at the speed
$v_{\mathrm{m}}$ through a 3D disordered medium (see SI Sec. S3E, Fig. S3
[46]). As a result, damage is expected to grow through bursts characterized by
scaling laws involving critical exponents reminiscent of the so-called
depinning transition. As a first test of our model, we compare the
theoretically predicted exponents with the one measured experimentally.
Investigating the avalanches dynamics of 2D interfaces with non-positive
interactions in the context of amorphous plasticity, Lin et al. [61, 65]
predicted the exponent values $\beta=1.51$, $d_{\mathrm{f}}=1.10$ and $z=0.57$
that agree reasonably well with the ones measured in our experiments.
Corrections to these predictions from the numerical solution of the evolution
equation (1) are provided at the end of our manuscript, in section VIII. They
improve further the agreement with the experimentally measured exponents.

Despite the ability of this approach to describe the scaling behavior of
precursors, we note that two important features of the damage evolution
equation (1) differ from standard models of driven elastic interfaces. First,
the long-range elastic interactions result in both reloading and unloading of
material elements in the aftermath of a damage event. Second, a subtler but
more important feature is that the evolution equation (1) describes a non-
stationary depinning scenario culminating in a bifurcation at the localization
threshold. In the following, we examine these aspects in details and discuss
their implications on the damage accumulation process preceding failure.

##  IV Atypical aspects of damage evolution as a driven disordered elastic
interface

###  IV.1 Non-positive elastic interactions

We describe below a methodology that we employ to characterize the elastic
interactions driving the cooperative dynamics of damage. This method takes
advantage of our full field characterization of damage precursors. Indeed, it
turns out that the spatial distribution of events within an avalanche encrypts
the range and the anisotropy of the elastic interactions, as we will see
below. The first step consists in computing the incremental damage field
$\delta d(\vec{r})$ during an avalanche. Its 2D auto-correlation function
$C(\vec{\delta r})=\langle\delta d(\vec{r})\cdot\delta d(\vec{r}+\vec{\delta
r})\rangle_{\mathrm{\vec{r}}}$ averaged over several avalanches is presented
in Fig. 3(b). The angular distribution (at a fixed distance) shows a clear
quadrupolar symmetry similar to the re-distribution pattern of the theoretical
interaction kernel (Fig. 3(a)).

Figure 3: (a) Angular distribution of the long-range interaction kernel
derived theoretically for the case of uni-axial compression. (b) 2D auto-
correlation map of the incremental damage field of precursors and angular
distribution of the correlations at a fixed distance as obtained from the 50
largest avalanches measured in one of our experiments and our numerical
simulations. (c) Variations of the correlation function of the incremental
damage field with distance along the horizontal $x$-axis and comparison with
the scaling $\psi\propto 1/\delta r^{2}$ of the theoretical interaction kernel
(dashed line). (d) Variations of the depth $\Delta d^{\ast}$ of the damage
cascade with the cascade size $S$. (e) Distribution of the local distance to
failure $\delta Y(\vec{x})$ close and far from failure.

Remarkably, the correlation along the horizontal axis where the reloading is
maximal decays as $C(\delta r)\propto 1/\delta r^{2}$ (see Fig. 3(c)), a
behavior also in line with the theoretical predictions of the elastic kernel,
Eq. (7). These observations support further the applicability of the non-local
damage mechanics to our experiments.

The presence of an unloading region in the interaction kernel has several
important implications. First at the scale of a material element, the approach
to failure is non-monotonic as the local driving force for damage can both
increase and decrease over time. As a result, the probability that an elements
damages more than an once during a cascade is low. In practice, the depth of
avalanches $\Delta d^{\ast}$ defined as the average damage increment of the
elements involved in the cascade is constant and does not vary with the
precursor size $S$, as shown in Fig. 3(d) (see also SI Sec. S4A [46]). This is
at odds with the behavior of driven elastic interfaces with positive
interactions for which the depth of avalanches scales with their size. Another
crucial, yet more subtle difference relates to the distribution of net driving
force $\delta Y(\vec{x})=Y_{\mathrm{c}}(\vec{x})-Y(\vec{x})$ that controls the
(marginal) stability of the specimen. $\delta Y(\vec{x})>0$ provides the
increment of driving force required for triggering damage. Its distribution,
computed over all the material elements, is expected to scale as $P(\delta
Y)\propto\delta Y^{\theta}$ [61, 65, 62, 66, 67]. Positive interactions lead
to $\theta=0$, pointing out the presence of a finite number of material
elements close to failure. On the contrary, the number of elements close to
failure vanishes for sign-changing interactions, leading to $\theta>0$. The
experimental determination of the exponent $\theta$ is quite challenging, as
it requires a priori the knowledge of the material disorder. In practice, we
circumvent this difficulty by computing the driving force $Y(\vec{x},t)$
(according to SI Eq.S4 [46]) at each time step for each material element and
determine $Y_{\mathrm{c}}[d(\vec{x},t)]$ retrospectively from the value of
$Y(\vec{x},t)$ when the material element damages (see Appendix A4, SI Sec. 2
and Fig. S2(a) [46] for details on the method). Figure 3(e) shows the
distribution $P(\delta Y/\langle\delta Y\rangle)$ close and far from
localization. In both cases, we measure a positive exponent $\theta>0$, a
particularly non-trivial property that comes in support of the proposed model.
Interestingly, $\theta$ increases as the specimen approaches failure, a
feature that possibly arises from the non-stationary nature of the evolution
equation (1). The exponent $\theta$ has been shown to be related to
$\beta=2-\frac{\theta}{\theta+1}\frac{d}{d_{\mathrm{f}}}$ [61] that
characterizes the distribution of precursor size $S$. This scaling relation
predicts a slight decrease of $\beta$ on approaching failure, a subtle effect
that is consistent with our experimental observations. A similar trend of an
increasing $\theta$ has been reported in direct simulations of sheared
amorphous solids that are also characterized by sign-changing interactions
[62, 68].

###  IV.2 Divergence of precursors

We now come back on the observations made in Fig. 1(d) and Fig. 2(a) of an
increase of the size and the spatial extent of precursors close to failure. As
shown in Fig. 4(a), the average precursor size increases as a power-law with
the distance to failure, $\langle S\rangle\sim 1/\delta^{-\alpha}$ where
$\alpha\simeq 0.57$. Following the scaling relations
$S\propto\xi^{d_{\mathrm{f}}}$ and $T\propto\xi^{z}$, the associated length
and time scales then also diverge on approaching failure. This is confirmed by
the variations of the precursor spatial extent directly measured from our
local analysis (see SI Sec.1D [46]) in the inset of Fig. 4(a).

Figure 4: Variation with distance to failure of (a) the average precursor size
$\langle S\rangle$, the average spatial extent of precursors
$\langle\xi\rangle$ in the inset, (b) the activity rate $dN_{\mathrm{S}}/dt$,
(c) the average size of load drop events $\langle A\rangle$ and (d) the event
activity rate $dN_{\mathrm{A}}/dt$. While $S$ and $dN_{\mathrm{S}}/dt$
characterize the intermittent damage activity under force controlled
conditions, $A$ and $dN_{\mathrm{A}}/dt$ are relevant for displacement imposed
conditions. The product of both quantities $dE_{d}/dt=\langle
S\rangle\,dN_{\mathrm{S}}/dt=\langle A\rangle\,dN_{\mathrm{A}}/dt$ provides
the dissipation rate that also diverges as $dE_{d}/dt\sim\delta^{-\alpha}$
with $\alpha=1/2$.

Notably, the activity rate $dN_{S}/dt$ i.e., the number of cascades per
interval of time is rather constant during damage accumulation, see Fig. 4(b).
This is in line with our previous observation of an exponential distribution
of waiting times, supporting further that precursors emerge from a (random)
Poisson process. As the dissipation rate $dE_{\mathrm{d}}/dt$ during the
intermittent damage evolution writes as the product of the average precursor
size with the precursor rate, $dE_{\mathrm{d}}/dt=\langle
S\rangle\,dN_{\mathrm{S}}/dt$, we obtain $dE_{\mathrm{d}}/dt\propto
1/\delta^{-\alpha}$. As argued in the next section, the divergence of the
dissipation rate on approaching peak load is reminiscent of damage
localization, a feature that results from the loss of stability of the
specimen. To further test this idea, we reanalyze our data considering the
actual displacement imposed conditions that also gives rise to damage
localization at peak load. The dissipation rate then writes as the product of
the average load drop size $\langle A\rangle$ with their activity rate
$dN_{A}/dt$ (refer to SI Sec.1B and Fig.S1(b) [46] for the definition of
precursors under displacement imposed loading conditions). We observe that
both load drops $\langle A\rangle\propto\delta^{-\alpha_{\mathrm{A}}}$ and
precursor rate $dN_{A}/dt\propto\delta^{-\alpha_{N_{A}}}$ diverge on
approaching localization, see Fig. 4(c) and (d). The exponent
$\alpha=\alpha_{A}+\alpha_{N_{A}}\simeq 0.44$ characterizing the divergence of
the dissipation rate under displacement imposed conditions is close to $1/2$,
as expected and accounts for the numerical observations of Girard et al. [33]
who reported $dE_{\mathrm{d}}/dt\propto\delta^{-0.4}$.

##  V Damage localization : depinning transition or standard bifurcation ?

The divergence of the precursor size, and hence of the characteristic
timescale and length scale of precursors on approaching failure supports a
priori the interpretation of compressive failure as a critical point [40, 26].
However, the comparison between the scaling exponents predicted by the
critical point scenario and the exponents measured in our experiments tells a
different story (see Table 1).

Table 1: Comparison between the exponents predicted by the critical point scenario that interprets failure as a depinning transition and the exponent measured in our experiments. The three exponents provided below describes the divergence $\sim\delta^{-\gamma}$ of the avalanche size $S$, its spatial extent $\xi$ and its duration $T$ on approaching failure. The theoretical exponents describe the behavior of a driven 2D interface with long-range non-positive elastic interactions approaching the depinning transition [61, 26]. | Critical point scenario | Experiments
---|---|---
| $\gamma$ | prediction - 2D | $\gamma$ |
$S$ | $\nu d_{\mathrm{f}}$ | 1.27 | $\alpha$ | 0.57
$\xi$ | $\nu$ | 1.16 | $\alpha/d_{\mathrm{f}}$ | 0.53
$T$ | $\nu z$ | 0.66 | $z\alpha/d_{\mathrm{f}}$ | 0.28

To explain the difference between the theoretical predictions and our
experimental observations, we further develop the model proposed in Sec. III.
We focus on the non-stationary aspects of the evolution equation (1), namely
the stiffness $\mathcal{K}$ and the driving speed $v_{\mathrm{m}}$ that are
given by the following expressions:

$\begin{split}\mathcal{K}(\sigma_{\circ})&=\left.\frac{\partial(Y_{\mathrm{c\circ}}-Y_{\circ})}{\partial
d_{\circ}}\right|_{\sigma_{\circ}}\\\
v_{\mathrm{m}}(\sigma_{\circ})&=v_{\mathrm{ext}}\left(\frac{\left.\partial
Y_{\circ}/\partial\sigma_{\circ}\right|_{d_{\circ}}}{\mathcal{K}(\sigma_{\circ})}\right).\end{split}$
(2)

$\mathcal{K}$ controls the stability of the damage evolution. Indeed, a
negative value of $\mathcal{K}$ implies that the net driving force
$Y-Y_{\mathrm{c}}$ increases with the damage level, leading to its unstable
growth and thus failure. It turns out that $\mathcal{K}$ goes to zero on
approaching peak load (see SI Sec. S3F [46]), in line with the stability
condition under force controlled conditions. The driving speed
$v_{\mathrm{m}}$ is inversely proportional to $\mathcal{K}$ and hence goes to
infinity. A linear expansion of the damage evolution equation close to peak
load $\sigma_{\circ}<\sigma_{\mathrm{c}}$ (see SI Sec. S3G [46]), provides the
asymptotic behavior of the damage growth rate $\dot{\Delta d}\sim
1/\sqrt{\delta}$. As $\Delta\dot{d}\propto v_{\mathrm{m}}$, owing to (2) we
obtain $v_{\mathrm{m}}\sim 1/\sqrt{\delta}$ and thus
$\mathcal{K}\sim\sqrt{\delta}$.

What are then the consequences of the divergence of the speed of the pseudo-
interface at peak load ? As the rate of dissipated energy is controlled by the
damage growth rate, one expects $dE_{\mathrm{d}}/dt\propto 1/\sqrt{\delta}$.
Considering the intermittency of damage evolution, the dissipation rate writes
as the product of the precursors’ size with the precursors’ rate,
$dE_{\mathrm{d}}/dt=\langle S\rangle\,dN_{\mathrm{S}}/dt$. As
$dN_{\mathrm{S}}/dt$ remains constant during the experiment (Fig. 4(b)), a
feature expected for disordered elastic interfaces, it follows that $\langle
S\rangle\sim 1/\sqrt{\delta}$, a prediction that accounts for our experimental
observations.

The implication of our model is clear: the divergence of the size of the
precursors close to failure, and hence the divergence of the length scale and
time scale of the fluctuations, result from the presence of a standard
bifurcation at peak load. The progressive loss of stability of the specimen on
approaching peak load is accompanied by a divergence of the damage growth rate
(and thus a divergence of the precursor size), a behavior that has nothing to
do with the presence of disorder. This mechanism is in stark contrast with the
divergence of fluctuations near a critical point, a feature that vanishes if
the disorder is shut down.

Figure 5: Interpretation of the two phenomena underlying compressive failure
of disordered solids. The left panel describes the behavior of the specimen as
it is driven to failure. The amplitude of the fluctuations of the damage field
(in red) and the dissipation rate (yellow) increase near failure ($\delta=0$).
The spatial variations $\langle(d(\vec{x})-\langle
d(\vec{x})\rangle_{\vec{\mathrm{x}}})^{2}\rangle_{\vec{\mathrm{x}}}$ of the
accumulated damage are essentially zero before peak load $(\delta<0)$, but
rapidly increase after localization $(\delta>0)$. The state of the specimen
bifurcates from a homogeneously growing damage to a localized damage growth.
This bifurcation is preceded by an increase
$v_{\mathrm{m}}\sim\delta^{-\alpha}$ of the dissipation rate, and thus an
increase of $\langle S\rangle\sim\delta^{-\alpha}$ of the fluctuation
amplitude. The right panel describes the process of damage growth at a fixed
distance to failure. By renormalizing the non-stationary features of the
evolution equation (see Eq.(3)), we realize that the damage field behaves as
an elastic interface driven in a disordered medium at constant speed.
$v_{\mathrm{m}0}\simeq v_{\mathrm{m}}\delta^{\alpha}$ at a fixed distance
above the critical point (here, the depinning transition) even when we reach
peak load ($\delta=0$) and thus bifurcation.

To further support this claim, we rewrite the damage evolution law (1) using
the new variable
$\underset{\widetilde{}}{d}(\vec{x})=d(\vec{x})\,\sqrt{\delta}$. The obtained
expression ensures a straightforward connection with standard (stationary)
depinning models [42, 44]:

$\begin{split}\underset{\widetilde{}}{\Delta\dot{d}}&\propto\mathcal{K}_{0}\left[v_{\mathrm{m}0}\,t-\Delta\underset{\widetilde{}}{d}(\vec{x},t)\right]+\\\
&~{}~{}~{}~{}{\psi}(d_{\circ}\sqrt{\delta})/\sqrt{\delta}\ast[\Delta\underset{\widetilde{}}{d}(\vec{x},t)-\langle\Delta\underset{\widetilde{}}{d}\rangle_{\mathrm{\vec{x}}}]-y_{c}[\vec{x},\underset{\widetilde{}}{d}(\vec{x},t)].\end{split}$
(3)

Under this form, the evolution equation displays both a constant stiffness
$\mathcal{K}_{0}=\mathcal{K}/\sqrt{\delta}$ and a constant driving speed
$v_{\mathrm{m}0}=v_{\mathrm{m}}\sqrt{\delta}$. The normalized damage field
$\underset{\widetilde{}}{d}$ behaves as an elastic interface driven at
constant finite speed during the whole damage accumulation regime. As a
result, the specimen remains at a fixed distance to the depinning transition
that corresponds to the limit $v_{\mathrm{m0}}\rightarrow 0$ and
$\mathcal{K}_{0}\rightarrow 0$. Our conclusions are in stark contrast with a
depinning scenario where the divergence of the precursory activity results
from the interpretation of compressive failure as a critical transition and
the evolution of the specimen towards this critical point. In other words, the
increasing applied load drives the specimen towards instability without
driving it towards more criticality. The intricate connection between damage
precursors and localization is illustrated in Fig. 5.

We now would like to highlight the strategic value of our findings to
structural health monitoring. In our theoretical description of compressive
failure, the evolution of precursors is described by robust scaling laws that
are independent of the material properties and the loading conditions. In
particular, as they apply for both force and displacement driven experiments,
they serve as early warning signals of impending failure.

##  VI Failure prediction from precursory activity

Figure 6: (a) Time series data of failure precursors obtained during the
compression test and data available for prediction (red) at
$t_{\mathrm{cur}}$. (b) Variation of the predicted remaining lifetime at
different instances $t_{\mathrm{cur}}$ represented as a fraction of the time
to failure. The error bars provide intervals with 90$\%$ confidence levels.
Inset : Error on the predicted remaining lifetime.

We now harness the scaling behavior of the precursors and bring an
experimental proof of concept of their predictive power by inferring the
residual life-time of our specimen. We perform a retrospective failure
prediction using the cascade size $S$ measured during the damage accumulation
regime. To do so, we follow the idea of Anifrani et al. [69] and use the
methodology proposed by Mayya et al. [70]. Considering a time-series of
measured precursory activity for the equivalent force-controlled experiment,
the normalized distance to failure here writes as
$\delta=(t_{\mathrm{c}}-t)/t_{\mathrm{c}}$ where $t=0$ corresponds to
$F=F_{\mathrm{el}}$ and $t=t_{\mathrm{c}}$ corresponds to $F=F_{\mathrm{c}}$.
We rewrite the scaling law for the cascade size variations as

$\centering\langle S\rangle=S_{0}/\sqrt{t_{\mathrm{c}}-t}\@add@centering$ (4)

where $S_{0}$ is a constant. Re-arranging the terms, we obtain $\langle
S\rangle^{2}\,t=\langle S\rangle^{2}\,t_{\mathrm{c}}+S_{0}$, an expression
that can be used for performing a linear-regression of our experimental data
set $(t,S)$ shown in Fig. 6(a). The average size of the precursors is obtained
over a non-overlapping time window of $10~{}\mathrm{s}$. The prediction is
made at time $t_{\mathrm{cur}}$ so that only the precursors recorded at time
$t<t_{\mathrm{cur}}$ can be used for the prediction. Note, however, that we
only use a short period (here $100~{}\mathrm{s}$) before $t_{\mathrm{cur}}$ to
make the prediction.

The linear regression provides $t_{\mathrm{c}}^{\mathrm{pred}}$ that is shown
in Fig. 6(b) as a function of $t_{\mathrm{cur}}$. As shown in inset, the error
on the predicted failure time reduces as the prediction is made closer to the
actual failure time $t_{\mathrm{c}}$. The prediction lies within $10~{}\%$
error when the prediction is made in the last $25~{}\%$ of the total lifetime.
Note that the same methodology can be implemented using the duration or the
rate of events (under displacement controlled conditions only), thus providing
several independent measurements to forecast final failure. Interestingly, the
predictions are conservative by providing shorter residual lifetime than the
actual one. Importantly, the proposed methodology does not require monitoring
from the beginning of the damage accumulation phase. The robustness of our
method strongly argues that the statistical analysis of failure precursors
based on the proposed scenario is an efficient quantitative tool of damage
monitoring and lifetime prediction of structures.

##  VII Acoustic emissions

We now verify the generality of our results and their applications to standard
brittle solids. To this effect, we analyze the acoustic emissions accompanying
damage growth in our experiments (see Appendix A1 for experimental details). A
typical acoustic emission time-series record is shown in Fig. 7(a). First, we
observe bursts of acoustic emissions separated by silent periods, a behavior
typical of compressive failure of brittle disordered materials. Then, we see
that the acoustic activity intensifies closer to peak load, near localization,
a behavior that is also observed close to failure in standard brittle solids,
and that we explore in more details below. The distribution of acoustic event
energy is shown in Fig. 7(b) and in the inset. It follows a power law with an
exponent $\beta_{\mathrm{tot}_{\mathrm{ae}}}\simeq 1.45$ that is slightly
smaller when we consider only the acoustic emissions close to failure
$\beta_{\mathrm{ae}}\simeq 1.21$. Tracking now the evolution of the acoustic
activity, we find that the energy $\langle A_{\mathrm{ae}}\rangle$ of the
acoustic events increases as $\langle
A_{\mathrm{ae}}\rangle\propto\delta^{-\alpha_{A_{\mathrm{ae}}}}$ on
approaching failure, see Fig. 7(c). The acoustic activity rate
$dN_{\mathrm{ae}}/dt\propto\delta^{-\alpha_{N_{\mathrm{ae}}}}$ also increases
as a power law on approaching failure, see Fig. 7 (d). As a result, the rate
of acoustic energy $dE_{\mathrm{ae}}/dt$ defined as the product of the average
acoustic event energy with the activity rate
$dE_{\mathrm{ae}}/dt\propto\langle A_{\mathrm{ae}}\rangle\cdot
dN_{\mathrm{ae}}/dt$ increases as
$dE_{\mathrm{ae}}\propto\delta^{-\alpha_{\mathrm{ae}}}$ where
$\alpha_{\mathrm{ae}}=\alpha_{A_{\mathrm{ae}}}+\alpha_{N_{\mathrm{ae}}}\simeq
1.35$. A comparison of the exponents measured in our experiments with the one
reported in the literature shows that damage spreading in our experimental
system shares strong similarities with damage spreading in standard brittle
solids. In particular, the exponents $\beta_{\mathrm{tot}_{\mathrm{ae}}}\simeq
1.4~{}\mathrm{and}~{}\alpha_{N_{\mathrm{ae}}}\simeq 0.7$ measured in our
experiments are similar to the one measured in rocks [23, 71, 72], coal [73]
and concrete [26]. Beyond confirming the applicability of our experimental and
theoretical findings to a broad range of brittle materials, these results
confirm versatility of acoustic emissions as a damage monitoring tool.

Figure 7: (a) Acoustic emissions recorded during a typical experiment of
compressive failure. (b) Distribution of acoustic event energy recorded during
the whole duration of experiments (inset) and close to failure (main panel).
Variations of the average value of (c) the acoustic energy $A_{\mathrm{ae}}$
and (d) their activity rate $dN_{\mathrm{ae}}/dt$ with distance to failure
$\delta$.

We would like to discuss the intricate connection between the acoustic
emissions [74] and the mechanical energy cascades that have been thoroughly
characterized in our experiments (see Secs. IIB and IIC). We clearly observe
here that acoustic events do not correspond to the damage cascades
characterizing the mechanical response of the specimen. First, we find that
their number (typically $10^{3}$) during an experiment is much larger than the
number of mechanical precursors (typically $10^{2}$). Then, we observe that
the exponent characterizing the increase of emitted acoustic energy on
approaching failure $\alpha_{\mathrm{AE}}$ is much larger than the exponent
$\alpha$ describing the increase of dissipated rate of mechanical energy.
Finally, we observe that the rate of acoustic events increase on approaching
failure, an observation in stark contrast with the constant activity rate
observed for damage cascades, see Fig. 4(b).
We infer that the acoustic emissions may rather be reminiscent of the
dissipative processes taking place at much lower spatio-temporal scales within
damage cascades. Such a scenario also agrees with the power-law distribution
of inter-event waiting time of acoustic events [23, 52] that contrasts with
the exponential distribution observed here for damage cascades, see Fig. 2(f).
A deeper analysis of our experimental data may be needed to understand the
connection between mechanical precursors and acoustic events. Nevertheless,
their increase on approaching failure may be harnessed for tracking damage
growth and to predict the residual lifetime of structures, after calibration
of material specific exponents [70].

##  VIII Numerical solution of the damage evolution equation

To conclusively validate our interpretation of the non-stationary avalanche
dynamics preceding compressive failure, we numerically solve the damage
evolution equation (1) using a 2D cellular automata (see Appendix B) for both
force (SI Sec. S4A [46]) and displacement (SI Sec. S4B [46]) control
conditions. We recover that intermittent damage accumulation culminates in a
bifurcation that manifests as the emergence of a localization band at peak
load (SI Figs. S4(a)-(b) & S5(a)-(b) [46]). The exponents characterizing the
damage cascades are measured numerically (SI Figs. S4(c)-(j) & S5(c)-(j) [46])
using the methods employed for analyzing the experimental data, thus allowing
for a systematic comparison with the statistical features of the precursors
measured experimentally. Numerical and experimental exponents are provided in
Table 2. We also proceed to a comparison with exponent values reported in the
literature (see SI Table S2 [46]). The good agreement supports our theoretical
framework as an adequate description. Importantly, the similarities between
the statistics in force and equivalent-force control validate our method of
reconstruction of the precursors from our experimental data set and
interpretation of failure built from the analyses. We also verify numerically
the method employed to characterize the elastic interactions from the
incremental damage field in Fig.3(b).

Table 2: Exponents from experiments and numerical model.

| Definition | Experiments | Simulations†
---|---|---|---
$d_{\mathrm{f}}$ | $S\propto\xi^{d_{\mathrm{f}}}$ | 1.07 $\pm$ 0.07 | 1.15
$z$ | $T\propto\xi^{z}$ | 0.53 $\pm$ 0.11 | 0.62
$\theta$ | $P(\delta Y)\propto\delta Y^{\theta}$ | 0.24 $\pm$ 0.03 | 0.35 (0.18)
$\beta$ | $P(S)\propto S^{-\beta}$ | 1.30 $\pm$ 0.11 | 1.36 (1.34)
$\alpha$ | $S\propto\delta^{-\alpha}$ | 0.57 $\pm$ 0.04 | 0.48 (0.60)
$z/d_{\mathrm{f}}$ | $T\propto S^{z/d_{\mathrm{f}}}$ | 0.49 $\pm$ 0.14 | 0.53 (0.64)
$\beta_{tot}$ | $\beta_{tot}=\beta+\frac{2-\beta}{\alpha}$ | 2.32 $\pm$ 0.18 | 2.2 (2.13)

†Values in brackets are from equivalent force control scenario.

##  IX Implications and conclusion

We would like to conclude by highlighting the implications of our findings.
First, our results may be relevant for any brittle solids such as rocks,
concrete and mortar, the generality of our approach being bonded to damage
mechanics that we used to describe our experiments (see Annexe 2). By showing
that compressive failure is not a critical phase transition, but instead a
standard bifurcation, we pave the way for reliable predictions of the failure
load of brittle solids and the use of precursory events to anticipate their
forthcoming failure. Indeed, this finding implies that the stability analysis
of a homogeneous specimen is sufficient to predict the failure load of
disordered solids. This approach, used in this study, captures the emergence
of the localization band at peak load and, more generally, the stress-strain
response of the specimen. Further experimental validation of the application
of these concepts to more complex materials are a part of on-going studies.
Another direct implication of our finding is the prediction of the forthcoming
failure from the statistics of precursors. The interpretation of compressive
failure as a standard bifurcation provides robust scaling laws that can be
harnessed to predict failure.

We also would like to highlight the implications of our work for a larger
class of phenomena. Compressive failure is an archetypical example of
seemingly critical phenomena, as it is preceded by precursors with diverging
size and duration on approaching failure. We show that this apparent
criticality is not bounded to the presence of a critical point. Instead, we
show that a standard bifurcation in presence of disorder and long-range
interactions results in an intermittent response with fluctuations that
diverge on approaching the bifurcation point. We expect these ideas to be
relevant for understanding the relationship between intermittent plastic flow
and shear banding in amorphous solids, where the nature of the yielding
transition and the localization has been vigorously debated these recent years
[61, 75, 68]. Our interpretation of the specimen being at a finite distance
above the critical point as it evolves towards failure aligns well with the
conclusion drawn by Lin et al. [62]. Using a discrete model of amorphous
plasticity, they show that avalanche size diverges on approaching failure.
Yet, this divergence is not reminiscent of a critical behavior. Instead, they
show that the avalanche size is controlled by the slope
$\frac{d\sigma}{d\epsilon}$ of the stress-strain response. The scaling law
$\langle S\rangle\propto\frac{1}{d\sigma/d\epsilon}$ is compatible with the
one obtained in our experiments $\langle
S\rangle\propto\frac{1}{\sqrt{\delta}}$ (see Fig.4(a)). This further supports
the relevance of the scenario proposed here for sheared amorphous media.

In summary, we investigated experimentally the precursors to localization
during compressive failure of 2D cohesionless soft cellular solids. This model
system was shown to behave like a wide range of standard (and somehow more
complex) brittle materials under compression. We used this simplified brittle
disordered material to characterize in depth (i) the highly correlated
sequence of spatially coherent clusters of failure events that compose
precursors, (ii) the divergence of the size, the spatial extent and the
duration of these precursors on approaching failure, (iii) the presence of the
long-range elastic interactions decaying as $\sim 1/r^{2}$ that drives this
cooperative dynamics, (iv) the power-law distributed acoustic bursts
accompanying these damage cascades and finally (v) the emergence of a
localization band at peak load. We characterized the non-stationary statistics
of the damage cascades observed during the damage accumulation regime prior
localization. We then derived from continuum damage mechanics an evolution
equation of the damage field that was shown to account quantitatively for all
the scaling properties measured experimentally. The avalanche-dynamics of
damage growth was thus shown to be reminiscent of a non-stationary depinning
scenario that reconciles the two competing approaches proposed to describe
compressive failure, namely standard bifurcation and critical transition.

Ultimately, failure precursors have been shown to be merely the by-products of
the progressive loss of stability of the specimen as it approaches
localization. Contrary to the critical transition scenario, specimens driven
towards failure stay at a fixed distance to criticality. Nevertheless, the
evolution of the statistical features of precursors can be harnessed to
anticipate and even predict the forthcoming failure.

##  Appendix A Experimental characterization of precursors

###  A.1 Experimental set-up

Our specimen consists of a hexagonal packing of about 1500 soft cylinders made
of polypropylene placed in a transparent Plexiglas box of dimensions
$205\mathrm{mm}\times 170\mathrm{mm}\times 30\mathrm{mm}$, as shown in Fig.
1(a). The cylinders are $25~{}\mathrm{mm}$ long with a $5\mathrm{mm}$
diameter. Displacement is applied to the upper layer through a Plexiglas beam
using an AG-X Shimadzu test machine. The force experienced by the specimen is
measured using a $10\mathrm{kN}$ load cell and sampled at a rate of
$100~{}\mathrm{Hz}$. During the test, images are recorded every $0.1$ seconds
from the lateral side of the box, using a Baumer HXC20 camera with a
resolution of $2048\times 1088$ pixels and binarized using the open source
Fiji software [1]. Precursory activity recorded during ten experiments with
loading rate of $2~{}\mathrm{mm/min}$ and two experiments with loading rate of
$1~{}\mathrm{mm/min}$ were analyzed to determine the statistical structure of
precursors. The peak load measured in our experiments is
$\mathrm{F}_{c}=228\pm 4~{}N$. Acoustic emissions were recorded during four of
the ten experiments using two low frequency sensors (type R3$\alpha$ from
Mistras group, Europhysical Acoustics) that were placed on the compression
platens. The signals crossing the fixed threshold 27 dB from the transducers
are preamplified (40 dB) and transferred to a PC1-2 acquisition system. A
detailed analysis in a study of acoustic emissions during dislocations showed
that the scaling exponents are robust even for the range of the timing
parameters that may result in superposition of events [77]. Still, we account
for the cellular nature of the packing and set the event detection time post
the first crossing of threshold as well as the lockout time for signals after
the event as $500\mu s$. Calibration of the timing parameters was performed by
recording signals from pencil lead-breaks on metallic rods placed within the
cylindrical cells.

###  A.2 Cohesionless soft cellular solids under compression as model brittle
solids

We justify here why packings of cohesionless soft cells employed in our
experiments can be considered as model materials to investigate the
intermittent response of brittle disordered solids under compression. First,
we show that our packings of soft cells behave like standard elasto-damageable
media, a description that has been shown to accurately describe the
compressive response of brittle solids [5]. Then, we explain how friction
between cells that control the non-elastic response of our packings leads to a
mechanical response that is similar to the ones of brittle solids under
compression.

First, we validate the applicability of damage mechanics to describe the
average mechanical response of our packings of soft cells. We remind that
elasto-damageable solids behave elastically until some critical load level is
reached. Then, damage (the level of which is described by an internal variable
noted $d_{\circ}$) increases. The impact of damage on the mechanical response
of the material is described by the decrease of the material Young’s modulus
$E(d_{\circ})$ and the variations of its Poisson’s ratio $\nu(d_{\circ})$
where $d_{\circ}$ is the accumulated damage. In our experiments, $d_{\circ}$
is obtained from the deviation of the cell geometry to cell circularity, a
definition that is justified below.

The Young’s modulus of our packings is obtained from the instantaneous
specimen stiffness $E=(1-\nu^{2})\sigma_{\circ}/\epsilon_{\circ}$ leading to
the following quadratic decrease $E\simeq E_{\circ}-a_{\circ}d_{\circ}^{2}$.
Its Poisson’s ratio $\nu\simeq 0.26$, obtained from the ratio
$-\langle\epsilon_{\mathrm{xx}}\rangle/\langle\epsilon_{\mathrm{yy}}\rangle$
computed in the central region of the specimen, see Fig. 8(c), is found to be
constant throughout our experiments and nearly independent of the damage level
$d_{\circ}$.

The force-displacement response predicted by elasto-damageable theory is shown
in Fig. 8(a). It results from the balance
$Y[\sigma_{\circ}(d_{\circ}),E(d_{\circ}),\nu]=Y_{\mathrm{c}}(d_{\circ})$
between the elastic energy release rate $Y$ and the damage energy
$Y_{\mathrm{c}}$ (see SI Section 3.A for the detailed calculation).The
predicted behavior captures very well our experimental data. The macroscopic
stress $\sigma_{\circ}$ and strain $\epsilon_{\circ}$ are computed here as
parametric functions of the average damage level $d_{\circ}$

$\left\\{\begin{array}[]{lcl}\vspace{3pt}\sigma_{\circ}(d_{\circ})&=&\displaystyle\frac{E(d_{\circ})}{\sqrt{1-\nu^{2}}}\,\sqrt{\frac{-2Y_{c\circ}(d_{\circ})}{E^{\prime}(d_{\circ})}},\\\
\epsilon_{\circ}(d_{\circ})&=&\displaystyle\sqrt{1-\nu^{2}}\sqrt{\frac{-2Y_{c\circ}(d_{\circ}}{E^{\prime}(d_{\circ})}}.\end{array}\right.$
(5)

The average damage resistance $Y_{c\circ}(d_{\circ})$ is inferred from the
displacement fields measured by tracking the cell positions, following the
procedure described in Appendix A4\. Note that beyond capturing the average
mechanical response of the packings of soft cells used in our experiments,
elasto-damageable theory captures the fluctuations around this average
response, and in particular the statistics of failure precursors in its most
subtle details. The application of the elasto-damageable theory to disordered
brittle solids is described in Section 3 of SI and the comparison of the
predictions derived from this model with our experimental observations is
summarized in Table LABEL:Table2. Note also that the long-range interactions
$\sim 1/r^{2}$, a unique feature of elasto-damageable solids, is retrieved in
our experiments (see Fig. 3(a)-(c)). Last but not least, elasto-damageable
theory predicts (see SI Section 3F) the emergence of localization at peak
load, a prediction in perfect agreement with our experimental observations
(see Appendix A.5. for the experimental determination of the onset of damage
localization).

Figure 8: Suitability of damage mechanics to describe the average mechanical
response of the cohesionless frictional cellular solids used in our
experiments. (a) Comparison of the mechanical response predicted by damage
mechanics with the experimental stress-strain curve. (b) Variations of the
elastic modulus and the average damage level during a representative
experiment as inferred from the global mechanical response and our tracking of
the cells position. (c) Determination of the Poisson’s ratio $\nu$ inferred
from the linear fit of the lateral strain vs the normal one after averaging
over the central gray region of the specimen shown in inset.
$\nu\simeq\nu_{\circ}(d_{\circ})=0.26$ is found to be nearly constant and
independent of $d_{\circ}$

The applicability of damage mechanics to our packings of soft cells support
the generality of our findings and their applicability to a broad range of
brittle disordered solids. We now would like to provide a microscopic
interpretation to this seemingly surprising behavior. The inelastic behavior
of our packings of soft cells is governed by the friction between neighboring
cells. This rearrangement of the cells explains the local softening of the
material that is appropriately descried by a local decrease of the Young’s
modulus $E(d_{\circ})$ within the theoretical framework of damage mechanics.
Meanwhile, cells rearrangements alter the geometry of the initially circular
cell. As a result, the deviation to cell circularity is an appropriate of
these rearrangements and the associated material softening, justifying the
definition of the damage level $d_{\circ}$ employed in our study.

Last but not least, friction is very often considered as the central damage
mechanism in standard brittle solids under compression. This is clearly
evidenced by the large success of the Mohr-Coulomb theory (a model that
directly derives from the assumption that the different regions of the
material can slide on each other) to describe the inelastic response of
brittle materials under compression.

Interestingly, Karimi et al. [48] previously noticed in their numerical work
that soft granular solids can be described by damage mechanics. They also
observed in their simulations that such systems failed through the emergence
of a localization band at peak load and even noticed that the size of the
failure precursors diverge as a power-law, two observations in agreement with
our experimental observations.

Overall, these various observations come in support to our central assumption:
packings of frictional soft particles behaves as elasto-damageable solids, a
conclusion already drawn by Houdoux  _et al._ [78, 79] from compression
experiments carried on standard granular solids.

###  A.3 Precursor size in an equivalent force control scenario

We explain here why the energy dissipated by damage during a precursory damage
event follows the expression
$S_{\mathrm{global}}=(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})F_{0}/2$. We
compute first the elastic energy stored in the specimen right before (for
$\Delta=\Delta_{\mathrm{ini}}$ ) and right after (for
$\Delta=\Delta_{\mathrm{end}}$) the precursory damage event. They follow
$E_{\mathrm{el,ini}}=\Delta_{\mathrm{ini}}F_{0}/2$ and
$E_{\mathrm{el,end}}=\Delta_{\mathrm{end}}F_{0}/2$, respectively. This leads
to the following expression $\Delta
E_{\mathrm{el}}=E_{\mathrm{el,end}}-E_{\mathrm{el,ini}}=(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})F_{0}/2$
for the variation of the total elastic energy stored in the specimen during
one precursory event. Note that this quantity is positive, meaning that the
stored elastic energy does increase during a precursory event. We now compute
the work of the external force during such an event. As the force is constant
during the considered force-controlled scenario, we obtain $\Delta
W=(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})F_{0}$. This quantity is
defined positive and corresponds to the energy injected by the loading machine
within the specimen during a precursory event. The determination of the energy
dissipated by damage during a precursor is now in order. Observing that the
energy injected by the loading machine during one precursor compensates
exactly the new elastic energy stored in the specimen and the one dissipated
by damage, one obtains $\Delta W=\Delta
E_{\mathrm{el}}+S_{\mathrm{global}}\Rightarrow S_{\mathrm{global}}=\Delta
W-\Delta
E_{\mathrm{el}}=(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})F_{0}/2$.
Interestingly, this balance of energy during one precursor can be understood
graphically from the force-displacement response of the specimen during the
event. As shown in Fig. 1(c), the work $\displaystyle\Delta
W=\int_{\Delta_{\mathrm{ini}}}^{\Delta_{\mathrm{end}}}F(\Delta)\mathrm{d}\Delta=F_{0}\int_{\Delta_{\mathrm{ini}}}^{\Delta_{\mathrm{end}}}\mathrm{d}\Delta=(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})F_{0}$
of the external force during one precursor corresponds to the area of the
rectangle located below the force plateau. This rectangle can be divided in
two other rectangles of similar size
$(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})F_{0}/2$, one of them
corresponding to the additional stored energy $\Delta E_{\mathrm{el}}$
(represented in grey in Fig. 1(c)) while the other one corresponding to the
dissipated energy $\Delta E_{\mathrm{d}}=S_{\mathrm{global}}$ (represented in
yellow). This equipartition of energy between the new stored elastic energy
and the dissipated one provides a simple way to track the size of failure
precursors during the phase of damage accumulation.

###  A.4 Tracking damage evolution

The continuous image acquisition allows for the tracking of both position and
circularity (used here to define damage) of individual cells. We then compute
the displacement fields $u_{x}(\vec{x},t)~{}\mathrm{and}~{}u_{y}(\vec{x},t)$
and the damage field $d(\vec{x},t)$ of the effective elasto-damageable medium
following the coarse-graining procedure described in Glasser and Goldhirsch
[80]. To resolve the damage cascades from the image stack, we refer to the
macroscopic response where an equivalent force control scenario is constructed
from the sequence of load drops. We thus obtain the start and end of each
damage cascade that are noted $t(\Delta_{\mathrm{ini}})$ and
$t(\Delta_{\mathrm{end}})$, respectively. The damage events belonging to the
image stack constituting the cascades are then clustered based on their
connectivity (26-connected neighborhood). This procedure provides a space-time
dissipation map of failure precursors composed of a series of highly
correlated clusters. We derive the strain field as
$\epsilon_{xx}(\vec{x})=\frac{du_{x}}{dx},~{}\epsilon_{yy}(\vec{x})=\frac{du_{y}}{dy}~{}\mathrm{and}~{}\epsilon_{xy}(\vec{x})=\frac{1}{2}\left(\frac{du_{x}}{dy}+\frac{du_{x}}{dy}\right)$.
The field of elastic energy stored in the specimen per unit volume is then
obtained from the relation,

$\begin{split}E_{\mathrm{el}}^{\mathrm{local}}(\vec{x},t)=&\frac{E_{\mathrm{local}}(\vec{x},t)}{2(1-\nu^{2})}\left[\epsilon_{\mathrm{xx}}^{2}(\vec{x},t)+\right.\\\
&\epsilon_{\mathrm{yy}}^{2}(\vec{x},t)+2\,\nu\epsilon_{\mathrm{xx}}(\vec{x},t)\epsilon_{\mathrm{yy}}(\vec{x},t)+\\\
&\left.2(1-\nu)\epsilon_{\mathrm{xy}}^{2}(\vec{x},t)\right].\end{split}$ (6)

The field of damage driving force is obtained as
$Y[d(\vec{x},t)]=\frac{dE_{\mathrm{el}}^{\mathrm{local}}}{dd}$. The field of
damage resistance $Y_{c}[d(\vec{x},t)]$ is inferred from the field of damage
driving force $Y[d(\vec{x},t)]$ computed in a retrospective manner. If damage
grows locally at $(\vec{x}_{0},t_{i})$, we assign
$Y_{c}[d(\vec{x}_{0}),t_{i}]=Y[d(\vec{x}_{0}),t_{i}]$ as well as for the
preceding time-steps $(t_{\mathrm{i-}})$. We find that the average damage
resistance $Y_{\mathrm{c\circ}}(d_{\circ})$ increases with the damage level
$d_{\circ}$. A linear approximation of the hardening behavior
$Y_{\mathrm{c\circ}}(d_{\circ})\sim\eta d_{\circ}$ with a hardening
coefficient $\eta\simeq 45$ provides a good description of our experimental
data.

###  A.5 Onset of damage localization

We provide here the methodology employed to determine the localization
threshold. The basic idea is to compute the field of incremental damage growth
at different distances to peak load as shown in Fig. 9(a) to determine whether
it is homogeneous or localized. The incremental damage field $\delta
d(\vec{x})=d(\vec{x},\delta)-d(\vec{x},\delta+\epsilon)$ where
$\epsilon\ll\delta$, is shown for three different values of $\delta$ : right
before peak load ($\delta=0.02$), at peak load ($\delta=0$) and right after
($\delta=-0.02$). We clearly see the localization of the damage growth in thin
band for $\delta=-0.02$ from which we infer that the localization starts at
peak load. The standard deviation $\sigma_{\delta d}(x)=\sqrt{\langle(\delta
d(x)-\langle\delta(\vec{x})\rangle_{\mathrm{y}})^{2}\rangle_{\mathrm{y}}}$ is
shown in Fig. 9(c) for several values of $\delta$. Here also, we see that
$\sigma_{\delta d}(x)$ increases above $0.01$, suggesting that localization
takes place at peak load.

Figure 9: The incremental damage field $\delta d(\vec{\mathrm{x}})$ is
computed for several locations $\delta$ near peak load as indicated by solid
circles on the force-displacement curve of panel (a). (b) $\delta
d(\vec{\mathrm{x}})$ is represented right before ($\delta=0.02$), at peak
($\delta=0$) and right after ($\delta=-0.02$). Note the localization of the
damage activity in the third panel suggesting that localization took place at
peak load. (c) This is further confirmed by the standard deviation
$\sigma_{\delta d}(x)$ of the incremental damage field that is shown to
increase for $\delta<0$.

##  Appendix B Numerical modeling of intermittent damage evolution and
localization

To solve the damage evolution equation under quasi-static loading conditions,
we adopt the procedure described in Berthier et al. [45]. We consider a grid
of size $L=21$ discretized into $L^{2}$ elements with periodic boundary
conditions. The interaction kernel is derived for the case of uni-axial
compression using the method described in Dansereau et al. [58] and writes as

$\psi(d_{\circ})=\left[\frac{E^{\prime}(d_{\circ})^{2}}{E(d_{\circ})^{3}}\right](1-\nu^{2})\sigma_{ext}^{2}\left[\frac{x^{4}-3y^{4}+6x^{2}y^{2}}{4\pi(x^{2}+y^{2})^{3}}\right].$
(7)

Interestingly, the kernel is independent of the type of loading configuration
(stress vs. strain imposed loading conditions). The amplitude of the kernel is
also non-stationary, in contrast with the non-positive interaction used to
describe amorphous plasticity [61, 62]. The variations of the elastic constant
with the damage level is found to follow
$E_{\circ}=E^{\circ}(1-d_{\circ}^{2})$ with $E^{\circ}=1.0~{}\mathrm{MPa}$. We
consider a heterogeneous field of damage resistance $Y_{c}[d(\vec{x},t)]$ that
evolves with the damage level $d(\vec{x},t)$ following the linear hardening
law $Y_{c\circ}=Y_{c}^{\circ}(1+\eta d_{\circ})$, where $Y_{c}^{\circ}\simeq
1.4~{}\mathrm{kJ/m^{3}}$ and $\eta\simeq 45$ are measured experimentally. The
stress (strain) is gradually increased until the damage criterion is fulfilled
for one of the elements $\vec{x}=\vec{x}_{0}$. The damage is increased locally
and the values of driving force $Y(\vec{x}_{0},t)$ and resistance
$Y_{c}(\vec{x}_{0})$ are updated. The non-local redistribution of driving
force given by the kernel $\psi(d_{\circ})$ may then trigger a cascade of
damage events. The cascade stops when all elements are stable following which,
we increase the stress (strain) again.

The reader is invited to refer to the Supplementary Information [46] for the
thorough description of the local analyses of experimental precursors, the
theoretical and numerical modeling of damage accumulation and the analogy with
driven disordered elastic interface.

###### Acknowledgements.

The authors gratefully acknowledge financial support from Sorbonne Université
through the Emergence grant for the research project, From damage spreading to
failure in quasi-brittle materials as well as CNRS and Satt-Lutech through the
tech transfer project, Development of a technology of predictive maintenance
for materials and structures under compression.

## References

  * Koivisto  _et al._ [2016] J. Koivisto, M. Ovaska, A. Miksic, L. Laurson, and M. J. Alava, Predicting sample lifetimes in creep fracture of heterogeneous materials, Physical Review E 94, 023002 (2016).
  * Kádár  _et al._ [2020] V. Kádár, G. Pál, and F. Kun, Record statistics of bursts signals the onset of acceleration towards failure, Scientific reports 10, 1 (2020).
  * Biswas  _et al._ [2020] S. Biswas, D. Fernandez Castellanos, and M. Zaiser, Prediction of creep failure time using machine learning, Scientific Reports 10, 1 (2020).
  * Debski  _et al._ [2021] W. Debski, S. Pradhan, and A. Hansen, Criterion for imminent failure during loading—discrete element method analysis, Frontiers in Physics 9, 223 (2021).
  * Ashby and Sammis [1990] M. Ashby and C. Sammis, The damage mechanics of brittle solids in compression, Pure and Applied Geophysics 133, 489 (1990).
  * Lockner  _et al._ [1991] D. Lockner, J. Byerlee, V. Kuksenko, A. Ponomarev, and A. Sidorin, Quasi-static fault growth and shear fracture energy in granite, Nature 350, 39 (1991).
  * Kachanov [1987] M. Kachanov, Elastic solids with many cracks: a simple method of analysis, International Journal of Solids and Structures 23, 23 (1987).
  * Mazars and Pijaudier-Cabot [1989] J. Mazars and G. Pijaudier-Cabot, Continuum damage theory—application to concrete, Journal of engineering mechanics 115, 345 (1989).
  * Kachanov [1993] M. Kachanov, Elastic solids with many cracks and related problems, Advances in applied mechanics 30, 259 (1993).
  * Fortin  _et al._ [2006] J. Fortin, S. Stanchits, G. Dresen, and Y. Guéguen, Acoustic emission and velocities associated with the formation of compaction bands in sandstone, J. Geophys. Res. Solid Earth 111 (2006).
  * Fortin  _et al._ [2009] J. Fortin, S. Stanchits, G. Dresen, and Y. Gueguen, Acoustic emissions monitoring during inelastic deformation of porous sandstone: comparison of three modes of deformation, Pure and Applied Geophysics 166, 823 (2009).
  * Manzato  _et al._ [2012] C. Manzato, A. Shekhawat, P. K. V. V. Nukala, M. J. Alava, J. P. Sethna, and S. Zapperi, Fracture strength of disordered media: Universality, interactions, and tail asymptotics, Phys. Rev. Lett. 108, 065504 (2012).
  * Manzato  _et al._ [2014] C. Manzato, M. J. Alava, and S. Zapperi, Damage accumulation in quasibrittle fracture, Phys. Rev. E 90, 012408 (2014).
  * Tal  _et al._ [2016] Y. Tal, B. Evans, and U. Mok, Direct observations of damage during unconfined brittle failure of carrara marble, J. Geophys. Res. Solid Earth 121, 1584 (2016).
  * Thilakarathna  _et al._ [2020] P. Thilakarathna, K. K. Baduge, P. Mendis, V. Vimonsatit, and H. Lee, Mesoscale modelling of concrete–a review of geometry generation, placing algorithms, constitutive relations and applications, Engineering Fracture Mechanics 231, 106974 (2020).
  * Rudnicki and Rice [1975] J. W. Rudnicki and J. Rice, Conditions for the localization of deformation in pressure-sensitive dilatant materials, Journal of the Mechanics and Physics of Solids 23, 371 (1975).
  * Olsson [1999] W. A. Olsson, Theoretical and experimental investigation of compaction bands in porous rock, J. Geophys. Res. Solid Earth 104, 7219 (1999).
  * Lockner [1993] D. Lockner, The role of acoustic emission in the study of rock fracture, in  _International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts_, Vol. 30 (Elsevier, 1993) pp. 883–899.
  * Petri  _et al._ [1994] A. Petri, G. Paparo, A. Vespignani, A. Alippi, and M. Costantini, Experimental evidence for critical dynamics in microfracturing processes, Phys. Rev. Lett. 73, 3423 (1994).
  * Deschanel  _et al._ [2006] S. Deschanel, L. Vanel, G. Vigier, N. Godin, and S. Ciliberto, Statistical properties of microcracking in polyurethane foams under tensile test, influence of temperature and density, Int. J. Frac. 140, 87 (2006).
  * Davidsen  _et al._ [2007] J. Davidsen, S. Stanchits, and G. Dresen, Scaling and universality in rock fracture, Phys. Rev. Lett. 98, 125502 (2007).
  * Rosti  _et al._ [2009] J. Rosti, X. Illa, J. Koivisto, and M. J. Alava, Crackling noise and its dynamics in fracture of disordered media, J. Phys. D: Appl. Phys. 42, 214013 (2009).
  * Baró  _et al._ [2013] J. Baró, Á. Corral, X. Illa, A. Planes, E. K. Salje, W. Schranz, D. E. Soto-Parra, and E. Vives, Statistical similarity between the compression of a porous material and earthquakes, Phys. Rev. Lett. 110, 088702 (2013).
  * Castillo-Villa  _et al._ [2013] P. Castillo-Villa, J. Baro, A. Planes, E. K. H. Salje, P. Sellappan, W. M. Kriven, and E. Vives, Crackling noise during failure of alumina under compression: The effect of porosity, J. Phys. Cond. Mat. 25, 292202 (2013).
  * Baró  _et al._ [2018] J. Baró, K. A. Dahmen, J. Davidsen, A. Planes, P. O. Castillo, G. F. Nataf, E. K. Salje, and E. Vives, Experimental evidence of accelerated seismic release without critical failure in acoustic emissions of compressed nanoporous materials, Phys. Rev. Lett. 120, 245501 (2018).
  * Vu  _et al._ [2019] C.-C. Vu, D. Amitrano, O. Plé, and J. Weiss, Compressive failure as a critical transition: Experimental evidence and mapping onto the universality class of depinning, Phys. Rev. Lett. 122, 015502 (2019).
  * Sornette [1994] D. Sornette, Sweeping of an instability: an alternative to self-organized criticality to get power laws without parameter tuning, J. Phys. I 4, 209 (1994).
  * Zapperi  _et al._ [1997a] S. Zapperi, P. Ray, H. E. Stanley, and A. Vespignani, First-order transition in the breakdown of disordered media, Phys. Rev. Lett. 78, 1408 (1997a).
  * Zapperi  _et al._ [1997b] S. Zapperi, A. Vespignani, and H. E. Stanley, Plasticity and avalanche behaviour in microfracturing phenomena, Nature 388, 658 (1997b).
  * Sornette [2002] D. Sornette, Predictability of catastrophic events: Material rupture, earthquakes, turbulence, financial crashes, and human birth, Proceedings of the National Academy of Sciences 99, 2522 (2002).
  * da Rocha and Truskinovsky [2020] H. B. da Rocha and L. Truskinovsky, Rigidity-controlled crossover: From spinodal to critical failure, Physical Review Letters 124, 015501 (2020).
  * Guarino  _et al._ [1998] A. Guarino, A. Garcimartin, and S. Ciliberto, An experimental test of the critical behaviour of fracture precursors, Eur. Phys. J. B 6, 13 (1998).
  * Girard  _et al._ [2010] L. Girard, D. Amitrano, and J. Weiss, Failure as a critical phenomenon in a progressive damage model, Journal of Statistical Mechanics: Theory and Experiment 2010, P01013 (2010).
  * Kun  _et al._ [2014] F. Kun, I. Varga, S. Lennartz-Sassinek, and I. G. Main, Rupture cascades in a discrete element model of a porous sedimentary rock, Phys. Rev. Lett. 112, 165501 (2014).
  * Kandula  _et al._ [2019] N. Kandula, B. Cordonnier, E. Boller, J. Weiss, D. K. Dysthe, and F. Renard, Dynamics of microscale precursors during brittle compressive failure in carrara marble, J. Geophys. Res. Solid Earth 124, 6121 (2019).
  * Roux  _et al._ [1988] S. Roux, A. Hansen, H. Herrmann, and E. Guyon, Rupture of heterogeneous media in the limit of infinite disorder, Journal of statistical physics 52, 237 (1988).
  * Delaplace  _et al._ [1996] A. Delaplace, G. Pijaudier-Cabot, and S. Roux, Progressive damage in discrete models and consequences on continuum modelling, Journal of the Mechanics and Physics of Solids 44, 99 (1996).
  * Garcimartin  _et al._ [1997] A. Garcimartin, A. Guarino, L. Bellon, and S. Ciliberto, Statistical properties of fracture precursors, Phys. Rev. Lett. 79, 3202 (1997).
  * Moreno  _et al._ [2000] Y. Moreno, J. B. Gòmez, and A. F. Pacheco, Fracture and second-order phase transitions, Phys. Rev. Lett. 85, 2865 (2000).
  * Weiss  _et al._ [2014] J. Weiss, L. Girard, F. Gimbert, D. Amitrano, and D. Vandembroucq, (finite) statistical size effects on compressive strength, Proceedings of the National Academy of Sciences 111, 6231 (2014).
  * Renard  _et al._ [2018] F. Renard, J. Weiss, J. Mathiesen, Y. Ben-Zion, N. Kandula, and B. Cordonnier, Critical evolution of damage toward system-size failure in crystalline rock, J. Geophys. Res. Solid Earth 123, 1969 (2018).
  * Barabási and Stanley [1995] A. L. Barabási and H. E. Stanley,  _Fractal concepts in surface growth_ (Cambridge University Press, 1995).
  * Sethna  _et al._ [2001] J. Sethna, K. Dahmen, and C. R. Myers, Crackling noise, Nature 410, 242 (2001).
  * Wiese [2021] K. J. Wiese, Theory and experiments for disordered elastic manifolds, depinning, avalanches, and sandpiles, arXiv preprint arXiv:2102.01215 . (2021).
  * Berthier  _et al._ [2022] E. Berthier, A. Mayya, and L. Ponson, Damage spreading in quasi-brittle disordered solids: Ii. what the statistics of precursors teach us about compressive failure, Journal of the Mechanics and Physics of Solids 162, 104826 (2022).
  * Mayya  _et al._ [2023] A. Mayya, E. Berthier, and L. Ponson, Supplementary information for “how criticality meets bifurcation in compressive failure of disordered solid” (2023).
  * Poirier  _et al._ [1992] C. Poirier, M. Ammi, D. Bideau, and J. P. Troadec, Experimental study of the geometrical effects in the localization of deformation, Phys. Rev. Lett. 68, 216 (1992).
  * Karimi  _et al._ [2019] K. Karimi, D. Amitrano, and J. Weiss, From plastic flow to brittle fracture: Role of microscopic friction in amorphous solids, Physical Review E 100, 012908 (2019).
  * Renard  _et al._ [2019] F. Renard, J. McBeck, N. Kandula, B. Cordonnier, P. Meakin, and Y. Ben-Zion, Volumetric and shear processes in crystalline rock approaching faulting, Proceedings of the National Academy of Sciences 116, 16234 (2019).
  * Cartwright-Taylor  _et al._ [2020] A. Cartwright-Taylor, I. G. Main, I. B. Butler, F. Fusseis, M. Flynn, and A. King, Catastrophic failure: How and when? insights from 4-d in situ x-ray microtomography, J. Geophys. Res. Solid Earth 125, e2020JB019642 (2020).
  * Lennartz-Sassinek  _et al._ [2014] S. Lennartz-Sassinek, I. Main, M. Zaiser, and C. Graham, Acceleration and localization of subcritical crack growth in a natural composite material, Phys. Rev. E 90, 052401 (2014).
  * Davidsen  _et al._ [2017] J. Davidsen, G. Kwiatek, E.-M. Charalampidou, T. Goebel, S. Stanchits, M. Rück, and G. Dresen, Triggering processes in rock fracture, Phys. Rev. Lett. 119, 068501 (2017).
  * Baud  _et al._ [2004] P. Baud, E. Klein, and T.-f. Wong, Compaction localization in porous sandstones: spatial evolution of damage and acoustic emission activity, Journal of Structural Geology 26, 603 (2004).
  * Lenoir  _et al._ [2007] N. Lenoir, M. Bornert, J. Desrues, P. Bésuelle, and G. Viggiani, Volumetric digital image correlation applied to x-ray microtomography images from triaxial compression tests on argillaceous rock, Strain 43, 193 (2007).
  * McBeck  _et al._ [2020] J. A. McBeck, J. M. Aiken, J. Mathiesen, Y. Ben-Zion, and F. Renard, Deformation precursors to catastrophic failure in rocks, Geophysical Research Letters 47, e2020GL090255 (2020).
  * Kandula  _et al._ [2022] N. Kandula, J. McBeck, B. Cordonnier, J. Weiss, D. K. Dysthe, and F. Renard, Synchrotron 4d x-ray imaging reveals strain localization at the onset of system-size failure in porous reservoir rocks, Pure and Applied Geophysics 179, 325 (2022).
  * Berthier  _et al._ [2017] E. Berthier, V. Démery, and L. Ponson, Damage spreading in quasi-brittle disordered solids: I. localization and failure, Journal of the Mechanics and Physics of Solids 102, 101 (2017).
  * Dansereau  _et al._ [2019] V. Dansereau, V. Démery, E. Berthier, J. Weiss, and L. Ponson, Collective damage growth controls fault orientation in quasibrittle compressive failure, Phys. Rev. Lett. 122, 085501 (2019).
  * Hentschel and Procaccia [1983] H. Hentschel and I. Procaccia, The infinite number of generalized dimensions of fractals and strange attractors, Physica D: Nonlinear Phenomena 8, 435 (1983).
  * Amitrano [2012] D. Amitrano, Variability in the power-law distributions of rupture events, Eur. Phys. J.: Spec. Top. 205, 199 (2012).
  * Lin  _et al._ [2014] J. Lin, E. Lerner, A. Rosso, and M. Wyart, Scaling description of the yielding transition in soft amorphous solids at zero temperature, Proceedings of the National Academy of Sciences 111, 14382 (2014).
  * Lin  _et al._ [2015] J. Lin, T. Gueudré, A. Rosso, and M. Wyart, Criticality in the approach to failure in amorphous solids, Phys. Rev. Lett. 115, 168001 (2015).
  * Lawn [1993] B. Lawn,  _Fracture of brittle solids_ (Cambridge University Press, 1993).
  * Berthier [2015] E. Berthier,  _Quasi-brittle failure of heterogeneous materials : Damage statistics and localization_ , Ph.D. thesis, Université Pierre et Marie Curie (2015).
  * Lin and Wyart [2016] J. Lin and M. Wyart, Mean-field description of plastic flow in amorphous solids, Phys. Rev. X 6, 011005 (2016).
  * Liu  _et al._ [2016] C. Liu, E. E. Ferrero, F. Puosi, J.-L. Barrat, and K. Martens, Driving rate dependence of avalanche statistics and shapes at the yielding transition, Phys. Rev. Lett. 116, 065501 (2016).
  * Budrikis  _et al._ [2017] Z. Budrikis, D. F. Castellanos, S. Sandfeld, M. Zaiser, and S. Zapperi, Universal features of amorphous plasticity, Nature communications 8, 1 (2017).
  * Ozawa  _et al._ [2018] M. Ozawa, L. Berthier, G. Biroli, A. Rosso, and G. Tarjus, Random critical point separates brittle and ductile yielding transitions in amorphous materials, Proceedings of the National Academy of Sciences 115, 6656 (2018).
  * Anifrani  _et al._ [1995] J.-C. Anifrani, C. Le Floc’h, D. Sornette, and B. Souillard, Universal log-periodic correction to renormalization group scaling for rupture stress prediction from acoustic emissions, Journal de Physique I 5, 631 (1995).
  * Mayya  _et al._ [2020] A. Mayya, E. Berthier, and L. Ponson, Procédé et dispositif d’analyse d’une structure. french patent application fr2002824 (2020).
  * Nataf  _et al._ [2014] G. F. Nataf, P. O. Castillo-Villa, J. Baró, X. Illa, E. Vives, A. Planes, and E. K. Salje, Avalanches in compressed porous si o 2-based materials, Phys. Rev. E 90, 022405 (2014).
  * Davidsen  _et al._ [2021] J. Davidsen, T. Goebel, G. Kwiatek, S. Stanchits, J. Baró, and G. Dresen, What controls the presence and characteristics of aftershocks in rock fracture in the lab?, Journal of Geophysical Research: Solid Earth 126, e2021JB022539 (2021).
  * Xu  _et al._ [2019] Y. Xu, A. G. Borrego, A. Planes, X. Ding, and E. Vives, Criticality in failure under compression: Acoustic emission study of coal and charcoal with different microstructures, Phys. Rev. E 99, 033001 (2019).
  * Zreihan  _et al._ [2019] N. Zreihan, E. Faran, E. Vives, A. Planes, and D. Shilo, Relations between stress drops and acoustic emission measured during mechanical loading, Phys. Rev. Mater. 3, 043603 (2019).
  * Parisi  _et al._ [2017] G. Parisi, I. Procaccia, C. Rainone, and M. Singh, Shear bands as manifestation of a criticality in yielding amorphous solids, Proceedings of the National Academy of Sciences 114, 5577 (2017).
  * Schindelin  _et al._ [2012] J. Schindelin, I. Arganda-Carreras, E. Frise, V. Kaynig, M. Longair, T. Pietzsch, S. Preibisch, C. Rueden, S. Saalfeld, B. Schmid, _et al._ , Fiji: an open-source platform for biological-image analysis, Nature methods 9, 676 (2012).
  * Lebyodkin  _et al._ [2013] M. Lebyodkin, I. Shashkov, T. Lebedkina, K. Mathis, P. Dobron, and F. Chmelik, Role of superposition of dislocation avalanches in the statistics of acoustic emission during plastic deformation, Physical Review E 88, 042402 (2013).
  * Houdoux  _et al._ [2018] D. Houdoux, T. B. Nguyen, A. Amon, and J. Crassous, Plastic flow and localization in an amorphous material: experimental interpretation of the fluidity, Physical Review E 98, 022905 (2018).
  * Houdoux  _et al._ [2021] D. Houdoux, A. Amon, D. Marsan, J. Weiss, and J. Crassous, Micro-slips in an experimental granular shear band replicate the spatiotemporal characteristics of natural earthquakes, Communications Earth & Environment 2, 90 (2021).
  * Glasser and Goldhirsch [2001] B. Glasser and I. Goldhirsch, Scale dependence, correlations, and fluctuations of stresses in rapid granular flows, Physics of Fluids 13, 407 (2001).

## Supplementary Information

Table 1: List of notations used in the manuscript and supplementary information $L,H,b$ | : Length, height and thickness of the specimen
---|---
$F$ | : Force experienced by the specimen as recorded at the macroscopic scale
$F_{\mathrm{el}},F_{\mathrm{c}}$ | : Elastic limit and peak force experienced by the specimen
$\Delta$ | : Applied displacement by the loading machine during the experiment
$E$ | : Elastic modulus of the specimen
$\nu$ | : Poisson’s ratio of the specimen
$\sigma_{\mathrm{ext}}$ | : Nominal stress experienced by the specimen during the experiment
$\epsilon_{\mathrm{ext}}$ | : Nominal strain experienced by the specimen solid during the experiment
$W$, $\Delta W$ | : Input work and incremental input work by the loading machine
$E_{el}$, $\Delta E_{el}$ | : Total elastic energy and incremental elastic energy in the specimen corresponding to input work
$E_{d}$, $\Delta E_{d}$ | : Total energy dissipated and incremental dissipation during damage growth in the specimen
${}_{\mathrm{global}},{}_{\mathrm{local}}$ | : Subscripts specifying the mode of analysis - macroscopic response (global) and at the cells (local)
$\Delta_{\mathrm{ini}}$, $\Delta_{\mathrm{end}}$ | : Values of displacement corresponding to the start and end of the damage cascade
in equivalent force control scenario
$\delta$ | : Normalized distance to failure
$\rho$ | : Dissipation energy density of the material element
$\rho_{\mathrm{th}}$ | : Thresholded value of the dissipation energy density of the material element
$S$ | : Size of the damage cascade
$\xi$ | : Spatial extent of the damage cascade
$T$ | : Duration of the damage cascade
$A$ | : Size of the load drops
$dN_{S}/dt$ | : Activity rate of the cascades
$dN_{A}/dt$ | : Activity rate of the load drops
$\tau_{w},\tau_{F}$ | : Waiting time between the cascades in terms of time and force
$\Delta d^{\ast}$ | : Characteristic damage value of the cascades
$d_{\mathrm{f}}$ | : Fractal dimension
$z$ | : Dynamic exponent
$\theta$ | : Exponent characterizing the marginal stability
$\beta$ | : Exponent characterizing the decay of distribution of cascade size, $S$ close to failure
$\beta_{tot}$ | : Exponent characterizing the decay of the stress-integrated distribution of cascade size, $S$
$\alpha$ | : Exponent characterizing the divergence of size of precursors on the approach to failure
$d_{\circ}$ | : Average damage level in the specimen
$\delta d(\vec{x})$ | : Perturbations to the mean damage level whose average over the field is zero
$d(\vec{x},t)$ | : Local damage level at any time $t$ during the experiment
$\Delta d(\vec{x},t)$ | : Incremental damage from $d_{\circ}$ at site $\vec{x}$ during the experiment
$Y[\vec{x},d(\vec{x},t),X_{0}]$ | : Local damage driving force determined for a fixed external driving parameter(stress or strain), $X_{0}$.
$\psi(d_{\circ})$ | : Non-local interaction kernel describing the redistribution of elastic energy
$Y_{c}[\vec{x},d(\vec{x},t)]$ | : Local damage resistance corresponding to the threshold for damage growth
$\eta$ | : Damage hardening parameter characterizing the variation of damage resistance with damage level
$y_{c}[\vec{x},d(\vec{x},t)]$ | : Fluctuations in the field of damage resistance visited by the damage like elastic interface
$\delta Y(\vec{x})$ | : Distance to local failure, i.e., incremental damage in line with energy based damage criterion
$\mathcal{K}$ | : Stiffness of springs driving the elastic interface describing the global stability of damage evolution
$v_{\mathrm{m}}$ | : Velocity of the rigid plate that is connected to the pseudo-interface
$\delta d_{\circ}$ | : Incremental damage that occurs whenever the damage criterion is satisfied
$\langle\,\rangle$ | : mean of the enclosed quantity

##  SI Appendix C Experimental precursors from an equivalent force control
scenario

The macroscopic response of the specimen under displacement controlled loading
conditions comprises intermittently occurring load drops during which damage
grows, followed by elastic reloading phases. Here, we construct an equivalent
force control scenario (cyan curve in Fig. 1(a)) from the displacement control
experiment (blue curve). As a result, damage evolution is interpreted as a
micro-instability manifesting as a jump of displacement from
$\Delta_{\mathrm{ini}}$ to $\Delta_{\mathrm{end}}$ at constant force, wherein
distinct damage events (force drops outlined in red) are separated by silent
times corresponding to elastic re-loading. The macroscopic response is thus
reconstructed as a sequence of force plateaus followed by elastic reloading.

###  C.1 Precursors size at the global scale

For the equivalent force control scenario, the incremental work done by the
external force during silent damage periods (for
$\Delta\leq\Delta_{\mathrm{ini}}$ and $\Delta_{\mathrm{end}}\leq\Delta$) is
converted into elastic energy, $dW=dE_{\mathrm{el}}=\frac{1}{2}F\Delta$.
During the equivalent force control cascade (for
$\Delta_{\mathrm{ini}}\leq\Delta\leq\Delta_{\mathrm{end}}$), only a part of
the incremental work $dW=F(\Delta_{\mathrm{end}}-\Delta_{\mathrm{ini}})$ is
stored in the material as elastic energy. The remaining part $\Delta
E_{\mathrm{d}}=\Delta W-\Delta E_{\mathrm{el}}=\Delta W/2$ is dissipated
through damage. This quantity defines the avalanche size $S_{\mathrm{global}}$
at the global (macroscopic) scale.

###  C.2 Precursors duration and damage events within a cascade

The duration of the force plateau is dependent on the cross-head velocity of
the loading machine. Thus, it cannot be used to define the intrinsic duration
$T$ of a precursor. To determine the duration of a cascade of damage events,
we focus on the load drops observed during the force-plateau (Fig. 1(b)). The
intermittent load drops are reminiscent of the individual events within a
cascade of damage growth during a force controlled experiment. Taking $\tau$
as the characteristic duration for these damage events, the duration of a
cascade writes as $T=N\tau$, where $N$ is the number of load drops during the
cascade. Following this framework, we consider the magnitude of the load
drops, denoted by $A$, as the size of the damage events during the cascade.
Only the load drops of magnitude greater than $0.1~{}\mathrm{N}$ are
considered as individual events. As damage events manifest as highly
correlated spatio-temporal clusters, load drops are expected to be macroscopic
manifestation of the clusters of correlated individual damage events observed
at the local scale, see Fig. 2(b) in the main article.

Figure 1:  (a) Schematic of a damage precursor (or avalanche) in an equivalent
force control scenario during the compression experiment in displacement
control. A schematic representation of the localization angle
$\theta_{\mathrm{loc}}$. (b) Global analysis of a precursor as measured in our
experiments: The damage events constituting the avalanche are revealed by
sudden load drops, identified by the red and black circles, denoting the
beginning and the end of the event, respectively. The vertical bars provide
the size $A$ of the damage event that is defined from the energy dissipated
during the event. (c) Binary format of the snapshots of the specimen taken at
different distances from localization. (d) The images are post-treated to
identify the position of the cell centers that are tracked during the
experiment. (e) The displacement of the cell centers with respect to their
initial position are coarse-grained to compute the displacement field
$u_{y}(\vec{x})$ at various distances to localization. (f) Evolution of the
different energies during the experiment as obtained from our global and local
analyses. The work $W$ of the loading machine is converted into elastic energy
$E_{\mathrm{el}}$ (computed both at the local and global scale) and dissipated
energy $E_{\mathrm{d}}$. The inset highlights the intermittent evolution of
the dissipated energy, a feature reminiscent of the avalanche dynamics of the
damage field. (g) Distribution of increments $\rho$ of local dissipation
energy density observed during an avalanche – the distribution is computed
here from the increments extracted from all the avalanches of a single
experiment. It follows an exponential decay
$P(\rho)\propto\mathrm{e}^{-\rho/\rho^{\star}}$, where the exponential cut-off
$\rho^{\star}\simeq 0.18~{}\mathrm{J/m}^{2}$ is used to threshold the
dissipation energy density maps of single avalanches. (h) Auto-correlation
functions $C(\delta r)$ of the thresholded maps of dissipation energy density
$\rho_{th}$ of $10$ precursors of different sizes. The spatial extent $\xi$ of
the precursors is obtained from the fit of the correlation function by an
exponential decay $C(\delta r)\propto\mathrm{-e}^{\delta r/\xi}$. Note that
larger precursors have smaller slopes and hence larger spatial extent.

###  C.3 Precursors size at the local scale

The continuous image acquisition during the experiment allows for the tracking
of both the location and the deformation of individual cells, which are then
used to analyze the spatio-temporal structure of damage cascades at the local
scale. The images are recorded at a rate of 10 frames per second such that
each frame corresponds to a duration comparable to the duration of the
smallest force plateaus. First, we obtain the binary formats of the images
processed using the open source Fiji software [1], Fig. 1(c). The center of
the cells and their relative circularity are then tracked to extract the local
displacement field and local damage field, respectively. The coordinates of
the cell centers recorded at different distances to failure during a typical
experiment is shown in Fig. 1(d). Using a coarse-graining technique [2] based
on mass conservation, we obtain the coarse-grained displacement field
$u_{x}(x,y)~{}\mathrm{and}~{}u_{y}(x,y)$ as well as the damage field $d(x,y)$
of an equivalent disordered continuum. Fig. 1(e) depicts the field
$u_{y}(x,y)$ during a typical experiment. Taking the gradient of the
displacement field, we obtain the strain fields
$\epsilon_{xx}(\vec{x})=\frac{du_{x}}{dx},~{}\epsilon_{yy}(\vec{x})=\frac{du_{y}}{dy}~{}\mathrm{and}~{}\epsilon_{xy}(\vec{x})=\frac{1}{2}\left(\frac{du_{x}}{dy}+\frac{du_{x}}{dy}\right)$.
The Poisson ratio $\nu$ is obtained from a linear fit of
$\epsilon_{xx}(\vec{x})$ vs. $\epsilon_{yy}(\vec{x})$ obtained over time and
space. To track the damage evolution at the local scale, we define the
relative circularity of the cells as our internal damage variable. It turns
out that the damage level does not affect the Poisson ratio ($\nu=0.26$) that
remains nearly constant all along the compression test. However, the Young’s
modulus ($E$), inferred from the global specimen stiffness
$\frac{\sigma_{yy}}{\epsilon_{yy}}=\frac{E}{(1-\nu^{2})}$, decays with the
average damage level $d_{\circ}$. In the following, we consider a damage
dependent elasticity (elasto-damageability) to describe the constitutive
response of the cellular solid. We use the relation
$E\sim\mathcal{O}(d_{\circ}^{2})$ inferred from the average elastic modulus of
the specimen $E$ and the average damage level $d_{\circ}$ to determine the
field $E_{\mathrm{local}}(\vec{x},t)$ of local Young’s moduli from the damage
field $d(\vec{x},t)$ where $\vec{x}\rightarrow(x,y)$. The field of elastic
energy stored in the specimen per unit volume is finally obtained from the
relation,

$E_{\mathrm{el}}^{\mathrm{local}}(\vec{x},t)=\frac{E_{\mathrm{local}}(\vec{x},t)}{2(1-\nu^{2})}\left[\epsilon_{\mathrm{xx}}^{2}(\vec{x},t)+\epsilon_{\mathrm{yy}}^{2}(\vec{x},t)+2\,\nu\epsilon_{\mathrm{xx}}(\vec{x},t)\epsilon_{\mathrm{yy}}(\vec{x},t)+2(1-\nu)\epsilon_{\mathrm{xy}}^{2}(\vec{x},t)\right].$
(8)

$\epsilon_{\mathrm{xx}}(\vec{x},t),\epsilon_{\mathrm{yy}}(\vec{x},t)~{}\mathrm{and}~{}\epsilon_{\mathrm{xy}}(\vec{x},t)$
are the local values of normal strains perpendicular and parallel to the
direction of loading and shear strains, respectively. At the macroscopic
scale, the damage cascades are evidenced from load drops at near constant
displacements. At the local scale, they manifest as a local decrease in stored
elastic energy of the damaging element. In practice, we consider the material
elements for which the elastic energy decreases as damaging elements and
compute the local increment of dissipation energy as $\delta
E_{\mathrm{d}}^{\mathrm{local}}(\vec{x})=-\delta
E_{\mathrm{el}}^{\mathrm{local}}(\vec{x})$. The evolution of the different
terms contributing significantly to the total energy of the isolated system
composed of the specimen and the loading machine, namely the input work, the
elastic energy and the dissipated energy as computed from the global analysis
is given in Fig. 1(f). Our measurement of the total elastic energy from the
local strains and elastic moduli are shown to be in good agreement with this
global analysis. As damage is accompanied by a decrease in local elastic
energy, we obtain the volumetric rate of dissipated energy as

$\rho(\vec{x},t_{i})=b[E_{\mathrm{el}}^{\mathrm{local}}(\vec{x},t_{i})-E_{\mathrm{el}}^{\mathrm{local}}(\vec{x},t_{i+1})].$
(9)

where $b$ is the specimen thickness and the time steps $t_{i}$ correspond to
the displacements
$\Delta(t_{i})\in\left[\Delta_{\mathrm{ini}},\Delta_{\mathrm{end}}\right]$.
Combining all the increments of dissipated energy belonging to the same damage
cascade provides the field $\rho(\vec{x})=\sum_{i}\rho(\vec{x},t_{i})$ of
energy dissipated during a given precursor. At the specimen level, we obtain
the precursor size byintegration over the whole specimen

$S_{\mathrm{local}}=\iint\rho(x,y)\,\mathrm{d}x\mathrm{d}y$ (10)

The distribution of $\rho(\vec{x})$ obtained from all cascades follows an
exponential decay $P(\rho)\propto\mathrm{e}^{-\rho/\rho^{\ast}}$ with
$\rho^{\ast}\simeq 0.18~{}\mathrm{J/m^{2}}$ (Fig. 1(g)). The value of
$\rho^{\ast}$ is used as a threshold to obtain binary formats of the
dissipation density maps $\rho_{\mathrm{th}}(\vec{x})$ of precursors (inset of
Fig. 2(a) in the main article).

###  C.4 Spatial extent of the precursors

To determine the spatial extent $\xi$ of the precursors, we compute the 2D
auto-correlation of the thresholded dissipation density maps
$\rho_{\mathrm{th}}(\vec{x})$. The correlation function $C(\delta
r)=\langle\rho_{\mathrm{th}}(\vec{x}).\rho_{\mathrm{th}}(\vec{x}+\delta\vec{x})\rangle_{\vec{x},|\delta\vec{x}|=\delta
r}/\langle\rho_{\mathrm{th}}(\vec{x})^{2}\rangle_{\vec{x}}$ is shown in Fig.
1(h) for precursors of different sizes $S_{\mathrm{local}}$. The function
decreases faster for smaller avalanche size, pointing out a smaller spatial
extent. The correlation function can be fitted by an exponential decay
$C(\delta r)\propto\mathrm{e}^{-\delta r/\xi}$, defining the correlation
length $\xi$. This quantity is plotted as a function of the precursors’ size
and duration in Fig. 2(c) and (d) of the main article.

##  SI Appendix D Determination of field of the damage resistance

###  D.1 Retrospective determination from the local damage driving force

One of the most appealing feature of our experiments is the full-field
measurement of the damage field and its evolution over time. This can be
harnessed to access the field of local damage resistance
$Y_{\mathrm{c}}(\vec{x})$, a quantity that has been hardly measured in the
literature, even at the specimen scale. The proposed methodology consists in
using the energy based damage criterion

$\left\\{\begin{array}[]{ll}Y[\vec{x},d(\vec{x},t)]&<Y_{c}[\vec{x},d(\vec{x},t)]~{}~{}~{}~{}\Rightarrow~{}~{}\delta
d(\vec{x},t)=0\\\
Y[\vec{x},d(\vec{x},t)]&=Y_{c}[\vec{x},d(\vec{x},t)]~{}~{}~{}~{}\Rightarrow~{}~{}\delta
d(\vec{x},t)>0\end{array}\right.$ (11)

that, as shown in the following, provides retrospectively the local value of
the damage resistance $Y_{c}[\vec{x},d(\vec{x},t)]$ each time a damage event
takes place. First, we compute the field of damage driving force by
differentiating Eq. 8 with respect to the damage parameter

$Y(\vec{x},t)=-\frac{E^{\prime}_{\mathrm{local}}(\vec{x},t)}{2(1-\nu^{2})}\left[\epsilon_{xx}^{2}(\vec{\mathrm{x}},t)+\epsilon_{yy}^{2}(\vec{\mathrm{x}},t)+2\,\nu\epsilon_{xx}(\vec{\mathrm{x}},t)\epsilon_{yy}(\vec{\mathrm{x}},t)+2(1-\nu)\epsilon_{xy}^{2}(\vec{\mathrm{x}},t)\right].$
(12)

where the derivative to the local elastic modulus
$E^{\prime}_{\mathrm{local}}\simeq\frac{dE}{dd_{\circ}}(d(\vec{x},t))$. Here,
we assume that the variations of the average macroscopic elastic modulus with
the average damage level $d_{\circ}$ hold at the local scale too. As
schematically illustrated in Fig. 2(a), the value of the local damage driving
force $Y(\vec{x},t)$ is assigned to the damage resistance
$Y_{\mathrm{c}}(\vec{x},d(\vec{x},t))$ each time a damage event takes place.
The damage driving force $Y(\vec{x},t_{i})$ at the onset of damage growth
$t_{i}$ is also assigned to the damage resistance
$Y_{\mathrm{c}}(\vec{x},d(\vec{x},\underset{\widetilde{}}{t}))=Y(\vec{x},t_{i})$
during the whole sequence of elastic reloading preceding $t_{i}$, i.e. for
$\underset{\widetilde{}}{t}\in[t_{i-k},t_{i-1}...t_{i}]$ where $t_{i-k}$ is
the time at which the previous damage event in the same material element
$\vec{x}$ took place. Now considering the average damage resistance
$Y_{\mathrm{c\circ}}(d_{\circ})=\langle
Y_{\mathrm{c}}(\vec{x},d(\vec{x},t))\rangle_{d(\vec{x},t)=d_{\circ}}$ at some
given damage level $d_{\circ}$, we see in Fig. 2(b) that, after an initial
transient regime during which damage evolution is dominated by randomly
distributed small localized events, $Y_{\mathrm{c\circ}}(d_{\circ})$ increases
with damage. This hardening can be described by a linear relation
$Y_{\mathrm{c\circ}}(d_{\circ})=Y_{\mathrm{c}}^{\circ}(1+\eta\,d_{\circ})$
where $Y_{\mathrm{c}}^{\circ}\simeq 1.4~{}\mathrm{kJ/m}^{3}$ and $\eta\simeq
44$.

Figure 2: (a) Schematic of the retrospective assignment procedure used to
determine the local damage resistance $Y_{\mathrm{c}}[\vec{x},d(\vec{x},t)]$
from the field of damage driving force $Y[\vec{x},d(\vec{x},t),t]$. (b)
Variations of the average damage resistance $Y_{\mathrm{c\circ}}$ with respect
to the damage level $d_{\circ}$. The increase of the damage resistance can be
described by a linear hardening
$Y_{\mathrm{c\circ}}(d_{\circ})=Y_{\mathrm{c}}^{\circ}(1+\eta\,d_{\circ})$
with a hardening coefficient $\eta\simeq 44$. (c) Variations of the average
damage resistance $Y_{c\circ}$ with the average damage level $d_{\circ}$
inferred from the variations with respect to $d-\circ$ of the accumulated
dissipated energy at the specimen level that can be described (see inset) by a
quadratic function (in red). (c) Comparison of the precursor size
$S_{\mathrm{local}}$ obtained our two independent methods: from the variations
of elastic energy in the abscissa and from the local increase of the damage
variable and the damage resistance $Y_{c\circ}(d_{\circ})$ in the ordinates.
(e) Dissipated energy density maps of a typical precursor as obtained from the
two methods.

###  D.2 Alternative methodology: determination from the accumulated
dissipated energy

As an alternative to the above retrospective determination of $Y_{c}$ at the
local scale, the average value of damage resistance can be estimated from the
variation of accumulated dissipated energy using

$Y_{c\circ}(d_{\circ})=\frac{1}{V_{\circ}}\frac{dE_{d}}{dd_{\circ}}$ (13)

where $V_{\circ}$ is the initial volume of the cellular solid and $E_{d}$ is
the accumulated dissipated energy obtained from the global balance of energy
(see inset of Fig. 2(c)). Assuming a linear variation $Y_{c\circ}\propto\eta
d_{\circ}$, the accumulated dissipated energy is fitted by a quadratic
function providing
$Y_{c}^{\circ}=1.1\mathrm{kJ/m^{3}~{}\mathrm{and}~{}\eta\simeq 45}$, both
values that are close to the one inferred from the previous method. Using this
linear relation, we then determine the field of damage resistance
$Y_{c}(\vec{x},t)=Y_{c\circ}(d(\vec{x},t))$ from the damage field. The energy
dissipated at each time step is then computed as
$\rho(\vec{x},t_{i})=\mathrm{max}(Y_{c}(\vec{x},t_{i})\delta
d(\vec{x},t_{i}),0)$, where $\delta d(\vec{x},t_{i})$ is the incremental
damage growth in $\vec{x}$. The dissipated energy density $\rho(\vec{x})$
during one precursor and the precursor size $S_{local}$ follow from (10). A
comparison with the first method is presented in Fig. 2(d) that shows a good
agreement. The maps of local dissipated energy density of a typical cascade as
obtained from the two approaches are compared in Fig. 2(e). The discrepancies
in the maps of the dissipation energy density are localized at the periphery
of the specimen where the rather large damage may be attributed to frictional
effects with the walls. Thus, we observe a good agreement between both
methods, thus validating the assumption of local energy balance made
previously.

##  SI Appendix E Theoretical modeling of the evolution of the damage field

We now detail the damage model used to describe our experiments. It intends to
account for the co-action of material disorder and long-range elastic
interactions in the aftermath of a damage event, in the spirit of the physics-
based non-local damage models proposed in Berthier et al. [3] and Dansereau et
al. [4]. In particular, we provide the redistribution kernel characterizing
the long-range interactions derived here for the particular case of a 2D
specimen compressed between two fixed lateral walls. The first step, detailed
in sections E.1 and E.2 for force and displacement imposed conditions,
respectively consists in calculating the distribution of damage driving force
resulting from a heterogeneous distribution of damage (see Eqs. (22) and (28)
that we write in a generic form valid for both loading conditions in (30)). We
then extend our calculations to the case of a monotonically increasing loading
amplitude in section E.3. In section E.4 we derive the damage evolution
equation from the damage criterion of (11). We then explain in section E.5 how
the evolution equation of the damage field belongs to the theoretical
framework of driven disordered elastic interfaces. The stability of the damage
spreading process is studied in section E.6. The localization threshold is
shown to coincide with peak load. Lastly, we show in section E.7 that the rate
of damage growth diverges on approaching localization.

We start by considering that the level of damage in each material element
located in $\vec{x}$ is $d(\vec{x})\geq 0$, where $d=0$ corresponds to the
intact initial material. The level of damage is assumed to affect the elastic
modulus $E[d(\vec{x})]$ of the cellular solid. Its Poisson’s ratio $\nu$ is
assumed to be independent of the damage level, as supported by our
experimental observations (Sec. S1C). The constitutive behavior of the
cellular solid is then the one of an elasto-damageable solid under plane
stress conditions

$\centering\bm{\epsilon}=\frac{(1+\nu)\bm{\sigma}-\nu\mathrm{Tr}(\bm{\sigma})}{E}.\@add@centering$
(14)

Due to the lateral confinement $(\epsilon_{xx}=0)$, the in-plane stress
components follow: $\sigma_{yy}=\epsilon_{yy}E/(1-\nu^{2})$ along the loading
direction and $\sigma_{xx}=\nu\sigma_{yy}$ along the lateral direction. In
absence of out-of-plane stress ($\sigma_{zz}=0$), Hooke’s law predicts
$\epsilon_{zz}=-\epsilon_{yy}\nu/(1-\nu)$. The elastic energy per unit volume
of the material,
$w(\bm{\epsilon},\bm{\sigma})=(1/2)~{}\bm{\epsilon}\colon\bm{\sigma}$, can
then be expressed as

$\centering
w(\bm{\epsilon},\bm{\sigma})=\frac{(1+\nu)\mathrm{Tr}(\bm{\sigma}^{2})-\nu{\mathrm{Tr}(\bm{\sigma})}^{2}}{2E}\@add@centering$
(15)

###  E.1 Damage driving force under force control conditions

For a homogeneous damage level $d_{\circ}$, the elastic energy release rate
that drives damage growth follows

$Y(d_{\circ},\bm{\sigma})=\left.\frac{\partial w}{\partial
d_{\circ}}\right|_{\bm{\sigma}}=\frac{\overline{\mu}^{\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma}^{2})}{4}+\frac{{\nu
E^{\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma})}^{2}}{2E(d_{\circ})^{2}}$ (16)

where $\mu=E/2(1+\nu)$ is the Lamé constant and $\overline{\mu}=1/\mu$. The
prime denotes the derivative with respect to the damage variable $d_{\circ}$.
To examine the effect of material disorder on the damage field and in
particular on the stress redistributions following individual damage events,
we introduce weak variations in the damage field $d(\vec{x})=d_{\circ}+\delta
d(\vec{x})$ where $d_{\circ}$ is the average damage level and $\delta
d(\vec{x})\ll d_{\circ}$ are its perturbations. The heterogeneities in the
damage field result from the spatial variations of the elastic modulus.
Consequently, the stress field is heterogeneous and writes as
$\bm{\sigma}(\vec{x})=\bm{\sigma_{\circ}}+\delta\bm{\sigma}(\vec{x})$ where

$\bm{\sigma_{\circ}}=\sigma_{\mathrm{ext}}\begin{bmatrix}\nu&0&0\\\ 0&1&0\\\
0&0&0\end{bmatrix}~{}\mathrm{and}~{}\delta\bm{\sigma}(\vec{x})=\begin{bmatrix}\delta\sigma_{xx}(\vec{x})&\delta\sigma_{xy}(\vec{x})&0\\\
\delta\sigma_{yx}(\vec{x})&\delta\sigma_{yy}(\vec{x})&0\\\
0&0&0\end{bmatrix}.$ (17)

Here, $\sigma_{\mathrm{ext}}=F/bL$ is the average stress imposed by the
loading machine to the upper surface $bL$ of the specimen. As a consequence,
the damage driving force provided in (16) is also heterogeneous and can be
decomposed as $Y[d(\vec{x})]=Y_{\circ}(d_{\circ})+\delta Y[d(\vec{x})]$ where

$\left\\{\begin{array}[]{ll}Y_{\circ}(d_{\circ},\bm{\sigma_{\circ}})&=\frac{\overline{\mu}^{\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma_{\circ}}^{2})}{4}+\frac{{\nu
E^{\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma_{\circ}})}^{2}}{2E(d_{\circ})^{2}}\\\
\delta
Y[d(\vec{x}),\bm{\sigma}_{\circ}]&=\Bigg{[}\frac{\overline{\mu}^{\prime\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma_{\circ}}^{2})}{4}-\frac{{\nu\mathrm{Tr}(\bm{\sigma_{\circ}})}^{2}}{2}\left(\frac{2E^{\prime}(d_{\circ})^{2}}{E(d_{\circ})^{3}}-\frac{E^{\prime\prime}(d_{\circ})}{E(d_{\circ})^{2}}\right)\Bigg{]}\delta
d(\vec{x})+\Bigg{[}\frac{\overline{\mu}^{\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma_{\circ}}\delta\bm{\sigma}(\vec{x}))}{2}+\frac{\nu
E^{\prime}(d_{\circ})\mathrm{Tr}(\bm{\sigma_{\circ}})\mathrm{Tr}(\delta\bm{\sigma}(\vec{x}))}{E(d_{\circ})^{2}}\Bigg{]}.\end{array}\right.$
(18)

The first term in the expression of $\delta Y[d(\vec{x})]$ is a local term,
namely its value in $\vec{x}$ depends only on the damage level in $\vec{x}$.
It can be obtained by considering a homogeneous damage level $d_{\circ}$ and
varying $Y[d_{\circ}=d(\vec{x})]$ with respect to it as follows

$\centering\left.\frac{\partial Y_{\circ}}{\partial
d_{\circ}}\right|_{\sigma_{\circ}}=\frac{(1-\nu^{2})}{2}\sigma_{ext}^{2}\left[\frac{2E^{\prime}(d_{\circ})^{2}}{E(d_{\circ})^{3}}-\frac{E^{\prime\prime}(d_{\circ})}{E(d_{\circ})^{2}}\right].\@add@centering$
(19)

To calculate the non-local contributions to $\delta Y[d(\vec{x})]$ (i.e., the
second term in (18)) we first need to evaluate the stress perturbations
$\delta\bm{\sigma}(\vec{x})$ emerging from the spatial variations of elastic
modulus, as performed in Dansereau et al. [4]. We then use a perturbative
approach in Fourier space and such calculations that gives

$\delta\bm{\tilde{\sigma}}(q_{x},q_{y})=\frac{\overline{\mu}^{\prime}(d_{\circ})}{\overline{\mu}(d_{\circ})}\left(\mathbf{O}\cdot\bm{\sigma_{\circ}}\cdot\mathbf{O}-\nu\left[\left(1-\mathbf{O}\right)\colon\bm{\sigma_{\circ}}\right]\mathbf{O}\right)\delta\tilde{d}(q_{x},q_{y}),$
(20)

where $\mathbf{O}$ stands for the Oseen tensor, $\mathbf{O}=1-\mathbf{Q}$ with
the tensor $\mathbf{Q_{ij}}=q_{i}q_{j}/|\vec{q}|^{2}$.
$\delta\bm{\tilde{\sigma}}(q_{x},q_{y})$ and $\delta\tilde{d}(q_{x},q_{y})$
are the 2D Fourier transform of the stress field perturbations and the damage
field perturbations, respectively. Substituting (20) into (18), we obtain the
non-local contributions to the damage driving force that writes as
$\tilde{\psi}\delta\tilde{d}$ in Fourier space, that defines the interaction
kernel $\tilde{\psi}(q_{x},q_{y})$. In real space and in polar coordinates, it