Andrea Cannizzo
email: andrea.cannizzo@iemn.fr

Stefano Giordano
email: stefano.giordano@univ-lille.fr

Thermal effects on fracture and brittle-to-ductile transition

The fracture behavior of brittle and ductile materials can be strongly influenced by thermal fluctuations, especially in micro-and nano-devices as well as in rubberlike and biological materials. However, temperature effects, in particular on the brittle-to-ductile transition, still require a deeper theoretical investigation. As a step in this direction we propose a theory, based on equilibrium statistical mechanics, able to describe the temperature dependent brittle fracture and brittle-to-ductile transition in prototypical discrete systems consisting in a lattice with breakable elements. Concerning the brittle behavior, we obtain closed form expressions for the temperature-dependent fracture stress and strain, representing a generalized Griffith criterion, ultimately describing the fracture as a genuine phase transition. With regard to the brittle-to-ductile transition, we obtain a complex critical scenario characterized by a threshold temperature between the two fracture regimes (brittle and ductile), an upper and a lower yield strength, and a critical temperature corresponding to the complete breakdown. To show the effectiveness of the proposed models in describing thermal fracture behaviors at small scales, we successfully compare our theoretical results with molecular dynamics simulations of Si and GaN nanowires.

I. INTRODUCTION

The mechanical degradation of a material typically results from the insurgence of cracks, from their geometric arrangement and interactions and, finally, from temperature. The classical Griffith energetic approach in fracture mechanics deduces that, under a homogeneous stress σ , a single slit crack with half-length L grows if σ > 2γ s E ′ /(πL), while if σ < 2γ s E ′ /(πL), it remains stable [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF][START_REF] Pugno | The centenary of Griffith's theory[END_REF]. Here, E ′ is the equivalent elastic modulus, equal to the Young modulus E in plane stress condition, and equal to E/(1ν 2 ) in plane strain condition, where ν is the Poisson ratio of the material. Moreover, γ s is the surface energy density, that is, the energy expended to debond a unit length crack. This stability criterion measures an energy competition between the free surface energy created by the fracture and the elastic energy stored in the deformable solid. The latter has been evaluated within the linear elasticity theory by Inglis [START_REF] Inglis | Stresses in a plate due to the presence of cracks and sharp corners[END_REF] and Kolosoff [START_REF] Kolosoff | Uber einige Eigenschaften des ebenen Problems der Elastizitatstheorie[END_REF] and used by Griffith to develop his criterion. The ingenious approach proposed by Griffith in his celebrated criterion has been largely and successfully tested in glass and other brittle materials containing cracks of controlled length [START_REF] Lawn | Fracture of Brittle Solids[END_REF][START_REF] Öchsner | Continuum Damage and Fracture Mechanics[END_REF], and also validated by atomistic simulations in ideal monocrystalline systems [START_REF] Mattoni | Atomic Scale Origin of Crack Resistance in Brittle Fracture[END_REF][START_REF] Bao | Molecular Dynamics Simulation of Nanocrack Propagation in Single-Layer MoS 2 Nanosheets[END_REF]. Its main limiting hypothesis is that the overall fracture energy coincides with the surface energy, i.e. with the energy needed to break the bonds between the two crack faces. Since Griffith's theory fails to apply to ductile materials (where fracture energy is much higher than the only surface energy [START_REF] Orowan | Fracture and strength of solids[END_REF][START_REF] Orowan | Energy Criteria of Fracture[END_REF]) it was generalized by Irwin to include plastic dissipation [START_REF] Irwin | Analysis of stresses and strains near the end of a crack traversing a plate[END_REF][START_REF] Irwin | Fracture[END_REF]. More advanced models for ductile fracture, taking into account an explicit description of the 'cohesive zone' where the plastic processes localize, have then been proposed by Dugdale and Barenblatt [START_REF] Dugdale | Yielding of steel sheets containing slits[END_REF][START_REF] Barenblatt | The mathematical theory of equilibrium cracks in brittle fracture[END_REF]. The theoretical relationship between the Griffith's theory, its Irwin modification, and the Dugdale-Barenblatt models was initially studied by Willis [START_REF] Willis | A comparison of the fracture criteria of Griffith and Barenblatt[END_REF] and further investigated by Rice through the concept of J-integral [START_REF] Rice | A path independent integral and the approximate analysis of strain concentration by notches and cracks[END_REF]. In addition to the single fracture study, an important topic in linear elastic fracture mechanics is represented by the collective degradation mechanism induced by populations of cracks that interact depending on their geometric arrangement [START_REF] Kachanov | Effective elastic properties of cracked solids: critical review of some basic concepts[END_REF][START_REF] Kachanov | Elastic Solids with Many Cracks and Related Problems[END_REF][START_REF] Giordano | Effects of the orientational distribution of cracks in solids[END_REF][START_REF] Giordano | Effects of the orientational distribution of cracks in isotropic solids[END_REF][START_REF] Giordano | Local elastic fields around cracks and their stress density of states[END_REF][START_REF] Giordano | Elastic properties of solids containing elliptic cracks[END_REF][START_REF] Dormieux | Stress-based estimates and bounds of effective elastic properties: The case of cracked media with unilateral effects[END_REF][START_REF] Giordano | Brittle fracture: from elasticity theory to atomistic simulations[END_REF][START_REF] Manca | Scaling Shift in Multicracked Fiber Bundles[END_REF][START_REF] Manca | Stochastic mechanical degradation of multi-cracked fiber bundles with elastic and viscous interactions[END_REF][START_REF] Dormieux | Micromechanics of Fracture and Damage[END_REF]. Due to the wide scientific and technological interest, the theoretical and experimental studies of fracture phenomena have been extremely extensive. Therefore, its history is long and complicated and, here, we refer the reader to the relevant literature [START_REF] Broberg | Critical review of some theories in fracture mechanics[END_REF][START_REF] Cherepanov | Fractal fracture mechanics -A review[END_REF][START_REF] Broberg | Cracks and Fracture[END_REF][START_REF] Erdogan | Fracture mechanics[END_REF][START_REF] Gross | Some Remarks on the History of Fracture Mechanics[END_REF]. We simply mention that current advanced researches concern the traction-separation relation in cohesive models [START_REF] Park | Cohesive Zone Models: A Critical Review of Traction-Separation Relationships Across Fracture Surfaces[END_REF][START_REF] Ferdousi | Characterize traction-separation relation and interfacial imperfections by data-driven machine learning models[END_REF], the instability in dynamic fracture [START_REF] Fineberg | Instability in dynamic fracture[END_REF][START_REF] Chen | Instability in dynamic fracture and the failure of the classical theory of cracks[END_REF], and the variational approach to fracture [START_REF] Bourdin | The Variational Approach to Fracture[END_REF][START_REF] Mesgarnejad | Validation simulations for the variational approach to fracture[END_REF]. Computational techniques for cracks propagation include the phase field method and the dual-horizon peridynamics formulation [START_REF] Miehe | A phase field model for rate-independent crack propagation: Robust algorithmic implementation based on operator splits[END_REF][START_REF] Ren | An explicit phase field method for brittle dynamic fracture[END_REF][START_REF] Agwai | Predicting crack propagation with peridynamics: a comparative study[END_REF][START_REF] Ren | Dual-horizon peridynamics[END_REF].

Previously described investigations are predominantly based on deterministic assumptions and theories. Of course, also statistical approaches have been widely applied to rupture phenomena [START_REF] Truskinovsky | Fracture as a phase transition[END_REF][START_REF] Zapperi | First-Order Transition in the Breakdown of Disordered Media[END_REF][START_REF] Moreno | Fracture and Second-Order Phase Transitions[END_REF][START_REF] Alava | Statistical models of fracture[END_REF][START_REF] Kawamura | Statistical physics of fracture, friction, and earthquakes[END_REF][START_REF] Gjerden | Universality Classes in Constrained Crack Growth[END_REF], and among others those based on the so-called fiber bundle model are particularly significant [START_REF] Sornette | Elasticity and failure of a set of elements loaded in parallel[END_REF][START_REF] Harlow | Approximations for the strength distribution and size effect in an idealized lattice model of material breakdown[END_REF][START_REF] Pradhan | Failure processes in elastic fiber bundles[END_REF][START_REF] Kovács | Brittle-to-ductile transition in a fiber bundle with strong heterogeneity[END_REF]. Importantly, the statistical analysis plays a crucial role for understanding the effect of disorder in failure processes [START_REF]Disorder and fracture[END_REF][START_REF] Pérez-Reche | Driving-Induced Crossover: From Classical Criticality to Self-Organized Criticality[END_REF][START_REF] Parisi | Shear bands as manifestation of a criticality in yielding amorphous solids[END_REF][START_REF] Ponson | Crack propagation through disordered materials as a depinning transition: A critical test of the theory[END_REF][START_REF] Niccolini | Dimensional transitions in creeping materials due to nonlinearity and microstructural disorder[END_REF][START_REF] Sinha | Phase transitions and correlations in fracture processes where disorder and stress compete[END_REF][START_REF] Lebihain | Effective toughness of disordered brittle solids: A homogenization framework[END_REF][START_REF] Borja Da Rocha | Mean field fracture in disordered solids: Statistics of fluctuations[END_REF]. Despite the wide diffusion of statistical techniques, the approaches that allow to study the effects of thermal fluctuations on the fracture are rather limited [START_REF] Ausloos | Phase transition theory approach to fracture of materials[END_REF][START_REF] Khitrin | Brittle fracture as a phase transition[END_REF][START_REF] Marder | Statistical mechanics of cracks[END_REF][START_REF] Ortiz | Plastic Yielding as a Phase Transition[END_REF][START_REF] Santucci | Thermal activation of rupture and slow crack growth in a model of homogeneous brittle materials[END_REF][START_REF] Guarino | The cooperative effect of load and disorder in thermally activated rupture of a two-dimensional random fuse network[END_REF][START_REF] Guarino | Thermally activated fracture of porous media[END_REF][START_REF] Vincent-Dospital | How heat controls fracture: the thermodynamics of creeping and avalanching cracks[END_REF][START_REF] Vincent-Dospital | Thermal weakening of cracks and brittleductile transition of matter: A phase model[END_REF]. In particular, the temperature dependence of crack stability criteria has not been studied explicitly. For this reason, we propose here two paradigmatic models able to evaluate the effects of thermal fluctuations on the quasi-static brittle fracture and on the brittle-to-ductile transition. These approaches make it possible to study how fracture stability is influenced by temperature changes and to determine the transition temperature between brittle and ductile behaviors. The proposed fracture models are based on equilibrium statistical mechanics and they are implemented by means of the spin variable approach, useful to deal with arbitrarily non-convex potential energies [START_REF] Giordano | Spin variable approach for the statistical mechanics of folding and unfolding chains[END_REF]. This method has been largely applied to several situations including the physics of muscles [START_REF] Caruel | Statistical mechanics of the Huxley-Simmons model[END_REF][START_REF] Caruel | Physics of muscle contraction[END_REF], the folding of macromolecules [START_REF] Benedito | Thermodynamics of small systems with conformational transitions: the case of two-state freely jointed chains with extensible units[END_REF][START_REF] Benedito | Isotensional and isometric forceextension response of chains with bistable units and Ising interactions[END_REF][START_REF] Florio | Unveiling the influence of device stiffness in single macromolecule unfolding[END_REF][START_REF] Bellino | The influence of device handles in single-molecule experiments[END_REF][START_REF] Marco | Thermodynamics of Extra-Toughness and Hidden-Length in Polymeric Materials with Sacrificial Bonds[END_REF], the adhesion or peeling processes [START_REF] Florio | Role of temperature in the decohesion of an elastic chain tethered to a substrate by onsite breakable links[END_REF][START_REF] Cannizzo | Temperature controlled decohesion regimes of an elastic chain adhering to a fixed substrate by softening and breakable bonds[END_REF], the phase transformations in solids [START_REF] Bellino | On the competition between interface energy and temperature in phase transition phenomena[END_REF][START_REF] Cannizzo | Thermal control of nucleation and propagation transition stresses in discrete lattices with non-local interactions and non-convex energy[END_REF], and the stick-slip on rigid substrate [START_REF] Giordano | Statistical mechanics of rate-independent stickslip on a corrugated surface composed of parabolic wells[END_REF]. This technique complements the more classical methods used to study the behaviour of physical systems with multiple stable and metastable states [START_REF] Prados | Sawtooth patterns in force-extension curves of biomolecules: An equilibriumstatistical-mechanics theory[END_REF][START_REF] Bonilla | Theory of forceextension curves for modular proteins and DNA hairpins[END_REF][START_REF] Tommasi | An energetic model for macromolecules unfolding in stretching experiments[END_REF][START_REF] Benichou | Structures undergoing discrete phase transformation[END_REF]. In the context of fracture, the prototypical models here proposed are discrete and based on quasi-one-dimensional lattices composed of breakable and unbreakable bonds. While the unbreakable springs serve to distribute the forces in the system, thus describing material elastic energy, the breakable springs are useful in mimicking the fracture propagation. It is important to note in this context that the role of discreteness in fracture models has already been highlighted in different studies [START_REF] Pugno | Quantized fracture mechanics[END_REF][START_REF] Ippolito | Role of lattice discreteness on brittle fracture: Atomistic simulations versus analytical models[END_REF][START_REF] Ippolito | Failure strength of brittle materials containing nanovoids[END_REF][START_REF] Pugno | Atomistic fracture: QFM vs[END_REF].

In the first model proposed here (elasto-fragile model), developed to describe temperature effects in brittle fracture, each breakable spring can be in one of two states, intact or broken, depending on its extension. Conversely, in the second proposed model (softening-fracture model), useful to describe temperature effects in brittle-to-ductile transitions, each breakable spring can be in one of three states, intact, softened or broken, depending again from the spring extension. In our model, the softened state represents an intermediate configuration where the elastic constant of the spring is smaller than that of the intact spring, but it is still not zero as instead assumed for the broken configuration. This intermediate state represents here the counterpart of the material behavior of the cohesive zone introduced in the Dugdale-Barenblatt model, having also the role of introducing an internal length scale. Both proposed models are approached by calculating the exact partition function, by an approximation obtained for large values of the number N of breakable springs, and, finally by the analysis of the thermodynamic limit. This multifaceted treatment allows us to state that both models exhibit a critical behavior with an associated phase transition, whose mechanical implications are thoroughly discussed. In our opinion the results are particularly useful for the interpretation of failure processes in micro-and nano-systems, where the effect of temperature is typically studied experimentally and through molecular dynamics simulations [START_REF] Dijkstra | Temperature-dependent fracture mechanisms in ultra-high strength polyethylene fibers[END_REF][START_REF] Kang | Brittle and ductile fracture of semiconductor nanowires -molecular dynamics simulations[END_REF][START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF][START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF][START_REF] Wu | Nanowire Failure: Long = Brittle and Short = Ductile[END_REF][START_REF] Wang | Near-ideal theoretical strength in gold nanowires containing angstrom scale twins[END_REF][START_REF] Ast | The brittle-ductile transition of tungsten single crystals at the micro-scale[END_REF][START_REF] Munshi | Insights into the mechanical properties and fracture mechanism of Cadmium Telluride nanowire[END_REF][START_REF] Cheng | In Situ Nano-thermomechanical Experiment Reveals Brittle to Ductile Transition in Silicon Nanowires[END_REF][START_REF] Puglisi | Multi-scale modelling of rubber-like materials and soft tissues: an appraisal[END_REF][START_REF] Florio | A coarsegrained mechanical model for folding and unfolding of tropoelastin with possible mutations[END_REF][START_REF] Florio | A predictive model for the thermomechanical melting transition of double stranded DNA[END_REF]. We want to emphasize that our models, being discrete and addressed to the study of thermal fluctuations, neglect important aspects related to the distribution of elastic fields around the fracture. This is consistent with the fact that they are not developed to replace classical models of continuum mechanics but rather to provide new elements to improve and complement them.

The paper is structured as follows. In Section II, we introduce the first model for the brittle fracture and we apply the tools of statistical mechanics to eventually obtain exact results. Then, in Section III, we obtain an approximate analytic solution for systems with a large number N of breakable springs, and in Section IV, we study the thermodynamic . Moreover, the first η nodes (i = 1, . . . , η) are also linked to the bottom layer (at y = 0) through η vertical breakable springs of elastic constant h (intact in blue or dark gray; broken in orange or intermediate gray). We underline that the first node (i = 0) and the last one (i = N +1) are anchored to the bottom and the top layer, respectively. Hence, the first and the last shear springs fix the direction of the crack propagation from the right to the left of the system. limit with N → ∞. This allows for a generalization of Griffith's criterion that takes temperature into account by means of a critical behavior. Concerning the model with the softening mechanism, we introduce its structure and we elaborate its formalism in Section V. Further, we obtain its asymptotic behavior for large values of N in Section VI, and we study the thermodynamic limit with N → ∞ in Section VII. Here we obtain the closed form expression for the brittle-to-ductile transition temperature and describe the corresponding complex critical scenario.

II. ELASTO-FRAGILE MODEL

We introduce here a discrete model that helps us to better understand the effect of thermal fluctuations on the brittle fracture processes in solid materials. As shown in Fig. 1, the model consists in a network of springs with different elastic constants arranged in a quasi-one-dimensional lattice, aimed at reproducing a mode I fracture geometry [START_REF] Giordano | Brittle fracture: from elasticity theory to atomistic simulations[END_REF]. Based on symmetry assumptions, we reduce the scheme in Fig. 1(a) to the one in Fig. 1(b). The structure in Fig. 1(b) connects a fixed substrate, at y = 0, with a rigid top layer that can be placed at different heights y = Y (isometric conditions within the Helmholtz ensemble). More in detail, see Fig. 1(b), this structure is composed of a series of N + 1 springs, with elastic constant k, linked together to form an horizontal chain (colored in black). The left end-side of the chain, at i = 0, is attached to the bottom fixed substrate at y = 0 while the other end, located at i = N + 1, is attached to the top layer at y = Y . The inner points of the chain, identified by i = 1, . . . , N, are individually linked to the top layer through N vertical springs with elastic constant l (colored in yellow or light gray), mimicking the elasticity of the upper half plane. Moreover, the first ones (i = 1, . . . , η), are linked to the bottom layer through η breakable springs (intact in blue or dark gray; broken in orange or intermediate gray). We assume that it exists an elongation threshold Y M after which the potential energy for a breakable spring is constant and the resulting elastic force is zero (broken state), see Fig. 2. The behavior of a single breakable spring corresponds to an elastic constant h when its elongation y i does not exceed the threshold Y M (see Fig. 2). When |y i | > Y M the potential energy is constant and then the resulting force is zero (see Fig. 2).

We remark that we adopted here simple piecewise linear constitutive equations for the springs of the system. We wish to mention that several important results have been obtained for the crack propagation in nonlinear materials described by power-law stress-strain behavior [START_REF] Hutchinson | Singular behaviour at the end of a tensile crack in a hardening material[END_REF][START_REF] Rice | Plane strain deformation near a crack tip in a power-law hardening material[END_REF][START_REF] Carpinteri | Asymptotic analysis in linear elasticity: from the pioneering studies by Wieghardt and Irwin until today[END_REF]. However, from one side it is difficult to combine the thermal analysis with nonlinear materials, and from the other side the nonlinear phenomena are rather limited in nanoscopic systems [START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF][START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF]. It is also worth mentioning that due to the discrete quasi-onedimensional structure of our model, it is not possible to find here the results concerning the stress singularities at the crack tip and the calculation of the corresponding stress intensity factor (this is true for both the elasto-fragile model and the softening-fracture model) [START_REF] Giordano | Local elastic fields around cracks and their stress density of states[END_REF][START_REF] Giordano | Brittle fracture: from elasticity theory to atomistic simulations[END_REF].

We analyze the fracture behavior of the proposed model in the framework of equilibrium statistical mechanics, introducing a temperature T of an embedding thermal bath. As previously anticipated, we investigate its behavior by adopting isometric conditions corresponding to the Helmholtz ensemble [START_REF] Manca | Elasticity of flexible and semiflexible polymers with extensible bonds in the Gibbs and Helmholtz ensembles[END_REF][START_REF] Gibbs | Elementary Principles in Statistical Mechanics[END_REF][START_REF] Weiner | Statistical Mechanics of Elasticity[END_REF]. We make the assumption that, during the system extension (i.e., increasing Y ), the system is composed of a segment with ξ intact elements on the left side of the system, and of a segment with ηξ broken elements on the right. As a result, the system evolution is characterized by the propagation of a single interface between intact and broken springs, regulated by the assigned traction conditions and by the temperature. This hypothesis (known as single domain wall assumption) simplifies the calculations and makes it possible to analytically derive the partition function and, thus, the important macroscopic physical quantities. The same hypothesis is considered in the classical continuum fracture models recalled in the introduction. The configurations previously described can be summarized by the relation 1 ≤ ξ ≤ η ≤ N, where ξ FIG. 2. Potential energy of a single breakable spring with elastic constant h (top panel) and corresponding force (bottom panel). The quantity Y M is the elongation after which the spring breaks, resulting in a force equal to zero.

represents the domain wall or interface position. We remark that the last Nη sites of the chain are always considered disconnected from the bottom layer in order to simulate a possible existing initial fractured domain. This is coherent with the assumption of an initial crack extension in the classical Griffith criterion, which is the milestone of the linear elastic fracture mechanics [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF][START_REF] Pugno | The centenary of Griffith's theory[END_REF].

Based on previous key premises, the total energy of the system is

Φ H (y 1 , ..., y N , ξ ) = N ∑ i=0 k 2 (y i+1 -y i ) 2 + N ∑ i=1 l 2 (Y -y i ) 2 + ξ ∑ i=1 h 2 y 2 i + η ∑ i=ξ +1 h 2 Y 2 M , (1)
with the boundary conditions y 0 = 0 and y N+1 = Y . Here, the variables y i represent the vertical coordinates of the lattice points while ξ assignes the discrete interface position. Observe that the first addend of Eq.( 1), proportional to the elastic constant k, is the energetic contribution of the shear unbreakable springs. We remark that the use of the shear springs with elastic constant k is an approximation valid under the small deformation assumption (in our models the central nodes can only move vertically). The second addend, proportional to the elastic constant l, is the contribution of the vertical unbreakable springs that connect the upper layer to the inner lattice points. Finally, the energetic contribution proportional to the elastic constant h, regarding the breakable springs, is split into two addends: the first one corresponds to the unbroken springs (going from i = 1 to i = ξ ) and the second one corresponds to the broken ones (from i = ξ + 1 to i = η). Again, the assumption of linear elastic springs, as typical in linear elastic fracture mechanics, allow for a proper description of fracture effects under the hypothesis of small strains. Such a simplification is crucial for the following analytical treatment. The energy contribution associated with broken bonds corresponds to the surface energy of the two exposed sides of the fracture, originally introduced by Griffith in the overall energy balance and eventually yielding his classical stability criterion [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF].

The energy in Eq.( 1) can be rewritten as

Φ H = k 2 N ∑ i=1 2 + l k y 2 i + ξ ∑ i=1 h k y 2 i - N-1 ∑ i=1 y i+1 y i - N-1 ∑ i=1 y i y i+1 +kY - N ∑ i=1 l k y i -y N + 1 2 kY 2 + 1 2 lNY 2 + 1 2 hY 2 M (η -ξ ).
(2) We can introduce the following N-component vectors

⃗ y = (y 1 , y 2 , . . . , y N ) , (3)
⃗ v = (β , . . . , β , 1 + β ) , (4)
and the tridiagonal matrix

A =          a 1 -1 0 • • • 0 -1 . . . . .
0 • • • 0 -1 a N          , (5)
where the diagonal components a i are defined as follows

a i = 2 + α, if 1 ≤ i ≤ ξ , 2 + β , if ξ + 1 ≤ i ≤ N, (6)
with

α = l + h k , β = l k , (7)
measuring extension versus shear springs stiffness of the lattice. Adopting the matrix A and the vectors ⃗ y and ⃗ v, we can write Eq.( 2) as

Φ H = 1 2 k⃗ y • A⃗ y -kY⃗ v •⃗ y + 1 2 kY 2 + 1 2 lNY 2 + 1 2 hY 2 M (η -ξ ). (8
) This new compact expression is more suitable to evaluate the partition function of the system, defined by

Z H (Y ) = η ∑ ξ =0 R N e - Φ H K B T d⃗ y. (9)
Here, we have integrated the continuous coordinates y i and we have summed over the discrete or spin variable ξ identifying the interface position. We have

Z H (Y ) = η ∑ ξ =0 exp - kY 2 2K B T - lNY 2 2K B T - hY 2 M (η -ξ ) 2K B T I ξ , (10) where
I ξ = R N exp - k 2K B T ⃗ y • A⃗ y + kY K B T ⃗ v •⃗ y d⃗ y. (11)
By using the classical Gaussian integral,

R N e -⃗ y•M⃗ y e ⃗ w•⃗ y d⃗ y = π N det M e 1 4 ⃗ w•M -1 ⃗ w , (12)
which is valid for a positive definite symmetric matrix M (as can be shown for the tridiagonal matrix A ) and for any vector ⃗ w, we get

I ξ = (2πK B T ) N k N det A exp kY 2 2K B T ⃗ v • A -1 ⃗ v . (13)
Thus, the partition function can be written as

Z H (Y ) = η ∑ ξ =0 exp - kY 2 2K B T - lNY 2 2K B T - hY 2 M (η -ξ ) 2K B T × (2πK B T ) N k N det A exp kY 2 2K B T ⃗ v • A -1 ⃗ v . (14)
Since A depends on ξ , see Eqs.( 5) and ( 6), both A -1 and det A depend on ξ in the sum of Eq.( 14). In Appendix A, we discuss an efficient method to determine A -1 and det A for a tridiagonal matrix. This method will be used to obtain asymptotic expressions, useful to study the system behavior for large values of N (and for the thermodynamic limit). By introducing the quantity

q = 1 + β N -⃗ v • A -1 ⃗ v, (15)
we can write the partition function in the form

Z H (Y ) = η ∑ ξ =0 (2πK B T ) N k N det A exp -hY 2 M (η -ξ ) -kY 2 q 2K B T . ( 16
)
We can now evaluate the expected values of macroscopic quantities. For example, expectation value of the force applied to the system is [START_REF] Manca | Elasticity of flexible and semiflexible polymers with extensible bonds in the Gibbs and Helmholtz ensembles[END_REF][START_REF] Gibbs | Elementary Principles in Statistical Mechanics[END_REF][START_REF] Weiner | Statistical Mechanics of Elasticity[END_REF]

⟨ f ⟩ = -K B T 1 Z H ∂ Z H ∂Y , (17)
resulting in

⟨ f ⟩ = ∑ η ξ =0 (det A ) -1 2 exp -hY 2 M (η-ξ )-kY 2 q 2K B T q ∑ η ξ =0 (det A ) -1 2 exp -hY 2 M (η-ξ )-kY 2 q 2K B T kY. (18)
This expression gives a physical interpretation to the expectation value ⟨q⟩ of the variable q. Since ⟨ f ⟩ = k⟨q⟩Y , the quantity k⟨q⟩ represents the effective stiffness ⟨ f ⟩/Y of the system. Using a similar analysis, it is possible to obtain the average value of the number of unbroken bonds, which is given by

⟨ξ ⟩ = ∑ η ξ =0 (det A ) -1 2 exp -hY 2 M (η-ξ )-kY 2 q 2K B T ξ ∑ η ξ =0 (det A ) -1 2 exp -hY 2 M (η-ξ )-kY 2 q 2K B T . (19)
In Fig. 3, we show the main effects on the fracture behavior of both temperature T and discreteness parameter N. The main effect that can be observed is related to temperature, which is able to shift the value of the extension corresponding to the fracture of the system. In particular, we may observe that, as typically experimentally observed, the higher the temperature, the lower the force and the extension required to induce fracture. The model then predicts a thermally activated fracture phenomenon. As described in the following, in the thermodynamic limit (i.e., for N → ∞), this behavior can be theoretically interpreted as a phase transition. We can observe since now (see Fig. 3) that, as the discreteness parameter N increases, the force-displacements curves become sharper, increasing the brittleness of the system. In summary, the model FIG. 4. Comparison between the approximated quantities given by Eqs.( 28) and (31) (dashed lines) and the corresponding exact results in Eqs.( 18) and ( 19) (continuous lines) for ⟨ f ⟩/(NhY M ) (top panel) and ⟨ξ ⟩/N (bottom panel) versus the dimensionless extension Y /Y M . The thermal to elastic energy ratio K B T /(hY 2 M ) is set to 0.5 (blue or dark gray curves) and 2 (orange or light gray curves). The total number of units is set to N = 50 and N = 200 (see arrows). The dimensionless quantities are set to l/h = 1 and k/h = 1. We also considered η = N. exhibits a temperature dependent brittle behavior, which can be thoroughly described by analytic expressions for large values of N, obtained in the following Section.

III. ASYMPTOTIC BEHAVIOR OF THE ELASTO-FRAGILE MODEL

Eqs.( 18) and (19) determine the expected value of ⟨ f ⟩, the force applied to the system, and the average value ⟨ξ ⟩ of the number of intact bonds as functions of both the assigned displacement Y and of the temperature T . Here, to give a clearer physical interpretation of such results, we obtain analytical approximated relations, effective in the case of large values of N. In particular, for large values of N, we have (see Appendix B)

q ∼ lh l + h ξ k + ε, (20)
where

ε = β 2 + 4β -β 2 > 0 (21) and det A ∼ τ ξ α τ N-ξ β , (22)
with

τ s = 2 + s + √ s 2 + 4s 2 , s = α, β . (23)
By using Eq.( 18), we obtain the following asymptotic expression of the force-extension relation

⟨ f ⟩ ∼ ∑ η ξ =0 exp -ξ 2 ln δ + hY 2 M ξ 2K B T -Y 2 2K B T lh l+h ξ q ∑ η ξ =0 exp -ξ 2 ln δ + hY 2 M ξ 2K B T -Y 2 2K B T lh l+h ξ kY, (24)
where we introduced δ = τ α /τ β . After defining

z = exp - 1 2 ln δ + hY 2 M 2K B T - Y 2 2K B T lh l + h , (25)
and adopting similar considerations for the quantity ⟨ξ ⟩, we may easily deduce

⟨ f ⟩ ∼ kε + lh l + h ∑ η ξ =0 ξ z ξ ∑ η ξ =0 z ξ Y, (26)
⟨ξ ⟩ ∼ ∑ η ξ =0 ξ z ξ ∑ η ξ =0 z ξ , (27)
for large values of N. Combining Eqs.( 26) and ( 27), we find

⟨ f ⟩ ∼ kε + lh l + h ⟨ξ ⟩ Y. (28)
So, if we calculate the expectation value of intact elements ⟨ξ ⟩, we also deduce the force-extension relation ⟨ f ⟩-Y . To do this in explicit form, we have to evaluate the sums that appear in the expression of ⟨ξ ⟩, i.e.

η ∑ ξ =0 ξ z ξ = z 1 -(η + 1)z η + ηz η+1 (1 -z) 2 , ( 29
) η ∑ ξ =0 z ξ = 1 -z η+1 1 -z , ( 30
)
where we adopted the variable z defined in Eq.( 25). To conclude, we write

⟨ξ ⟩ ∼ 1 -(η + 1)z η + ηz η+1 1 -z η+1 z 1 -z , (31)
These results represent the approximated expressions for the response of the system under isometric condition for large values of N. As shown in Fig. 4 (for η = N), we can observe that the approximations given in Eqs.( 28) and (31) (dashed lines in the figure) are in perfect agreement with the exact results (continuous lines), previously obtained in Eqs.( 18) and [START_REF] Giordano | Effects of the orientational distribution of cracks in solids[END_REF]. We remark that the agreement is very good for different temperatures and for different (large) values of N. These approximated results are particularly useful to study the thermodynamic limit or, equivalently, to study the limiting case for N → ∞, as discussed below.

IV. THERMODYNAMIC LIMIT OF THE ELASTO-FRAGILE MODEL

We perform now the limit for N → ∞. Since previous asymptotic results in Eqs.( 28) and ( 31) depend on powers of z, we study the inequality z > 1. To verify this condition we need to set the argument of the exponential in Eq.( 25) larger than zero, obtaining the following condition on Y

|Y | ≤ l + h lh hY 2 M -K B T ln δ ≜ Y s , (32)
where we introduced the critical extension Y s , the physical interpretation of which will be given below. We also define a critical temperature T c for the system through the condition hY 2 M -K B T c ln δ = 0 (see Eq.( 32)) that, once solved for T c , gives

T c = hY 2 M K B ln   2 + l+h k + l+h k 2 + 4 l+h k 2 + l k + l k 2 + 4 l k   , (33)
where we used δ = τ α /τ β . This is a specific value of the temperature, depending on the main material parameters of the system, which corresponds to a phase transition, as discussed below. For the moment, we can write

z > 1 ⇐⇒ |Y | ≤ l + h l Y 2 M 1 - T T c ≜ Y s (34)
so that Y s = 0 when T = T c .

We first analyze the limit for N → ∞ of the average fraction of intact elements

⟨ξ ⟩ N = 1 -(η + 1)z η + ηz η+1 1 -z η+1 z 1 -z 1 N = 1 -z η+1 + (η + 1)(z η+1 -z η ) 1 -z η+1 z 1 -z 1 N = 1 N + η + 1 N 1 -z z -z -η z 1 -z . (35)
Now, we consider η = N -M, where M = φ N is the number of initially absent breakable springs (initial fractured domain).

Here, φ is the percentage of initially absent breakable springs over the total number N of elements. The fraction of intact elements is determined as follows lim M ) = 2 (in orange or light gray) while the total number of units for the large N approximation is set to N = 500. The dimensionless quantities are set to l/h = 1 and k/h = 1. Moreover, we adopted η = N corresponding to φ = 0.

N→∞ ⟨ξ ⟩ N = 1 -φ if z > 1 or |Y | < Y s , 0 if z < 1 or |Y | > Y s . (36)
This means that all elements are broken simultaneously (brittle fracture) when Y = Y s .

For the stress ⟨ f ⟩/N (density of force over the number of bonds), we get lim

N→∞ ⟨ f ⟩ N = lh l+h (1 -φ )Y if z > 1 or |Y | < Y s , 0 if z < 1 or |Y | > Y s . (37)
Thus (see Fig. 5), after an initial linear behavior, the stress collapses to zero at the extension threshold Y s . Together with the extension threshold Y s , we can therefore introduce the stress threshold σ s as follows

Y s =Y M l + h l 1 - T T c , (38)
σ s =hY M (1 -φ ) l l + h 1 - T T c . (39)
In other words, σ s is the value of ⟨ f ⟩/N in correspondence of We can observe that both quantities present a critical behavior corresponding to a phase transition for T = T c . We adopted l/h = 1, k/h = 1, and φ = 0. Y = Y s . We can say that Y s is the fracture extension while σ s is the fracture strength inducing the breaking process. The behavior of the system is shown in Fig. 5, where the dimensionless quantities ⟨ f ⟩/(NhY M ) and ⟨ξ ⟩/N are represented versus the dimensionless extension Y /Y M . We compared here the response for a large value of N and the thermodynamic limit. In this limit, the breaking is fully brittle with temperature dependent fracture extension and stress. Thus, as typical in collective phenomena of complex systems, although each breakable spring has a temperature-independent breakage threshold, the overall system exhibits a temperaturedependent fracture point due to the interactions between the springs and the thermal bath. In particular, the system undergoes a phase transition for T = T c , with both fracture extension and stress decreasing to zero at T = T c when the system breaks without any external mechanical actions. This is described in Fig. 6, where we plot the dimensionless fracture extension Y s /Y M and the dimensionless fracture strength σ s /(hY M ) versus the temperature ratio T /T c . This critical behavior corresponds to a classical second order phase transition. We also remark that Eq.( 39) represents an extension of the Griffith criterion of the linear elastic fracture mechanics, accounting for the additional effects of temperature. From Fig. 2, we see that the energy necessary to break an element is given by hY 2 M /2 and then the Griffith surface energy density γ s is proportional to hY 2 M /2. Equivalently, Y M is proportional to √ γ s and therefore it is easily seen that the fracture strength given in Eq.( 39) is proportional to √ γ s , exactly as in Griffith's criterion [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF]. Moreover, at constant temperature it is well seen that the breaking strength decreases if φ increases, which is exactly what the Griffith's criterion states [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF]. This means that, if the initial system is degraded, a smaller force is required to continue its mechanical degradation. Of course, our version is quantitatively different from the original one because of the simplified geometry we used. In particular, we do not consider the exact elastic energy distributed over the deformed continuum due to fracture. Our model, however, in-FIG. 7. Tensile strength of [START_REF] Jia | A Dugdale-Barenblatt Model for a Plane Stress Semi-Infinite Crack under Mixed Mode Concentrated Forces[END_REF]-oriented silicon nanowires as a function of temperature for wires with different diameters: comparison between molecular dynamics simulations results [START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF] (symbols) and our theory given by Eq.( 39) (continuous lines). The parameters used are reported in the main text. troduces thermal effects in brittle fracture and, in particular, shows the phase transition at the critical temperature T c given in Eq. [START_REF] Park | Cohesive Zone Models: A Critical Review of Traction-Separation Relationships Across Fracture Surfaces[END_REF].

To show the effectiveness of the obtained results, we analyze the temperature dependent fracture behavior of [START_REF] Jia | A Dugdale-Barenblatt Model for a Plane Stress Semi-Infinite Crack under Mixed Mode Concentrated Forces[END_REF]oriented silicon nanowires [START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF]. In Fig. 7, we compare the theoretical fracture force given by Eq.( 39) with the molecular dynamics results discussed in Ref. [START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF]. We observe that the theory well predicts the brittle fracture behavior of the nanowires both temperature-and diameter-wise (see Fig. 4(b) of Ref. [START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF]). In the figure the theoretical force in Eq. [START_REF] Miehe | A phase field model for rate-independent crack propagation: Robust algorithmic implementation based on operator splits[END_REF] has been divided by the area S pertaining to each breakable spring in order to obtain the stress σ = ⟨ f ⟩ /(NS) = σ s /S. The strain has been determined as ε = Y /ℓ, where ℓ is the characteristic lengthscale induced by the crystal structure. For all curves we adopted the parameters Y M = 1.78 × 10 -11 m, S = 2.27 × 10 -21 m 2 , ℓ = 1.82 × 10 -10 m, k = 88.4N/m, and K B = 1.38 × 10 -23 J/K. Moreover, for the blue curve (or dark gray, D = 5nm) we used l = 9.07N/m and h = 2.00N/m; for the yellow curve (or light gray, D = 6nm) we used l = 9.54N/m and h = 2.05N/m; for the orange curve (or intermediate gray, D = 7nm) we used l = 9.94N/m and h = 2.08N/m. While most of geometrical parameters were available in the original paper dealing with molecular dynamics simulations, the other physical parameters (in particular the elastic constants), were fitted to correctly reproduce the results. The elastic constants take effective values pertinent to the springs of our lattices and therefore cannot be directly obtained from the data available in the above papers. Interestingly, all the obtained (fitted) values are reasonable and consistent with the underlying physics of the system. In particular, the fact that h and l increase with the diameter is consistent with the results of Ref. [START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF], providing evidence that the nanowires Young modulus E increases with diameter (scale effect). This coherence is also quantitative since in our case we have E = ℓ/[(1/l + 1/h)S], which assumes the values 130GPa, 135GPa and 138GPa, for the three diameters 5nm, 6nm and 7nm, in agreement with Fig. 4(a . . , η) are also linked to the bottom layer (at y = 0) by η vertical softenable and breakable springs with elastic constant h when intact (colored in blue or dark gray), or p when softened (colored in green, springs with less coils). The broken elements are represented in orange and identified by a rupture in the springs. We remark that the first node (i = 0) and the last one (i = N + 1) are anchored to the bottom and the top layers, respectively.

Ref. [START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF]. The good agreement between theory and simulations makes us confident on the applicability of our theory to micro-and nanoscopic systems.

V. SOFTENING-FRACTURE MODEL

The previous model, useful to describe brittle fracture, is further generalized here to introduce a material ductile behavior of the elements possibly resulting in a brittle-to-ductile transition. Specifically, in the same spirit of the Dugdale-Barenblatt model [START_REF] Dugdale | Yielding of steel sheets containing slits[END_REF][START_REF] Barenblatt | The mathematical theory of equilibrium cracks in brittle fracture[END_REF], we introduce a cohesive zone between the elastic and fractured domains of the breakable springs, characterized by two different states before the broken configuration, depending on their extension y i (see Fig. 8). More precisely, each breakable spring presents an elastic constant h when its extension is less than the softening point Y p and a lower elastic constant p < h for larger extensions, until the breaking point corresponding to the extension Y b is attained and the link is broken (see Fig. 9). As we can see, each breakable element behaves as a spring of elastic constant h when -Y p ≤ y i ≤ Y p . Then, the spring is softened with an elastic constant p < h when Y p ≤ |y i | ≤ Y b . After the breaking point Y b , the potential energy is constant and therefore the resulting force is zero. Thus, the potential energy of a breakable spring is

U(y i ) =      1 2 hy 2 i if |y i | ≤ Y p , 1 2 py 2 i + ∆E if Y p ≤ |y i | ≤ Y b , 1 2 pY 2 b + ∆E if |y i | ≥ Y b . (40)
FIG. 9. Potential energy of a single softenable and breakable spring of elastic constants h and p (top panel) and corresponding force (bottom panel). We see that Y p is the elongation after which the spring is weakened or softened, and Y b is the elongation after which the spring is broken.

Given the two elastic moduli h and p (with p < h) and the energy gap ∆E > 0 , we obtain the softening point

Y p = 2∆E h -p , (41)
which must always satisfy the condition Y p < Y b . Thus ∆E + pY 2 p /2 = hY 2 p /2 is the energy necessary to weaken one breakable element of the system, and

∆E + pY 2 b /2 = hY 2 p /2 + p(Y 2 b -Y 2 p )/2
is the energy necessary to break the element. This reproduces in the discrete context considered here the Irwin generalization of the Griffith's criterion [START_REF] Irwin | Analysis of stresses and strains near the end of a crack traversing a plate[END_REF][START_REF] Irwin | Fracture[END_REF].

The total energy Φ H (y 1 , ..., y N , ξ , χ) of the system is

Φ H = N ∑ i=0 k 2 (y i+1 -y i ) 2 + N ∑ i=1 l 2 (Y -y i ) 2 + ξ ∑ i=1 h 2 y 2 i + χ ∑ i=ξ +1 p 2 y 2 i + ∆E + η ∑ i=χ+1 p 2 Y 2 b + ∆E . ( 42
)
Here we introduced the position ξ of the interface between intact and softened elements, the position χ of the interface between softened and fully broken elements and, finally, the position η of the interface between fully broken and initially absent elements. The value of η corresponds to the initial state of the system and is therefore fixed. The two interfaces at ξ and χ can move as a function of temperature and mechanical actions on the system. The aim of this section is to study the (quasi-static) evolution of these interfaces determining the fracture propagation phenomenon. The region between ξ and χ, characterized by softened elements, identifies the cohesive zone of the rupture phenomenon. When the cohesive zone is absent or negligible, the fracture is brittle; on the other hand, when the cohesive zone is not negligible, the fracture becomes ductile. Therefore, as we show in the following, this model allows to describe the brittle-to-ductile transitions. As in Section II, we remark that the use of the shear springs with elastic constant k is an approximation valid under the small deformation assumption (in our models the central nodes can only move vertically). We can rewrite Eq. [START_REF] Ren | Dual-horizon peridynamics[END_REF] as

Φ H = k 2 N ∑ i=1 2 + l k y 2 i + ξ ∑ i=1 h k y 2 i + χ ∑ i=ξ +1 p k y 2 i -2 N-1 ∑ i=1 y i+1 y i + kY - N ∑ i=1 l k y i -y N + 1 2 kY 2 + 1 2 lNY 2 + p 2 Y 2 b + ∆E (η -χ) + (χ -ξ )∆E. ( 43
)
As before, to simplify the mathematical structure of the energy function, we introduce the vectors in Eqs.( 3) and ( 4), and the tridigonal matrix in Eq.( 5), where the diagonal elements a i are now defined as follows

a i =      2 + α if 1 ≤ i ≤ ξ , 2 + γ if ξ + 1 ≤ i ≤ χ, 2 + β if χ + 1 ≤ i ≤ N, (44)
with the parameters

α = l + h k , β = l k , γ = l + p k , (45)
satisfying the condition β < γ < α. By introducing the matrix A and the vectors⃗ y and⃗ v, we are able to write the total energy as

Φ H = k 2 ⃗ y • A⃗ y -kY⃗ v •⃗ y + k 2 Y 2 + l 2 NY 2 + p 2 Y 2 b (η -χ) + (η -ξ )∆E. (46)
We suppose to embed the system in a thermal bath at temperature T and, assuming to be not far from the thermodynamic equilibrium, we can evaluate the partition function

Z H (Y ) = η ∑ χ=0 χ ∑ ξ =0 R N e - Φ H K B T d⃗ y. (47)
By using Eq.( 46), it can be evaluated as

Z H (Y ) = η ∑ χ=0 χ ∑ ξ =0 I ξ ,χ e λ ξ ,χ , (48)
where

λ ξ ,χ = - kY 2 2K B T - lNY 2 2K B T - ∆E K B T (η -ξ ) - pY 2 b 2K B T (η -χ), (49)
and

I ξ ,χ = R N exp - k 2K B T ⃗ y • A⃗ y + kY K B T ⃗ v •⃗ y d⃗ y. (50)
Using Eq.( 12) we get

I ξ ,χ = (2πK B T ) N k N det A exp kY 2 2K B T ⃗ v • A -1 ⃗ v . (51)
Summing up, we obtain the partition function as

Z H (Y ) = η ∑ χ=0 χ ∑ ξ =0 (2πK B T ) N k N det A e -kY 2 2K B T (1+βN-⃗ v•A -1 ⃗ v) × exp - ∆E K B T (η -ξ ) - pY 2 b 2K B T (η -χ) . (52)
In this case, by using Eq.( 17), we get that the expected value of the applied force is

⟨ f ⟩ = η ∑ χ=0 χ ∑ ξ =0 (det A ) -1 2 exp 2∆Eξ + pY 2 b χ -kY 2 q 2K B T q η ∑ χ=0 χ ∑ ξ =0 (det A ) -1 2 exp 2∆Eξ + pY 2 b χ -kY 2 q 2K B T kY, (53)
where we used the definition of q in Eq.( 15). This is the expression for the average value of the force necessary to impose the extension Y to the system. Similarly, we obtain the interfaces positions

⟨ξ ⟩ = η ∑ χ=0 χ ∑ ξ =0 (det A ) -1 2 exp 2∆Eξ + pY 2 b χ -kY 2 q 2K B T ξ η ∑ χ=0 χ ∑ ξ =0 (det A ) -1 2 exp 2∆Eξ + pY 2 b χ -kY 2 q 2K B T , ( 54
) ⟨χ⟩ = η ∑ χ=0 χ ∑ ξ =0 (det A ) -1 2 exp 2∆Eξ + pY 2 b χ -kY 2 q 2K B T χ η ∑ χ=0 χ ∑ ξ =0 (det A ) -1 2 exp 2∆Eξ + pY 2 b χ -kY 2 q 2K B T . ( 55
)
These results allow us to fully describe the fracture behavior for a ductile material. In particular, based on Eqs.( 54) and [START_REF] Parisi | Shear bands as manifestation of a criticality in yielding amorphous solids[END_REF], we are able to determine when the fracture is brittle, without the region of softened elements, or when the fracture is ductile, i.e., with a non negligible fraction of softened elements, representing the cohesive region.

We mentioned the Dugdale and Barenblatt models since historically they are the most important approaches to introduce a process zone in fracture phenomena. It is useful to remember that the original Dugdale model has been developed for plane strain conditions. Other approaches have been developed successively to consider plane stress conditions [START_REF] Hahan | Local yielding and extension of a crack under plane stress[END_REF][START_REF] Jia | A Dugdale-Barenblatt Model for a Plane Stress Semi-Infinite Crack under Mixed Mode Concentrated Forces[END_REF]. However, our model is composed of a quasi one-dimensional lattice of springs that does not allow the access to realistic elastic fields in the structure. Hence, it is difficult to quantitatively compare our results with elastic models in both plane stress and plane strain. Moreover, Dugdale model does not account for hardening phenomena, such as our approach, which is completely linear. In spite of these limitations, Dugdale model has been generalized for strain hardening materials [START_REF] Chen | Dugdale model for strain hardening materials[END_REF][START_REF] Neimitz | Modification of Dugdale model to include the work hardening and in-and out-of-plane constraints[END_REF]. It is also important to remember that the process zone in real situations extends beyond the fracture growth plane, a point neglected in both Dugdale's original model and ours. In real fractures, the actual deformation is represented by a complicated three-dimensional field, completely disregarded in our one-dimensional analysis. To conclude, the purpose of our models is not to improve aspects related to continuous elastic fields but rather to introduce the effects of temperature into a simplified model. With this in mind, our approaches are not created to replace classical ones but only to inform them of how temperature acts in fracture phenomena.

It is also important to discuss the physical meaning of the softened state of breakable springs. The ductility in metallic materials is related to a population of dislocations originated by the moving crack, generating a damaged zone near the crack tip with degraded elastic properties [START_REF] Rice | Ductile versus brittle behaviour of crystals[END_REF][114][START_REF] Rice | Dislocation nucleation from a crack tip: An analysis based on the Peierls concept[END_REF][START_REF] Thompson | Predicting whether a material is ductile or brittle[END_REF]. Since, we do not have the possibility to consider realistic dislocations in our model, we introduced the weakened state for the breakable springs, corresponding to the degraded elastic properties of the damaged zone. In the realistic case, the brittle-to-ductile transition is controlled by the competition between continuing the fracture (as in the brittle case) or using an amount of energy to generate dislocations that degrade the material. In our model, we have similar competition between the intact/broken switching (brittle regime) or the intact/softened switching (ductile regime). This competition is strongly influenced by temperature and our model explains this effect in detail. Besides metals, a similar damaged zone, describing the physical state of the material between the intact and the fully broken conditions, has been also observed in different systems including concrete [START_REF] Hu | Fracture energy and fracture pro-cess zone[END_REF], soft materials [START_REF] Creton | Fracture and adhesion of soft materials: a review[END_REF][START_REF] Trentadue | A predictive micromechanically-based model for damage and permanent deformations in copolymer sutures[END_REF], polymeric networks [START_REF] Awaja | Cracks, microcracks and fracture in polymer structures: Formation, detection, autonomic repair[END_REF][START_REF] Buche | Chain breaking in the statistical mechanical constitutive theory of polymer networks[END_REF], bones [START_REF] Kirchner | Ductility and brittleness of bone[END_REF].

VI. ASYMPTOTIC BEHAVIOR OF THE SOFTENING-FRACTURE MODEL

Once again to obtain clearer analytic results, we consider the behavior of systems with large values of N. We have (see Appendix B)

q ∼ β 2 ξ γ - ξ α - χ γ + χ β + ε, ( 56
)
where ε is given in Eq.( 21), and

ln det A N ∼ ln τ β + ξ N ln τ α τ γ + χ N ln τ γ τ β , (57)
where τ α , τ β and τ γ are given in Eq.( 23), for s = α, β , γ. Thus, we can write the expressions for the average force and the average interface positions as it follows

⟨ f ⟩ ∼ ∑ η χ=0 ∑ χ ξ =0 β 2 ξ γ -ξ α -χ γ + χ β + ε e q ξ ξ +q χ χ ∑ η χ=0 ∑ χ ξ =0 e q ξ ξ +q χ χ kY, (58)
⟨ξ ⟩ ∼ ∑ η χ=0 ∑ χ ξ =0 ξ e q ξ ξ +q χ χ ∑ η χ=0 ∑ χ ξ =0 e q ξ ξ +q χ χ , (59) ⟨χ⟩ ∼ ∑ η χ=0 ∑ χ ξ =0 χe q ξ ξ +q χ χ ∑ η χ=0 ∑ χ ξ =0 e q ξ ξ +q χ χ , (60)
where we introduced the quantities

q ξ = - 1 2 ln τ α τ γ - kY 2 β 2 2K B T 1 γ - 1 α + ∆E K B T , (61)
q χ = - 1 2 ln τ γ τ β - kY 2 β 2 2K B T 1
β - 1 γ + pY 2 b 2K B T . ( 62
)
Using the expressions for ⟨ξ ⟩ and ⟨χ⟩, we can rewrite ⟨ f ⟩ in the simpler form

⟨ f ⟩ ∼ ε + β 2 1 γ - 1 α ⟨ξ ⟩ + β 2 1 β - 1 γ ⟨χ⟩ kY. (63)
Hence, once we know the expected values ⟨ξ ⟩ and ⟨χ⟩ of the interfaces positions, we also know the force required to impose the extension Y . To simplify the notation, we introduce

w = e q ξ , (64)
z = e q χ . (65)
By Eqs.( 59) and (60), using Eqs.( 29) and ( 30),after long but straightforward calculations, we obtain

⟨ξ ⟩ = N ξ D , (66)
⟨χ⟩ = N χ D , (67)
where, for the sake of readability, we introduced

N ξ = w 1 -w 1 -z η+1 1 -z - 1 -(wz) η+1 1 -wz - (1 -w)wz (1 -wz) 2 1 -(wz) η (1 + η) + η(wz) η+1 , (68)
N χ = z (1 -z) 2 1 -z η (1 + η) + ηz η+1 - w 2 z (1 -wz) 2 1 -(wz) η (1 + η) + η(wz) η+1 , (69)
D = 1 -z η+1 1 -z -w 1 -(wz) η+1 1 -wz . ( 70
)
These results approximate the behavior of the fracture process in the presence of the softening phenomenon for large values of N. In particular, we can determine the limit for N → ∞ of the main observables, in order to provide a precise physical interpretation of the brittle-to-ductile transition. The obtained expressions depend on w η , z η and (wz) η and, since in our model η = N(1φ ) where φ is the fraction of initially absent breakable springs, they present an exponent going to infinity when N → ∞. We know that, when N → ∞, a generic power x N tends to infinity if x > 1 and tends to zero if |x| < 1, hence we study the three inequalities w > 1, z > 1 and wz > 1, which will be useful to better understand the system behavior. These inequalities are equivalent to study the positive character of their exponents q ξ , q χ , and q ξ + q χ . We start by setting the exponent of w larger than zero

q ξ = - 1 2 ln τ α τ γ - kY 2 β 2 2K B T 1 γ - 1 α + ∆E K B T > 0. (71)
In terms of Y , this inequality gives

|Y | < 1 kβ 2 1 γ -1 α (h -p)Y 2 p -K B T ln τ α τ γ ≜ Y ξ , (72)
where we introduced a first extension threshold Y ξ . We observe that kβ 2 1 γ -1 α is always positive because p < h by definition and, then, the argument of the square root is positive when the temperature T is smaller than the critical temperature T ξ defined as

T ξ = (h -p)Y 2 p K B ln τ α τ γ . ( 73
)
The meaning of Y ξ and T ξ will be clarified later. By setting the exponent of z greater than zero, we define the inequality It can be solved with respect to Y , eventually giving the result

q χ = - 1 2 ln τ γ τ β - kY 2 β 2 2K B T 1 β - 1 γ + pY 2 b 2K B T > 0. (74)
|Y | < 1 kβ 2 1 β -1 γ pY 2 b -K B T ln τ γ τ β ≜ Y χ , (75)
where we introduced a second extension threshold Y χ . As before, the quantity kβ 2 1 β -1 γ is always positive and therefore the whole square root argument is positive for values of the temperature below the critical temperature T χ defined as

T χ = pY 2 b K B ln τ γ τ β . (76)
As before, Y χ and T χ will be physically interpreted in the following. Finally, we set the exponent of wz greater than zero, which corresponds to q ξ + q χ > 0. We obtain the inequality

|Y | < 1 kβ 2 1 β -1 α pY 2 b + (h -p)Y 2 p -K B T ln τ α τ β ≜ Y ξ χ , (77)
where we introduced a third extension threshold Y ξ χ . Being

kβ 2 1 β -1 α
always positive, the square root has a positive argument when T < T ξ χ , where

T ξ χ = pY 2 b + (h -p)Y 2 p K B ln τ α τ β . (78)
Again, we will discuss in the following the physical meaning of Y ξ χ and T ξ χ . Summarizing these results, we can write

w > 1 q ξ > 0 ⇔ |Y | < (l + p)(l + h)Y 2 p l 2 1 - T T ξ , (79)
z > 1 q χ > 0 ⇔ |Y | < (l + p)Y 2 b l 1 - T T χ , (80)
wz > 1

q ξ + q χ > 0 ⇔ |Y | < pY 2 b + (h -p)Y 2 p lh l+h 1 - T T ξ χ , (81)
where we used the definition for the critical temperatures T ξ , T χ and T ξ χ previously introduced.

The three dimensionless extension thresholds

Y ξ /Y b , Y χ /Y b and Y ξ χ /Y b are plotted versus the dimensionless temperature K B T /(hY 2
b ) in Fig. 10. We can already anticipate that brittle or ductile behavior depends on the sign of Y χ -Y ξ . Indeed, we have a brittle fracture if Y χ < Y ξ , and a ductile fracture if Y χ > Y ξ (see Fig. 10). In the first brittle case, the rupture occurs for Y = Y ξ χ . In the second ductile case, the softening occurs for Y = Y ξ and rupture for Y = Y χ . We observe therefore that it exists a brittle-to-ductile transition temperature T * that separates the brittle behavior from the ductile one (see the yellow point in Fig. 10). This temperature is defined by equating Y ξ and Y χ , as follows When solved, this equation gives the value of T * as

(h -p)Y 2 p kβ 2 1 γ -1 α 1 - T * T ξ = pY 2 b kβ 2 1 β -1 γ 1 - T * T χ . (82)
K B T * hY 2 b = 1 α-β α Y p Y b 2 -β α α-γ ln τ α τ γ -β γ-β ln τ γ τ β . (83)
Interestingly, this quantity can be also explicitly written in terms of the elastic constants of the system

T * = (l + h)Y 2 p -lY 2 b K B l+h h-p ln τ α τ γ -l p ln τ γ τ β . ( 84
)
To justify the introduction of all these quantities and notations, we use now Eqs.( 63), ( 66) and ( 67) to observe the behavior of the system with different values of N and temperature T . In particular, we consider Fig. 10 and we show the system behavior for three values of the temperature belonging to the regions 0 < T < T * (brittle response, Fig. 11), T * < T < T ξ (ductile response, Fig. 12), and T ξ < T < T χ (over-ductile response, Fig. 13). We do not consider values of the temperature larger than T χ since, in this case, all elements are broken due to the only thermal effects, without the application of mechanical actions.

In Fig. 11, we can find the dimensionless force given by Eq.( 63) in the first panel, and the three quantities ⟨ξ ⟩/N, ⟨χξ ⟩/N and ⟨N(1φ ) -χ⟩/N representing the fraction of intact, softened and broken elements, calculated through Eqs.( 66) and [START_REF] Guarino | Thermally activated fracture of porous media[END_REF], respectively, in the second panel. For simplicity, we always considered φ = 0. Note that the dimensionless force is divided by N so as to be consistent with the definition of mechanical stress. We can see that, with a temperature in the range 0 < T < T * , the force drops to zero in correspondence to the extension threshold Y ξ χ , describing the simultaneous rupture of all elements. Indeed, it can be seen in the second panel that the elements change almost completely from the intact to the broken state, with a fraction of softened elements that is negligible. The response becomes increasingly sharp as the value of N increases. In particular, the fraction of softened elements decreases to zero for N growing. A direct transition from intact to broken elements without an intermediate phase is therefore observed. This confirms that the response is brittle for T < T * and the rupture of the system occurs in this case at the applied extension Y ξ χ .

In Fig. 12, we represent the same functions for a temperature in the range T * < T < T ξ . In this case, we observe an almost simultaneous transition of all breakable elements from the intact to the softened state at the extension Y ξ and a subsequent transition from the softened to the broken state at the threshold Y χ . This behavior reproduces a ductile fracture, and the intermediate phase, characterized by the softened elements, mimics the cohesive phase of the fracture phenomenon. Also, in the force-extension diagram we see a first peak in correspondence to the softening of the elements, and a second peak describing the actual rupture. This curve is sharper for high values of N and smoother for low values. The comparison of Figs. 11 and12 shows the transition from a brittle to a ductile fracture as temperature increases, T * being the threshold temperature between the two regimes. This transition is described by the bifurcation at T * exhibited in Fig. 10 (see the yellow point), which gives rise to the intermediate region with softened elements.

To complete the picture on the system behavior, we also show in Fig. 13 the extreme situation when temperature is in the range T ξ < T < T χ (over-ductile response). Since the temperature is larger than T ξ , at the beginning of the traction almost all elements are already in the softened state even without an applied mechanical action. As a result, we can observe only one transition between the softened state and the FIG. 13. Over-ductile response of the fracture phenomenon (T ξ < T < T χ ). Left panel: dimensionless force versus dimensionless extension for different values of N = 500, 1000, 1500 (as indicated by the arrow). Right panel: fraction of intact, softened and broken elements for different values of N = 500, 1000, 1500 (as indicated by arrows). Inset: same plot as Fig. 10, where the temperature used here is indicated by the black point. We adopted the parameters α = 7/5, β = 2/5, γ = 9/10, ∆E/(hY 2 b ) = 1/10, φ = 0 (i.e., η = N), and K B T /(hY 2 b ) = 6/5.

FIG. 14. Dimensionless force versus dimensionless extension for

different values of the thermal to elastic energy ratio K B T /(hY 2 b ) = {0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1} (as indicated by the arrow). We observe that the behavior of the model changes from brittle, at low values of temperature, ductile for intermediate temperatures, to over ductile at high temperatures. We adopted the parameters N = 1000, α = 7/5, β = 2/5, γ = 9/10, φ = 0, and ∆E/(hY 2 b ) = 1/10.

broken state at the extension threshold Y χ . Consequently, in this temperature range, we observe a brittle transition between thermally-softened and broken elements. Since this response is observed only after the classical ductile behavior, we called it over-ductile response.

To better visualize the transitions between the different fracture regimes, we show in Fig. 14 some force-extension curves corresponding to different temperatures, spanning over brittle, ductile, and over-ductile regimes. It is interesting to remark that, within the ductile fracture, the shape of the forceextension curve is smoother in correspondence to the system softening/breaking since the cohesive phase is able to absorb an amount of energy before the final rupture. We also note that, within the ductile regime, the softening stress is higher than the failure stress for lower temperatures and conversely the softening stress becomes lower than the failure stress for higher temperatures. This point will be further discussed below. It is important to underline that all curves seen in Figs.11, 12, 13, and 14 have been obtained through Eqs.( 63), ( 66) and ( 67) with a large, but finite value of N. In the following we also describe the thermodynamic limit N → ∞.

What has been described so far represents the modeling of brittle-to-ductile transition induced by thermal fluctuations. Our model also allows us to describe a parametric brittle-toductile transition, i.e., intrinsic to the structure of the system. This means that there can be systems that exhibit only brittle or ductile behavior, regardless of temperature. On the one hand, an example of always brittle system is shown in Fig. 15, left panel, where we represent the three dimensionless extension thresholds

Y ξ /Y b , Y χ /Y b and Y ξ χ /Y b versus the dimensionless temperature K B T /(hY 2 b ). We can see that Y χ < Y ξ χ < Y ξ
for any value of the temperature. This means that there is no temperature high enough to induce a ductile fracture (T * > T ξ ). On the other hand, an example of always ductile system is shown in Fig. 15, right panel. In this case, Y χ > Y ξ χ > Y ξ for any value of the temperature so that there is no temperature low enough to induce a brittle fracture (T * < 0). These two situations describe materials that are always brittle or always ductile, regardless of the considered temperature.

To conclude, we obtained two types of brittle-to-ductile transitions: a thermal transition, induced by the effects of thermal fluctuations, and a parametric transition, induced by the values of the elastic parameters. It is worth to point out that this rather rich fracture behavior has been obtained based on a minimal system depending on the competition between elastic, entropic and fracture energy terms, regulated by the temperature and material parameters.

VII. THERMODYNAMIC LIMIT OF THE SOFTENING-FRACTURE MODEL

To give an even clearer physical description, we deduce here analytic results in the thermodynamic limit, N → ∞. We start the analysis by examining the average value of the number of intact, softened and broken elements of the system. If we consider the brittle behavior, with 0 < T < T * , the thermodynamic limit gives lim

N→∞ ⟨ξ ⟩ N = 1 -φ if Y < Y ξ χ , 0 if Y > Y ξ χ , (85) lim
N→∞ ⟨χ -ξ ⟩ N = 0 for all Y, (86) lim
N→∞ ⟨N(1 -φ ) -χ⟩ N = 0 if Y < Y ξ χ , 1 -φ if Y > Y ξ χ . (87)
In this case, we observe a direct transition between intact and broken elements without going through the softened state. However, if we take into account the ductile behavior with

T * < T < T ξ , we obtain for N → ∞ lim N→∞ ⟨ξ ⟩ N = 1 -φ if Y < Y ξ , 0 if Y > Y ξ , (88) lim
N→∞ ⟨χ -ξ ⟩ N =    0 if Y < Y ξ , 1 -φ if Y ξ < Y < Y χ , 0 if Y > Y χ , (89) lim
N→∞ ⟨N(1 -φ ) -χ⟩ N = 0 if Y < Y χ , 1 -φ if Y > Y χ . (90)
In this case, we observe the emergence of a region with softened elements, corresponding to the cohesive zone. Finally, the over-ductile regime, characterized by

T ξ < T < T χ , for N → ∞ leads to lim N→∞ ⟨ξ ⟩ N = 0 for all Y, (91) lim
N→∞ ⟨χ -ξ ⟩ N = 1 -φ if Y < Y χ , 0 if Y > Y χ , (92) lim
N→∞ ⟨N(1 -φ ) -χ⟩ N = 0 if Y < Y χ , 1 -φ if Y > Y χ . (93)
In this regime, all elements are initially softened and therefore the single transition corresponds to their complete breaking.

Concerning the expected value of the force, from Eq.( 63), we can write

lim N→∞ ⟨ f ⟩ N = β 2 1 γ - 1 α ⟨ξ ⟩ N + β 2 1 β - 1 γ ⟨χ⟩ N kY, (94)
where we can substitute the values of ⟨ξ ⟩/N and ⟨χ⟩/N pertinent to each fracture regime. We remark that in Eq.( 94), we have canceled out the first term shown in Eq.( 63) since N → ∞. For the brittle behavior (0 < T < T * ), we have

lim N→∞ ⟨ f ⟩ N = β 2 (1 -φ ) 1 β -1 α kY if Y < Y ξ χ , 0 if Y > Y ξ χ . (95)
For the ductile behavior (T * < T < T ξ ), we have

lim N→∞ ⟨ f ⟩ N =        β 2 (1 -φ ) 1 β -1 α kY if Y < Y ξ , β 2 (1 -φ ) 1 β -1 γ kY if Y ξ < Y < Y χ , 0 if Y > Y χ . (96)
Finally, for the over-ductile behavior (T ξ < T < T χ ), we get lim

N→∞ ⟨ f ⟩ N = β 2 (1 -φ ) 1 β -1 γ kY if Y < Y χ , 0 if Y > Y χ . (97)
The thermodynamic limit behavior (N → ∞) of the intact, softened and broken elements together with the value of the stress ⟨ f ⟩/N is exhibited in Fig. 16, where all the three fracture regimes brittle, ductile and over-ductile are considered. We remark that the resulting overall picture is coherent with the plots in Figs.11, 12, 13, where the same quantities were represented for large, but finite values of N.

These results allow us to identify the values of stress corresponding to the behavioral transitions. In the brittle regime (0 < T < T * ) we identify the fracture or breaking stress corresponding to σ B = lim N→∞ ⟨ f ⟩/N, for Y = Y ξ χ , which assumes the value

σ B = (1 -φ ) kβ 2 1 β - 1 α pY 2 b + (h -p)Y 2 p 1 - T T ξ χ = (1 -φ ) lh l + h pY 2 b + (h -p)Y 2 p 1 - T T ξ χ , (98)
depending on the critical temperature T ξ χ (see fourth row, first panel, of Fig. 16). In the ductile regime (T * < T < T ξ ), we have a first transition coinciding with the softening of all elements. It represents the beginning of the plastic regime. Two different values of stress describe this transition: the upper yield strength σ + S = lim N→∞ ⟨ f ⟩/N (for Y = Y - ξ , i.e. on the left of Y ξ ), and the lower yield strength σ - S = lim N→∞ ⟨ f ⟩/N (for Y = Y + ξ , i.e. on the right of Y ξ ), given by

σ + S = (1 -φ ) kβ 2 1 β -1 α 2 1 γ -1 α (h -p)Y 2 p 1 - T T ξ = (1 -φ ) h 2 l + p l + h Y 2 p 1 - T T ξ , (99)
σ - S = (1 -φ ) kβ 2 1 β -1 γ 2 1 γ -1 α (h -p)Y 2 p 1 - T T ξ = (1 -φ ) p 2 l + h l + p Y 2 p 1 - T T ξ , (100)
which depend on the critical temperature T ξ (see fourth row, second panel, of Fig. 16). These two values are useful to calculate the stress jump corresponding to the softening mechanism (yielding) within the ductile regime

σ + S -σ - S = (1 -φ )Y p l(h -p) (l + h)(l + p) 1 - T T ξ , (101)
which is always positive since h > p. Still in the ductile regime (T * < T < T ξ ), we observe the second transition describing the complete failure of all the elements for a stress σ F = lim N→∞ ⟨ f ⟩/N, for Y = Y χ , assuming the value

σ F = (1 -φ ) k 1 β - 1 γ p 1 - T T χ β Y b = (1 -φ ) l l + p 1 - T T χ pY 2 b , (102)
depending on the critical temperature T χ (see fourth row, second panel, of Fig. 16). In the over-ductile regime (T ξ < T < T χ ), the complete breaking of the system occurs at the same stress σ F given in Eq.( 102) and shown in the fourth row, third panel, of Fig. 16. The behavior of these transition stresses is summarized in Fig. 17, where they are plotted versus the temperature T . In addition, different stress-extension curves are shown at different temperatures of interest. In the first two cases, A) and B), we observe a brittle behavior characterized by the breaking of the system when the stress reaches the value σ B and the extension the value Y ξ χ . While the case B) corresponds to a temperature slightly smaller than T * (brittle), the case C) represents a temperature slightly larger than T * , being therefore in the ductile region. We see here both the softening transition at Y ξ and the failure transition at Y χ . In this case C), the stresses satisfy the relationship σ - S < σ F < σ + S and then the softening peak (upper yield strength) is larger than the failure peak. We can now continue to increase the temperature until σ - S < σ F = σ + S , that is, until the softening peak is equal to the failure peak. This condition is fulfilled in the panel D) of Fig. 17, and it corresponds to the temperature T * * , defined as

T * * = h 2 l+h Y 2 p -p 2 l (p+l) 2 Y 2 b K B h 2 (l+h)(h-p) ln τ α τ γ -pl (p+l) 2 ln τ γ τ β . ( 103
)
In the stress-temperature plot, this temperature value T * * represents the intersection of the two curves σ F and σ + S versus T . Increasing the temperature further, we enter the region T * * < T < T ξ (once again ductile), represented in the panel E), where the failure peak is larger than the softening peak (upper yield strength), σ - S < σ + S < σ F . Finally, for values of temperature in the range T ξ < T < T χ , we are in the over-ductile regime and the softening peak disappears, remaining only the failure peak σ F for the overall system, as shown in panel F) of Fig. 17. In panels A), B), C), D), E), and F) of Fig. 17 we also represent the stress-extension response for a finite (large) value of N in order to show the good agreement between the approximated expressions and the thermodynamic limit. FIG. 16. Response of the system in the thermodynamic limit within the three temperature regimes 0 < T < T * (brittle, first column), T * < T < T ξ (ductile, second column), and T ξ < T < T χ (over-ductile, third column). We plotted the average number of intact (first row), softened (second row) and broken elements (third row), and the stress ⟨ f ⟩/N (fourth row), for N → ∞. To compact the notation, we defined φ = 1φ , and we introduced the characteristic stresses σ B , σ + S , σ - S , and σ F , as defined in Eqs.( 98), ( 99), ( 100) and [START_REF] Florio | A predictive model for the thermomechanical melting transition of double stranded DNA[END_REF].

As a conclusion to this discussion, we would like to point out that the strength (rupture stress) of the system as a function of temperature is finally represented by a discontinuous curve formed by the branch σ B for 0 < T < T * (brittle) and by the branch σ F for T * < T < T χ (ductile), as one can see in the panel σ -T of Fig. 17. This discontinuity can be easily explained by observing that the brittle-to-ductile transition involves the phenomenon of softening and thus the synchronized lowering of the elastic constant of all breakable elements. Since we are applying a stretching to the system controlled by the extension, the reduction of the overall elastic constant produces a consequent reduction in stress (which is therefore discontinuous). We further remark that the strength (σ B or σ F , depending on the temperature) is proportional to the factor 1φ , which represents the fraction of initially present elements (φ is in fact the fraction of initially absent elements). This is reminiscent of the Griffith criterion, stating that the stress at fracture is lower if the initial crack opening is larger [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF]. In our case, the initial crack opening is proportional to φ and, therefore, the Griffith criterion is respected. However, we add here the temperature dependent nature of this criterion, which is described by the classical term 1 -T T ξ χ in σ B , or 1 -T T χ in σ F , which represents the critical behavior eventually resulting in a genuine phase transition.

The brittle-to-ductile transition has been observed in GaN nanowires through molecular dynamics simulations [START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF], and Fig. 18 shows the comparison with our theoretical results. We considered a GaN nanowire oriented in the direction [0001], with a diameter of 1.92 nm and a length of 6.12 nm, as reported in Fig. 2(d) of Ref. [START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF]. The lateral facets of this system are oriented along the 11 20 side planes, as shown in Fig. 1(b) of Ref. [START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF]. In Fig. 18, the blue curves (or dark gray) represent the brittle behavior whereas the orange (or light gray) ones describe the ductile behavior. We remark that Eqs.( 95) and ( 96) define a relation ⟨ f ⟩ /N = F (Y ) where F is a given function. We have to introduce the real stress σ = ⟨ f ⟩ /(NS), where S is the area pertaining to each breakable spring, and the real strain ε = Y /ℓ, where ℓ is the characteristic lengthscale induced by the crystal structure. The stressstrain relation can be therefore written as σ = F (εℓ)/S, where F is the relationship defined by Eqs.( 95) and [START_REF] Wang | Near-ideal theoretical strength in gold nanowires containing angstrom scale twins[END_REF]. In Fig. 18, we adopted the parameters h = 84), the temperature T * * corresponds to σ + S = σ F , see Eq. [START_REF] Hutchinson | Singular behaviour at the end of a tensile crack in a hardening material[END_REF]. In panels A), B), C), D), E), and F) we also show the stress-extension response for N = 1000 (dim colors or light gray curves). FIG. 18. Tensile stress-strain curves for a GaN nanowire oriented in the direction [0001] (diameter of 1.92 nm, length of 6.12 nm). The lateral facets are oriented along the 11 20 side planes. Comparison between molecular dynamics simulations results [START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF] (dashed lines) and our theory given by Eqs.( 95) and ( 96) (continuous lines). The parameters used are reported in the main text. parameters were available in the original paper dealing with molecular dynamics simulations and the others were fitted to correctly reproduce the results. We then plotted the stressstrain curves for the two temperatures T = 300K (brittle behavior) and T = 1800K (ductile behavior). It is interesting to note that the structural parameters used in our model are able to predict the correct brittle-to-ductile transition as obtained through molecular dynamics simulations. Moreover, also the upper and lower yield stresses of the ductile behavior are in quite good agreement with simulations. We remark that in our model there is a single softening process and therefore we can see only one failure peak after the softening peak in the stress-strain relation. As discussed in the conclusions, the model could be generalized with more softening steps to describe real damage such as that of the nanowires studied here.

VIII. CONCLUSIONS

We proposed and studied two prototypical models able to describe temperature effects in fracture processes. The first is aimed at explaining the temperature dependent behavior for brittle systems, and the second at showing the complex damage processes occurring in the presence of ductile breakable links with a possibility of a brittle-to-ductile transition regulated by thermal effects. Both models are based on a simple lattice structure built through unbreakable and breakable springs. The system is confined between two layers, one fixed and one movable, and is designed in such a way that lifting the top layer results in a force experienced by the system, being able to generate fracture propagation. This structure is supposed to be embedded into a thermal bath at fixed temperature. Hence, the models are developed within the equilibrium statistical mechanics formalism. The difference between the two models lies in the behavior at rupture of the breakable springs. In the first model for brittle fracture, each breakable spring can be in two states, namely elastic or broken, depending on the extension applied to the element itself. The state transition occurs through the absorption of an amount of energy that corresponds to the typical surface energy of the Griffith criterion [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF]. In the second model, each breakable spring can be in three different states representing the elastic, softened, and broken regimes. The intermediate softened state is introduced to reproduce the possible ductile regime of the fracture process. In this case, we have a 'yielding' point between the elastic and the ductile regimes, followed by a final failure point corresponding to fracture. The transition between intact and softened states occurs after the yielding energy is absorbed, and the transition between softened and broken occurs through the absorption of another amount of energy corresponding to fracture. Thus, in both models the energy balance is similar to what is typically assumed in linear elastic fracture mechanics since Griffith's and Irwin's pioneering works [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF][START_REF] Irwin | Analysis of stresses and strains near the end of a crack traversing a plate[END_REF][START_REF] Irwin | Fracture[END_REF], with the fracture phenomenon regulated by elastic, damage and fracture (surface) energy. However, including thermal fluctuations sensibly modifies the results and add important features to the system description.

As for the model for brittle fracture, we obtain a temperature-dependent fracture stress and a corresponding fracture strain, representing a 'genuine' phase transition. Thus, we obtain a critical temperature at which both fracture stress and fracture strain are zero and therefore the material is always broken for supercritical temperatures. Interestingly, the obtained temperature-dependent strength is in good agreement with several experiments and molecular dynamics simulations as demonstrated previously. It is interesting to note that although breakable springs have a temperatureindependent breaking behavior, the overall system exhibits a breaking point that is highly dependent on temperature. This a typical case of a complex system with collective behavior, giving rise to a critical phenomenon. We argue that this effect is relevant in the case of weak links, such as hydrogen bonds in biological materials, or in rubber, where the elasticity has an entropic character, or in small size metallic or semiconductor systems, such as the considered nanowires [START_REF] Dijkstra | Temperature-dependent fracture mechanisms in ultra-high strength polyethylene fibers[END_REF][START_REF] Kang | Brittle and ductile fracture of semiconductor nanowires -molecular dynamics simulations[END_REF][START_REF] Wang | Atomistic simulations of the size, orientation, and temperature dependence of tensile behavior in GaN nanowires[END_REF][START_REF] Kang | Size and temperature effects on the fracture mechanisms of silicon nanowires: Molecular dynamics simulations[END_REF][START_REF] Wu | Nanowire Failure: Long = Brittle and Short = Ductile[END_REF][START_REF] Wang | Near-ideal theoretical strength in gold nanowires containing angstrom scale twins[END_REF][START_REF] Ast | The brittle-ductile transition of tungsten single crystals at the micro-scale[END_REF][START_REF] Munshi | Insights into the mechanical properties and fracture mechanism of Cadmium Telluride nanowire[END_REF][START_REF] Cheng | In Situ Nano-thermomechanical Experiment Reveals Brittle to Ductile Transition in Silicon Nanowires[END_REF][START_REF] Puglisi | Multi-scale modelling of rubber-like materials and soft tissues: an appraisal[END_REF][START_REF] Florio | A coarsegrained mechanical model for folding and unfolding of tropoelastin with possible mutations[END_REF][START_REF] Florio | A predictive model for the thermomechanical melting transition of double stranded DNA[END_REF].

The model with ductile breakable elements exhibits an even richer behavior. In this case it is the full response of the system that is temperature dependent. Indeed, we have demonstrated the existence of a brittle-to-ductile transition temperature T * (whose expression is obtained in closed form) that regulates the behavior of the fracture process. On the one hand, for temperatures lower than T * , we observe a brittle behavior characterized by a direct transition of the springs from the elastic to the broken state, without passing through the intermediate softened state. On the other hand, for temperatures higher than T * , we see that, as the extension applied to the system increases, first the springs soften (yielding point), and then they switch from the softened to the broken state (failure point). The intermediate softened region reproduces in this discrete context the cohesive zone of the classical Dugdale-Barenblatt model of the ductile fracture [START_REF] Dugdale | Yielding of steel sheets containing slits[END_REF][START_REF] Barenblatt | The mathematical theory of equilibrium cracks in brittle fracture[END_REF]. Of course, both yielding and failure point depend on temperature and are again characterized by phase transitions. In fact, importantly, we are not only able to calculate the brittle-to-ductile transition temperature, but also to predict the critical behavior of the upper and lower yield strengths, and the thermal properties of the fracture strength. Our model is also able to predict the existence of a special fracture regime, here called over-ductile, in which the temperature is high enough to damage all elements without mechanical action. In this situation, as the extension of the system increases, we observe the only transition from the softened to the broken regime.

From the methodological point of view, to elaborate the closed form expression of the partition function in both proposed models, we adopted specific techniques particularly suitable for calculating the determinant and inverse of tridiagonal matrices [START_REF] Usmani | Inversion of a tridiagonal Jacobi matrix[END_REF][START_REF] Usmani | Inversion of Jacobi's tridiagonal matrix[END_REF]. These approaches allow the derivation of exact solutions as shown in Appendix A, but also asymptotic approximations as discussed in Appendix B. Although these mathematical developments are relegated to the appendices, they are of crucial importance for obtaining the physical results on fracture processes.

We point out that even if the models here presented clarify fundamental aspects of thermally activated rupture phenomena, they should be generalized to take into account the complex reality of these processes. We want to mention here at least four points that partially limit the applicability of these models to real situations. The first issue concerns the spatial homogeneity of the adopted models. We have always considered all springs of discrete systems having the same mechanical behavior (in terms of elastic constants, failure thresholds, etc.). In reality, this is true only for perfect monocrystalline structures that are quite rare. It would be interesting to study these phenomena in disordered systems that, on the one hand, are more similar to several real structures, and on the other hand, may generate even more interesting critical behaviors typical of complex systems with quenched disorder [START_REF]Disorder and fracture[END_REF][START_REF] Parisi | Shear bands as manifestation of a criticality in yielding amorphous solids[END_REF][START_REF] Ponson | Crack propagation through disordered materials as a depinning transition: A critical test of the theory[END_REF][START_REF] Niccolini | Dimensional transitions in creeping materials due to nonlinearity and microstructural disorder[END_REF][START_REF] Sinha | Phase transitions and correlations in fracture processes where disorder and stress compete[END_REF][START_REF] Lebihain | Effective toughness of disordered brittle solids: A homogenization framework[END_REF][START_REF] Borja Da Rocha | Mean field fracture in disordered solids: Statistics of fluctuations[END_REF]. The second point to be explored is the kinetic of rupture processes. Here we have considered only quasistatic phenomena studied by means of equilibrium statistical mechanics. In real experiments, traction can be applied at different tensile velocities, and the response obviously depends on these traction rates [START_REF] Chen | Instability in dynamic fracture and the failure of the classical theory of cracks[END_REF][START_REF] Vincent-Dospital | How heat controls fracture: the thermodynamics of creeping and avalanching cracks[END_REF]. To model these phe-nomena one would have to adopt out-of-equilibrium statistical mechanics and then base the analysis on Langevin or Fokker-Planck methodologies [START_REF] Risken | The Fokker-Planck equation[END_REF][START_REF] Coffey | The Langevin equation[END_REF][START_REF] Benichou | Rate Dependent Response of Nanoscale Structures Having a Multiwell Energy Landscape[END_REF][START_REF] Benichou | The rate dependent response of a bistable chain at finite temperature[END_REF]. To conclude, the third point that could be improved concerns the fact that the softening process is restricted to a single step of reduction of the elastic constant of the breakable springs. In order to be more adherent to the physical reality of the yielding process one would have to imagine a series of steps where several reductions of the elastic constant take place progressively. In this sense, the yielding point would be implemented through a multi-softening process, more similar to what happens in real nonlinear materials. The fourth and final point concerns the too-simple geometry of our model, which should be improved (with 2D or 3D lattices) in order to be able to represent real elastic fields to be compared with models from continuum mechanics.

FIG. 19. The quantities S 2 (ξ )/N, S 1 (ξ ) and S 0 (ξ ) are obtained using Usmani relations (continuous curves) for different values of N = {10, . . . , 300} (with step of 10). We observe that, as N increases, the quantities approach their relative approximations obtained for large N (dashed lines).

We can therefore write S 2 (0) as

(τ -N β -τ N β ) + 2 + 2τ β (τ β -1) 2 (1 + τ β -τ N β -τ 1-N β ) ,
(B7) where we used several times the geometric sum. The expression for S 2 (0) given in Eq.(B7) is not transparent but, in the limit of large N, it can be approximated by

S 2 (0) ∼ N β - β 2 + 4β -β β 2 = N β +C(0), (B8)
where

C(0) = - β 2 + 4β -β β 2 . ( B9
)
This term represents the zeroth order correction (with respect to N), previously discussed. When ξ = N, we can obtain the result by simply substituting β with α, eventually obtaining

S 2 (N) ∼ N α - √ α 2 + 4α -α α 2 = N α +C(N), (B10)
where

C(N) = - √ α 2 + 4α -α α 2 . (B11)
Finally, we can write the general approximation for S 2 (ξ ) in the limit of large N as

S 2 (ξ ) ∼ N β + N α - N β ξ N +C(ξ ), if α, β ̸ = 0, (B12)
where C(0) is given in Eq.(B9), C(N) is given in Eq.(B11), and C(ξ ) assumes a constant value for ξ ∈ {1, . . . , N -1} (for large N), which is always in the range between C(0) and C(N).

We do not determine here this value since is not relevant for our analysis. Indeed, although the zeroth order term of S 2 (ξ ) is represented by three different values of the constant depending on ξ , in the application to the fracture problem we adopt the value C(0) in all calculations. It is not difficult to realize that this is the only value playing a role in our model since it describes the behavior of the system when ξ = 0, i.e. when all the breakable springs are fractured. In this condition, only one spring links together the two layers of the system and the constant C(0) is able to describe the exact stiffness of the resulting spring network. The other values C(ξ ), for ξ ̸ = 0, are negligible when N → ∞.

The approach used to find the approximation of S 2 (ξ ) for large N, can be also applied for S 1 (ξ ) and S 0 (ξ ). Concerning S 1 (ξ ), as shown in the center panel of Fig. 19, its value numerically obtained with Eq.(A6) approaches the constant value S 1 (0) as N increases. This value is therefore the approximation of S 1 (ξ ) for large N. Through previous definitions, we can write S 1 (0) as

S 1 (0) = N ∑ i=1 (A -1 ) i,N = N ∑ i=1 θ i-1 θ N = N ∑ i=1 G (β , i) G (β , N + 1)
, (B13) which, in the limit of large N, leads to

S 1 (ξ ) ∼ S 1 (0) ∼ 2 β + ∆ β . ( B14
)
We can observe that the exact values of S 0 (ξ ) approach the constant value S 0 (0) for large N, as one can see in the right panel of Fig. 19. We have the exact expression S 0 (0) = (A -1 ) N,N = θ N-1 θ N , (B15) which, in the limit of large N, gives S 0 (ξ ) ∼ S 0 (0) ∼ 1

τ β = 2 2 + β + β 2 + 4β , ( B16
)
as shown in the right panel of Fig. 19. Now we determine the value of the quadratic form in Eq.(B1), for large N, as

⃗ v • A -1 ⃗ v ∼ β 2 N β + 1 α - 1 β ξ +C(0) (B17) + 2β 2 β + β 2 + 4β + 2 2 + β + β 2 + 4β
.

To conclude, we recall the definition of q, stated in Eq.( 15), and we get

q = 1 + β N -⃗ v • A -1 ⃗ v ∼ lh l + h ξ k + ε, (B18)
where we introduced (B20)

ε = β 2 + 4β -β 2 . ( B19
Adopting the results of Appendix A, it is possible to evaluate the determinant of A for different ξ = 0, . . . , N. The results can be found in Fig. 20, from which we realize that the value of ln(det A )/N is approximated by a straight line that links together the values of ln det A /N when ξ = 0, and when ξ = N, in the limit of large N. It is simple to prove that for ξ = 0 we have ln θ N /N ≃ ln τ β when N → ∞, and similarly for

which proves Eq.( 22) of the main text.

We discuss now the same results for the model with the softening mechanism. In this case, the quadratic form ⃗ v • A -1 ⃗ v, can be written as in Eq.(B1) ⃗ v • A -1 ⃗ v = β 2 S 2 (ξ , χ) + 2β S 1 (ξ , χ) + S 0 (ξ , χ), (B23) where S 2 , S 1 and S 0 depends now on both interface positions ξ and χ. We start by analyzing the behavior of S 2 . Since S 2 (0, 0) ≃ N/β , S 2 (0, N) ≃ N/γ and S 2 (N, N) ≃ N/α, we get for large values of N the following expression

S 2 (ξ , χ) ≃ N β + N α - N γ ξ N + N γ - N β χ N +C, (B24)
where considerations similar to the previous ones confirm that the zero-th order term C assumes the same value in Eq.(B9). Similarly, we can prove that for large values of N the quantities S 1 and S 0 assume the same values obtained for the purely brittle fracture model, i.e. which are independent of ξ andχ. Adopting these approximations in the expression for q, we obtain for N → ∞

q ∼ β 2 ξ γ - ξ α - χ γ + χ β + ε, (B27)

FIG. 1 .

 1 FIG.1. Panel (a): scheme of a crack propagating within an arbitrary crystal lattice. Panel (b): reduced scheme of the elasto-fragile model, based on symmetry assumption. The central horizontal chain (colored in black) is composed of N + 1 linear elastic springs of elastic constant k. The nodes of this chain are connected to the top layer of the system (at y = Y ) with N vertical linear elastic springs of elastic constant l (colored in yellow or light gray). Moreover, the first η nodes (i = 1, . . . , η) are also linked to the bottom layer (at y = 0) through η vertical breakable springs of elastic constant h (intact in blue or dark gray; broken in orange or intermediate gray). We underline that the first node (i = 0) and the last one (i = N +1) are anchored to the bottom and the top layer, respectively. Hence, the first and the last shear springs fix the direction of the crack propagation from the right to the left of the system.

FIG. 3 .

 3 FIG. 3. Behavior of the brittle-fracture model with a variable number of units (N = {50, 125, 200} as indicated by arrows) and a variable thermal to elastic energy ratio K B T /(hY 2 M ) = 0.5 (blue or dark gray), and K B T /(hY 2 M ) = 2 (orange or light gray). The dimensionless quantities ⟨ f ⟩/(NhY M ) (top panel) and ⟨ξ ⟩/N (bottom panel) are represented versus the dimensionless parameter Y /Y M , where, in both cases, l/h = k/h = 1. We also fixed η = N, which means there are no missing or broken elements in the initial configuration.

FIG. 5 .

 5 FIG. 5. Comparison between the thermodynamic limit obtained for N → ∞ and the approximations obtained with large values of N (see arrows). We plotted ⟨ f ⟩/(NhY M ) (top panel) and ⟨ξ ⟩/N (bottom panel) versus the dimensionless extension Y /Y M . The thermal to elastic energy ratio is set to K B T /(hY 2 M ) = 0.5 (in blue or dark gray) and to K B T /(hY 2 M ) = 2 (in orange or light gray) while the total number of units for the large N approximation is set to N = 500. The dimensionless quantities are set to l/h = 1 and k/h = 1. Moreover, we adopted η = N corresponding to φ = 0.

FIG. 6 .

 6 FIG.6. Fracture extension Y s /Y M (in orange or light gray) and fracture strength σ s /(hY M ) (in blue or dark gray) versus the reduced temperature T /T c . All quantities are written in dimensionless form. We can observe that both quantities present a critical behavior corresponding to a phase transition for T = T c . We adopted l/h = 1, k/h = 1, and φ = 0.

FIG. 8 .

 8 FIG.8. Scheme of the fracture model with the softening mechanism. The central horizontal chain (colored in black) is composed by N + 1 linear springs with elastic constant k. The nodes of this chain are connected to the top layer (at y = Y ) by N vertical linear springs with elastic constant l (colored in yellow or light gray). The first η nodes (i = 1, . . . , η) are also linked to the bottom layer (at y = 0) by η vertical softenable and breakable springs with elastic constant h when intact (colored in blue or dark gray), or p when softened (colored in green, springs with less coils). The broken elements are represented in orange and identified by a rupture in the springs. We remark that the first node (i = 0) and the last one (i = N + 1) are anchored to the bottom and the top layers, respectively.

FIG. 10 .

 10 FIG. 10. Behavior of the dimensionless extension thresholds Y ξ /Y b , Y χ /Y b and Y ξ χ /Y b versus the dimensionless temperature K B T /(hY 2 b ) of the system. From the physical point of view, Y ξ χ (yellow or light gray curve) describes the brittle fracture below the transition temperature T * , and the couple Y ξ (blue or dark gray), Y χ (orange or intermediate gray) describes the ductile fracture above the temperature T * . The three curves Y ξ /Y b , Y χ /Y b and Y ξ χ /Y b versus K B T /(hY 2 b ) intersect at the bifurcation black point, characterized by T * . We adopted the parameters α = 7/5, β = 2/5, γ = 9/10, and ∆E/(hY 2 b ) = 1/10.

FIG. 11 .

 11 FIG. 11. Brittle response of the fracture phenomenon (0 < T < T * ). Left panel: dimensionless force versus dimensionless extension for different values of N = 100, 250, 500 (as indicated by the arrow). Right panel: fraction of intact, softened and broken elements for different values of N = 100, 250, 500 (as indicated by arrows). Inset: same plot as Fig.10, where the temperature used here is indicated by the black point. We adopted the parameters α = 7/5, β = 2/5, γ = 9/10, ∆E/(hY 2 b ) = 1/10, φ = 0 (i.e., η = N), and K B T /(hY 2 b ) = 3/10.

FIG. 12 .

 12 FIG.12. Ductile response of the fracture phenomenon (T * < T < T ξ ). Left panel: dimensionless force versus dimensionless extension for different values of N = 500, 1000, 1500 (as indicated by the arrow). Right panel: fraction of intact, softened and broken elements for different values of N = 500, 1000, 1500 (as indicated by arrows). Inset: same plot as Fig.10, where the temperature used here is indicated by the black points. We adopted the parameters α = 7/5, β = 2/5, γ = 9/10, ∆E/(hY 2 b ) = 1/10, φ = 0 (i.e., η = N), and K B T /(hY 2 b ) = 4/5.

FIG. 15 .

 15 FIG. 15. Behavior of the dimensionless extension thresholds Y ξ /Y b , Y χ /Y b and Y ξ χ /Y b versus the dimensionless temperature K B T /(hY 2 b ) for a system always brittle (left panel), and for a system always ductile (right panel). In the left panel we have Y χ < Y ξ χ < Y ξ and we adopted the parameters α = 13/10, β = 3/10, γ = 4/5, and ∆E/(hY 2 b ) = 1/5. In the right panel we have Y χ > Y ξ χ > Y ξ and we adopted the parameters α = 9/5, β = 4/5, γ = 13/10, and ∆E/(hY 2 b ) = 1/10.

  FIG. 17. Behavior of the characteristic stresses σ B , σ + S , σ - S , and σ F versus the temperature T , as defined in Eqs.(98), (99), (100) and (102) with φ = 0. Moreover, some stress-extension curves are plotted in correspondence of the following values of the temperature: A) 0 < T < T * ; B) T = (T * ) -(on the left of T * ); C) T = (T * ) + (on the right of T * ); D) T = T * * ; E) T * * < T < T ξ ; F) T ξ < T < T χ . While the temperature T * indicates the switching between brittle and ductile behavior, see Eq.(84), the temperature T * * corresponds to σ + S = σ F , see Eq.[START_REF] Hutchinson | Singular behaviour at the end of a tensile crack in a hardening material[END_REF]. In panels A), B), C), D), E), and F) we also show the stress-extension response for N = 1000 (dim colors or light gray curves).

S 2 (GG 2 (

 22 (β , i)G (β , N + 1i) G (β , N + 1) (β , i)G (β , N + 1j) G (β , N + 1) . (B6)Using the definition of G (γ, z) in Eq.(A12), and the properties of τ γ and ρ γ introduced in Eq.(A15), we get S

)Finally,

  Eqs.(B18) and (B19) prove Eqs.[START_REF] Giordano | Effects of the orientational distribution of cracks in isotropic solids[END_REF] and[START_REF] Giordano | Local elastic fields around cracks and their stress density of states[END_REF] of the main text.To complete this part, we study the approximation of det A for large values of N. Thanks to Usmani theory[START_REF] Usmani | Inversion of a tridiagonal Jacobi matrix[END_REF][START_REF] Usmani | Inversion of Jacobi's tridiagonal matrix[END_REF], we have det A = θ N .

FIG. 20 .

 20 FIG.[START_REF] Giordano | Effects of the orientational distribution of cracks in isotropic solids[END_REF]. The quantity (ln det A )/N is obtained using Usmani relations (continuous curves) for different values of N = {10, . . . , 300} (with steps of 10). We observe that, as N increases, (ln det A )/N approaches its approximation for large N (dashed line).

S 1 (ξ , χ) ∼ 2 β + β 2

 122

ACKNOWLEDGMENTS

The authors have been supported by Central Lille and Region Hauts-de-France under project MiBaMs. Moreover, the authors acknowledge stimulating scientific discussions with Giuseppe Florio and Giuseppe Puglisi.

Appendix A: Exact results for tridiagonal matrices Since the matrix A defined in Eq.( 5) is tridiagonal, we can analytically evaluate the inverse A -1 and the determinant det A [START_REF] Usmani | Inversion of a tridiagonal Jacobi matrix[END_REF][START_REF] Usmani | Inversion of Jacobi's tridiagonal matrix[END_REF]. We consider a generic tridiagonal matrix M , and we define its elements as M i,i = b i (main diagonal), M i,i-1 = a i (lower diagonal), and M i,i+1 = c i (upper diagonal). All other elements are zero. We can introduce the quantities θ i by means of the following recurrence relation

where, in particular, θ N = det M . Furthermore, it is possible to define the quantities φ i through the recurrence formula

where φ 1 = θ N = det M . These definitions can be used to determine the elements of the inverse matrix M -1 [START_REF] Usmani | Inversion of a tridiagonal Jacobi matrix[END_REF][START_REF] Usmani | Inversion of Jacobi's tridiagonal matrix[END_REF], as follows

(A3) By considering our particular case, the elements in the main diagonal of A are defined as a i = 2 + α for 1 ≤ i ≤ ξ , and a i = 2 + β for ξ + 1 ≤ i ≤ N. Moreover, we have that A i,i+1 = A i+1,i = -1 for the upper and lower diagonals. For this special situation, θ i and φ i are defined by the rules

and

Consequently, the elements of the inverse matrix A -1 are given by

Hence, we need to find θ i and φ i in order to obtain the inverse matrix elements. We start by evaluating θ i for i ≤ ξ . In this case, a i = 2 + α, and Eq.(A4) becomes

To find a solution, we substitute θ i = λ i in the last equation, and we obtain a second degree algebraic equation with solutions

Then, a generic solution for θ i , with i ≤ ξ , is given by the following linear combination

where we introduced ∆ α = α 2 + 4α. We can obtain the two coefficients A and B by the initial conditions in Eq.(A4). We obtain

Therefore, the final solution for θ i , when i ≤ ξ , is

where we introduced the function

A12) with ∆ γ = γ 2 + 4γ. If we introduce the parameters τ γ and ρ γ as follows

the function G (γ, z) can be written as

We note that

Now, we evaluate θ i when i ≥ ξ + 1. In this case, a i = 2 + β and Eq.(A4) becomes

As before, we find that the general solution is given by the linear combination

where ∆ β = β 2 + 4β . To find the coefficients C and D, we exploit the initial conditions θ ξ = G (α, ξ + 1) and θ ξ -1 = G (α, ξ ). Straightforward calculations lead to the solution for θ i , when i ≥ ξ + 1, in the form

A18) where we used the function defined in Eq.(A12). We consider Eq.(A5) and we proceed with the evaluation of φ i . We start with the case where i ≥ ξ + 1. In this condition, the recurrent equation becomes

which must be combined with the initial conditions in Eq.(A5). Eventually, we obtain φ i for i ≥ ξ + 1 in the form

We can find φ i when i ≤ ξ by using the two conditions

After straightforward calculations, we get for i ≤ ξ

The obtained values of θ i and φ i allow the calculation of A -1 and det A , useful to implement the determination of the partition function in Eq.( 16) and the quantities in Eqs.( 18) and [START_REF] Giordano | Effects of the orientational distribution of cracks in solids[END_REF]. Moreover, these results are useful to develop some asymptotic expressions in Appendix B.

Appendix B: Asymptotic analysis

Considering the brittle model with a large number N of units, it is possible to derive approximations to simplify the partition function and the main average quantities. More specifically, we can find approximations for ⃗ v • A -1 ⃗ v and for det A . We start our analysis by expanding the first quadratic form as follows

(B1) where we introduced S 2 (ξ ) as the sum over all the elements of the inverse matrix, S 1 (ξ ) as the sum over all the elements of the N-th column of the inverse matrix, and S 0 (ξ ) as the element (N, N) of the inverse matrix. We observe that these three quantities are in general function of ξ .

Exploiting the symmetry of the inverse matrix A -1 , we write S 2 (ξ ) as

The evaluation of S 2 (ξ ) for a matrix A that shows heterogeneous diagonal elements a i (ξ ̸ = {0, N}) can be done but is not straightforward. Fortunately, it is possible to observe that, in the limit of large N, the form of S 2 (ξ ) is a linear combination of the two values S 2 (0) and S 2 (N), each corresponding to a matrix with homogeneous diagonal. In fact, when ξ = 0, the diagonal components of A are all equal to 2 + β , and when ξ = N, the diagonal components of A are all equal to 2 + α.

With the help of the left panel of Fig. 19, we can observe that, as N increases, the form of S 2 (ξ ) numerically obtained through Eq.(A6) (continuous curves), approach the straight line joining S 2 (0) and S 2 (N) (dashed line). To improve the approximation, we observe that in the linear solution for S 2 (ξ ) we could also add a zeroth order term (with respect to N), represented by an additional small quantity C(ξ ), possibly dependent on ξ , but independent of N. Now, we analytically evaluate S 2 (0). By means of this value, it is also easy to obtain S 2 (N) simply substituting β with α. When ξ = 0, we have from Appendix A

where i = 1, . . . , N. In addition, Eq.(A6) gives

where ε is given by

This result corresponds to Eq.( 56) of the main text. Finally, we can also find an approximated expression for det A (ξ , χ) that now depends on both ξ and χ. Since previous approximation for ln det A (ξ ) was a linear function of ξ linking the values obtained for the two homogeneous matrices at ξ = 0 and ξ = N, we can now assume that the approximation for ln det A (ξ , χ) is a linear function in ξ and χ passing through the three points identified by (ξ , χ) = (0, 0), (N, 0) and (N, N). Hence, we assume that

where a, b and c are coefficients, which can be found with the assumptions a = ln τ β , a + cN = ln τ γ , and a + cN + bN = ln τ α , with