{'fraction_non_alphanumeric': 0.04107845556435798, 'fraction_numerical': 0.01418200236366706, 'mean_word_length': 4.663847002574476, 'pattern_counts': {'":': 0, '<': 35, '<?xml version=': 0, '>': 28, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Thermal control of nucleation and propagation transition stresses in discrete lattices with non-local interactions and non-convex energy

Non-local and non-convex energies represent fundamental interacting effects regulating the complex behavior of many systems in biophysics and materials science. We study one dimensional, prototypical schemes able to represent the behavior of several biomacromolecules and the phase transformation phenomena in solid mechanics. To elucidate the effects of thermal fluctuations on the non-convex non-local behavior of such systems, we consider three models of different complexity relying on thermodynamics and statistical mechanics: (i) an Ising-type scheme with an arbitrary temperature dependent number of interfaces between different domains, (ii) a zipper model with a single interface between two evolving domains, and (iii) an approximation based on the stationary phase method. In all three cases, we study the system under both isometric condition (prescribed extension, matching with the Helmholtz ensemble of the statistical mechanics) and isotensional condition (applied force, matching with the Gibbs ensemble). Interestingly, in the Helmholtz ensemble the analysis shows the possibility of interpreting the experimentally observed thermal effects with the theoretical force-extension relation characterized by a temperature dependent force plateau (Maxwell stress) and a force peak (nucleation stress). We obtain explicit relations for the configurational properties of the system as well (expected values of the phase fractions and number of interfaces). Moreover, we are able to prove the equivalence of the two thermodynamic ensembles in the thermodynamic limit. We finally discuss the comparison with data from the literature showing the efficiency of the proposed model in describing known experimental effects.


Introduction
Several natural and artificial systems, as typically encountered in biology and modern nanotechnology, exhibit a combination of fundamental effects of nonlocality and non-convexity, resulting in a variety of rather complex physical responses [1][2][3][4]. A deep understanding of these phenomena is therefore essential for the analysis and design of such systems. The non-convexity feature is related to the possibility that the potential energy of the system may have different basins leading to wiggly energy landscapes and possibly many competing metastable states. Often, these systems are composed of several units (which may be identical, or inhomogeneous) and each unit is characterized by a bistability (or in general a multistability). This assumption describes the possibility for each unit to be in two (or more) distinct states. In this framework, non-locality describes the possibility that the units can be in strong interaction with each other and thus the state in which one unit is found affects the state of the others, particularly those that are spatially closer. In the biological context this form of interaction is sometimes referred to as cooperativity [5] and may be fundamental in many crucial biological and medical phenomena such as protein folding/unfolding, DNA degradation and resulting diseases [6,7]. The complex response of such systems is the result of the energy competition among the many metastable configurations, regulated by the temperature controlled exploration of the overall energy landscape. We show how three different statistical mechanics approaches, taking care of the elastic properties of the different configurations, can give important insights for many observed physical and biophysical phenomena.
In this paper, we are interested in the observed thermal effects in previous examples, that become increasingly important as the size of the system decreases. The analysis is therefore relevant in systems of nanoscopic dimensions [53] or when, as in biological and polymeric soft matter, the competing contributions are of entropic type with small energy differences and low barriers [54]. For this reason, in order to have a correct physical description of the static and dynamic features of these systems, it is not sufficient to rely on classical discrete or continuous mechanics, but we must take into account the equilibrium or non-equilibrium statistical mechanics. So doing, we can obtain important information on the thermodynamic picture of the folding/unfolding process in a macromolecular chain or a detailed description of the thermal effects on the microstructure evolution in a two-phase solid material.
Concerning materials science, and in particular the mechanics of solids for multiphase materials, many approaches can be found in the literature to describe the microstructure evolution. The Ericksen pioneering work proposed a variational energetic approach in the context of non-linear continuum elasticity theory with non-convex energy densities [36]. This methodology has been further generalized to describe phase transformations at the microstructure level [37][38][39]. Nevertheless, these continuum variational approaches neglect interfacial energy effects and non-local interactions, which are crucial contributions for the description of the realistic microstructure evolutions [40]. As a matter of fact, the minimization of the non-convex elastic energy without nonlocal interactions cannot completely describe the nucleation and propagation of finite domains [41,42]. Therefore, surface energy contributions have been introduced by means of higher gradient energy terms [43,44], and by non local interactions [45]. A similar research line, in the context of the discrete mechanics, has been developed from the pioneering work of Müller and Villaggio [46]. The basic model is composed of a one-dimensional lattice of units with a nonconvex potential energy and an intrinsic length-scale [47]. This scheme allows the description of energy barriers, metastable states, quasi-plastic and pseudoelastic behaviors [48,49]. Also in this context, the model has been extended with non local energy terms able to capture the different features of phase nucleation and propagation [50]. Further generalizations consider the influence of boundary conditions, enabling a more detailed identification of the internal and boundary phase nucleations [51,52]. Recent models investigate the austenite to martensite phase transitions in wires, eventually describing the shape memory alloy behavior under uniaxial tension [55][56][57].
Many similar methods have been elaborated to describe the configurational transitions in biological macromolecules (mainly proteins), undergoing folding/unfolding processes. In particular, those theories are able to explain the saw-tooth-like force-extension response observed in several experiments. A model has been proposed for macromolecules unfolded in atomic force microscopes, and validated for titin and RNA hairpins [58]. An approach based on the equilibrium statistical mechanics is based on a Landau-like free energy and predicts a sequence of first-order phase transitions in correspondence to the unfolding processes [59,60]. Further investigations are based on the energy minimization of a bistable system and agree with the pattern observed in titin experiments [61]. Also the Monte-Carlo implementation of a two-state theory for single-molecule stretching experiments has been proposed [62]. Finally, the mechanical unfolding of proteins has been also modeled through domains interactions described by the Ising model [63,64].
In order to properly introduce the thermodynamics on non-convex, nonlocal systems, we adopt the method of the spin variables, by extending classical one-dimensional schemes to consider the fundamental effect of stiffness and elasticity of the different states. The first theoretical ideas introducing this technique can be found in the early models of the bio-mechanical response of skeletal muscles [19,20]. This technique has been further generalized to study different multistable systems [22-25, 54, 65, 66], and macromolecular chains [67][68][69][70][71][72]. This approach is based on the introduction of a series of discrete variables (the so-called spin variables), which are able to identify the state of the units. In other words we associate to each unit a sort of "bit", identifying the folded or unfolded state of the unit. So doing, we can consider two separated and different quadratic functions representing the wells of the potential energy, instead of the more complicated original bistable function. The introduction of the spin variables strongly simplifies the calculation of the partition functions and, consequently, the analysis of the macroscopic thermodynamic quantities. Indeed, in order to calculate these partition functions, we sum over the spin variables and we integrate the classical continuous variables. Since the separated wells are represented by quadratic terms, the integration can be performed straightforwardly since it acts on Gaussian functions. This theoretical approach is therefore able to yield closed form results useful to better understand the underlying physics. Moreover, by means of this technique we can study different ensembles of the statistical mechanics, corresponding to different mechanical boundary conditions. From one hand, we can analyze the behavior of the system under an applied force (isotensional condition), corresponding to the Gibbs ensemble. On the other hand, we can also investigate the features of the system with prescribed extension (isometric condition), corresponding to the Helmholtz ensemble. As described in Ref. [67], these boundary Fig. 1 Scheme of the discrete system with variable phase configuration. Top panel: isotensional loading (Gibbs ensemble) of the chain with applied force f . Bottom panel: isometric loading (Helmholtz ensemble) of the chain with the prescribed extension x N . Each chain unit, in both ensembles, feels an Ising interaction energy J due to the nearest neighbors. As an example, in both Gibbs and Helmholtz schemes, we show two interfaces between folded and unfolded domains.
conditions can be considered as limiting configurations of the more realistic case of an elastic device interacting with the system. Moreover, from a theoretical point of view, the comparison of the force-extension response within different ensembles is an important point useful to understand the concept of ensembles equivalence in the thermodynamic limit, as discussed below. The spin variable method has been successfully used to study the nanomechanics of macromolecules [67][68][69][70][71][72], the denaturation of macromolecules [73,74], the non-local or cooperative effects [54,64], skeletal muscles [22][23][24][25], and systems with transitions between unbroken and broken states [65,66,[75][76][77].
Here, we analyze a discrete chain with N bistable units (see Fig.1), where the two potential energy wells of each unit are characterized by different elastic constants k (pristine or folded state) and αk (extended or unfolded state), where α > 0 (see Fig.2). In addition, the two states are separated by an energy jump ∆E, representing the (Helmholtz) transition energy, and have their equilibrium lengths equal to ℓ and χℓ, with a length rise of ∆x = ℓ(χ − 1), where χ > 1 (see Fig.2). The assumption of considering two different elastic constants of the wells, already considered in the purely mechanical case in Ref. [48], involves important novelties compared to the case with identical constants. Indeed, when the folded and unfolded elastic constants are equal (α = 1), the conformational transitions correspond to a temperature independent average plateau force [54,[67][68][69][70][71][72]. This result can be simply explained in the framework of the Bell relation f = ∆E/∆x, discovered in the context of cell adhesion [62,78,79]. This plateau force, which depends neither on the spring constant nor on the temperature T , can be explained as follows. We consider two poten- to the folded and unfolded states of the unit under force f when α = 1. In both cases, the equilibrium lengths are defined by ∂U f /∂x = 0 and ∂U u /∂x = 0 and we get x f = ℓ + f /k and x u = χℓ + f /k. Finally, the unfolded configuration is more favorable than the folded one when U u (x u ) < U f (x f ), which corresponds to f > ∆E/∆x. The quantity f M = ∆E/∆x is the so-called Maxwell force for the case with α = 1. This approach can be easily generalized to the case with α = 1 and T = 0 (purely mechanical behavior). The same analysis yields in fact the following quadratic equation for the transition force where ∆x = ℓ(χ − 1), as before. This value of force has two important properties (both valid for α < 1 and α > 1): (i) on the plane (U, x) it is represented by an inclined straight line that is the common tangent to the two parabolas of the wells (see solid green lines in Fig.2, top panels); (ii) on the plane (dU/dx, x) is represented by a horizontal straight line (see solid green lines in Fig.2, bottom panels) that makes equal the two areas of the indicated triangles (see shaded regions in Fig.2, bottom panels). This last property gives the name Maxwell force to the transition force because the equality of the two triangles is reminiscent of Maxwell's construction on the isothermal Van der Waals curves in the pressure-volume plane of a real gas [80]. The situation becomes much more complicated when α = 1 and T > 0 and the transition force f M , based on statistical mechanics analysis, is temperature dependent. Indeed, the asymmetry of energy wells makes entropic contributions crucial. It means that a temperature dependent term must be added to Eq.(1) when α = 1 and T > 0. This point is the focus of this paper and it is extensively examined in the following development.
Other two special values of the force are of interest in this papers and represent the maximum (roof) f R and the minimum (floor) f F values of the force, such that for f ∈ (f F , f R ) both the first phase and the second phase can exist. They can be simply obtained as where x * represents the point of intersection of the two parabolas (where the top of the energy barrier is reached, see Fig.2). We also introduce in the discrete lattice of non-convex elements (see Fig.1) a non-local interaction described through an Ising scheme. This feature is crucial to consider a form of cooperativity in the biological context or, equivalently, an interface energy between folded and unfolded domains in the materials science context. The Ising interaction energy J can be considered positive (cooperative case) when adjacent units prefer to be in the same state, and negative (anti-cooperative case) when they prefer to be in two different states. To focus mainly on the interesting temperature dependent force plateau behavior, in this study we will discuss only the ferromagnetic-like interactions, with J > 0. The important theoretical novelty with respect to classical spin models in physics is the fundamental role of elasticity and, in particular, the effect generated by the different elastic behavior of the two phases. Furthermore, for these The folded well corresponds to S i = −1 and it is characterized by L(−1) = ℓ, Q(−1) = 0, K(−1) = k. The unfolded well, identified by S i = +1, is defined by L(+1) = χℓ, Q(+1) = ∆E, K(+1) = αk. Bottom panels: we show the quantity dU/dx versus x, exhibiting the characteristic force jump. Both the cases with α < 1 and α > 1 are shown on the left and on the right, respectively. In all plots, the solid green lines represent the transition paths (at T = 0), discussed in the main text.
systems, we can determine the analytic expression of the partition function in both Gibbs and Helmholtz ensembles (isotensional and isometric conditions, as shown in Fig.1). This can be done by means of an ad hoc implementation of the transfer matrix technique for the Gibbs case [81], and by using the Laplace transform relationship between the partition functions of conjugated ensembles for the Helmholtz case [82].
The first analysis, within the Gibbs ensemble, provides evidence that the force plateau is temperature dependent and we obtain its expression in the limit of large values of N and large (positive) values of the ratio J/K B T (K B being the Boltzmann constant), which represents the strongly ferromagnetic case in the thermodynamic limit. Defining β = 1/K B T , we can introduce the quantityβ = Jβ that represents the competition between the Ising interaction energy and the entropic contributions. Moreover, for isometric loading, in the Helmholtz ensemble, we find a peak force at the beginning of the plateau, representing the nucleation of a new domain with unfolded units. This is an important feature, typically observed in experimental measurements [83][84][85][86][87][88], and in molecular dynamics simulations of the microstructure evolution in nanosystems [53,[89][90][91][92][93][94][95][96][97][98][99][100][101]. In this Ising model, the microstructure evolution of folded and unfolded domains is free, and regulated by the competition between interface energies and entropic contributions and, as we will show, the stress peak is an effect induced by the introduction of interface energy terms. The number of interfaces may vary in the whole range between 0 and N − 1 and typically it increases with the temperature [54]. The knowledge of the partition function allows a full analysis of the configurational properties of the system leading to the quantitative evaluation of the average number of unfolded units and the average number of interfaces. In particular, this allows us to observe that the microstructure evolution, under increasing extension, is characterized by a single moving domain wall between folded and unfolded regions only wheñ β ≫ 1, i.e. when the system is strongly ferromagnetic, whereas forβ ≪ 1, i.e. when temperature increases, entropic energy terms favor solutions with an increasing number of interfaces.
An interesting point concerning these systems is the equivalence of the ensembles in the thermodynamic limit (i.e., for N → ∞) [102][103][104][105][106]. Two conjugated ensembles are said to be equivalent when the macroscopic behavior described by the force-extension relation is the same for N → ∞. In general, it is difficult to prove for a given system if two statistical ensembles are equivalent. Although there are some particular rules, there are no general criteria or theorems for determining whether a system satisfies such an equivalence [105]. Several examples of non-equivalence are well known in the literature [107][108][109][110][111][112]. In our case, the analysis of the equivalence is rather difficult, mainly due to the overly complicated mathematical form of the Helmholtz partition function. But, since the systems with positive and sufficiently intense Ising interaction (β ≫ 1) are the most interesting for practical applications, we can limit our analysis to the case of a system where the number of interfaces can take only the values 0 and 1. This observation is at the origin of the second approach proposed in the paper, called the zipper model, previously adopted in other statistical mechanics investigations [113][114][115][116]. The main assumption is the analysis of solutions with none or one interface, which correspond to the previous Ising scheme only whenβ ≫ 1. This simplification makes the thermodynamic limit analysis more transparent and, in this zipper case, we can explicitly prove the equivalence of the isotensional and isometric ensembles for N → ∞. In addition, the result obtained under the zipper assumption, within the Helmholtz ensemble, can be further simplified by means of the stationary phase method (large values of N ), leading to the third approach here discussed to describe non-convex discrete systems with non-local interactions. This final approximation is particularly useful since yields a compact mathematical expression for the force-extension curve and it allows the explicit calculation of the force peak in the Helmholtz ensemble, which is a crucial quantity in several experiments and numerical simulations, as discussed above. From one side, it leads to draw some comparisons between our theoretical results and data from the literature; from the other side, it suggests both the possibility of designing new materials with required transition properties, and the possibility of controlling them through external thermal fields.
The structure of the paper is the following. In Section II we introduce the Ising scheme and we discuss in detail the proposed solutions for the Gibbs and Helmholtz ensembles. In Section III we explain how to obtain the thermodynamic limit within the Gibbs ensemble and we discuss the important effect of temperature dependent force plateaux. In Section IV we introduce the zipper model: we discuss the Gibbs and Helmholtz ensembles and we draw a comparison with the previous Ising scheme. In Section V we study the approximation based on the stationary phase and we perform a detailed analysis of the first peak force within the Helmholtz ensemble. Finally, in Section VI, to support the obtained analytical results, we discuss some explicit, quantitative comparisons with data from the literature concerning the behavior of nanowires with microstructural evolution.

Two-state chain with Ising non-local interactions
Consider a discrete chain of N two-state elements, each described by a bistable potential energy (see Figs.1 and 2), that interact also non-locally. We distinguish the two phases using the spin variable S i assuming values in {−1, +1}.
In particular S i = −1 corresponds to the first well (folded element), whereas S i = +1 corresponds to the second well (unfolded element). While in a previous work the authors considered identical wells [54], the main hypothesis of this work is that the two phases are characterized by two different elastic constants K(S i ), associated with two natural lengths L(S i ) and two basal energies Q(S i ). Specifically, we assume, without loss of generality, that L(−1) = ℓ, L(+1) = χℓ, Q(−1) = 0, Q(+1) = ∆E, K(−1) = k, K(+1) = αk, where α > 0, χ > 1 and ∆E is the energy jump between the states. In these hypotheses, and introducing Ising non-local interaction terms, the overall Hamiltonian assumes the compact form Here, the non dimensional parameter λ i = (x i − x i−1 )/ℓ is the ratio between the i-th spring length and the folded rest length L(−1) = ℓ, i.e. the spring stretch and λ 0 (S i ) = L(S i )/ℓ, i.e. the natural (zero-force) spring stretch. The parameter J measures the non-local interaction strength. We remark that J > 0 corresponds to the ferromagnetic case, favoring phases coalescence. It is useful to take into account the adimensional HamiltonianH, obtained dividing H by the interface energy J whereQ(S i ) = Q(S i )/J,K(S i ) = K(S i )ℓ 2 /J and ̟ = 1 is a constant that we inserted in order to be able, later on, to calculate the average number of interfaces. We also introducek = kℓ 2 /J so thatK(−1) =k andK(+1) = αk. We remark that the spin variables approach can be adopted only when we work not far from the thermodynamic equilibrium [68,71]. Indeed when rate effects are considered, the relaxation times of the system strongly depend on the energy barriers between the potential wells, which are neglected within our approach (see Ref. [117] and recent generalizations in Refs. [118,119]). On the other hand, in the rate independent regime considered here, this approach allows us to describe non-convexity through the discrete parameters S i , in the sense that a for fixed phases configuration (S i , i=1,...,N ), the energy is convex with respect to other (stretch) variables. This energy structure also ensures that the equilibrium solutions represent local elastic energy minima (metastable equilibrium states). It is also important to remark that the application of the spin variable approach is correct only when the energy barrier between the two states is sufficiently larger than the thermal energy K B T . For completeness, in Appendix A, we prove that the Hamiltonian function given in Eq.(4), based on Ising interactions among spins, can be obtained through an approximation of the Hamiltonian function defined by non-local next-nearest-neighbor (NNN) elastic interactions. Observe that this is coherent with the results in Refs. [51,52], where the author shows that for small values of J the non-local interaction energy is, as in the case of the Ising model, proportional to the number of interfaces plus possible higher order boundary energy terms (that are neglected in this paper). This is an important point since some previous investigations on thermal effects for multi-stable lattices considered NNN interactions [54]. The resulting analysis was complex enough to prevent a fully analytical solution, whereas the scheme with Ising interactions proposed here represents an important step forward because the results are obtained in closed form.
Here we analyze separately the two cases of assigned force (Gibbs ensemble) and assigned displacement (Helmholtz ensemble). We remark that detailed computations are reported in the Appendix B.

Ising model within the Gibbs ensemble
The statistical mechanics in the case of assigned force f , within the Gibbs ensemble (see Fig.1, top panel), can be introduced by calculating the canonical partition function Due to the fact that it is more useful to consider dimensionless parameters, in Eq. (6) we substitute x i with the stretch λ i and introduce the dimensionless forcef = f ℓ/J and energyH = H/J together with the main non dimensional parameter of the paperβ = βJ. We end up with the following partition function expression The sums over {S i } are to be considered extended to the values +1 and −1 for each spin variable (i = 1, . . . , N ). Moreover, we have that where, without loss of generality, to avoid rigid motions, we have assumed x 0 = 0). Following the calculations in Appendix B, we obtain . This result is similar to the one obtained in Ref. [64], where however the elastic constants were considered equal (α = 1), and where a three-dimensional structure was studied to deal with polymeric cooperative systems. The knowledge of the Gibbs partition function allows us to calculate the expected value x N of the chain length (i.e., the average value of the last position x N , that from now on we rename x for simplicity of notation), the average number of unfolded units n u , and the average number of interfaces ι between folded and unfolded units. Since using the dimensionless parameters introduced before, we get wherex = x/ℓ and ∆Ẽ = ∆E/J. In Fig.3, the behavior of the mechanical quantities defined in Eqs. (10), (11), and (12) is shown by varying the values ofβ (different colors) and α (different rows). More precisely, each row of Fig.3 corresponds to different values of α (namely, 1/3, 1 and 3 in the first, second and third row, respectively), whileβ = J/K B T varies in each single plot (assuming the values 1/2, 1 and 3/2 for the blue, yellow and red curves, respectively). The first column shows the average chain length x , the second one the average number of unfolded units n u , and the third one the average number of interfaces ι . To begin, considering the force-extension figure in the case with identical wells α = 1, we can notice that at the higher dimensionless Ising parameter β = 3/2 (red curve) there is a force plateau that corresponds to a cooperative unfolding of the chain units, a typical behavior characterizing the Gibbs ensemble [51]. At the lower valueβ = 1/2 (blue curve), when interaction terms decrease, this plateau is less sharp due to the fact that thermal fluctuations induce the system to explore a larger part of the wiggly energy landscape, eventually smoothing the force plateau. Again for α = 1, the cooperative behavior is further confirmed by the expectation value of the number of unfolded units n u that shows a transition from 0 to N at the same threshold of the forceextension curve. In the plot of the average number of interfaces ι (always with α = 1), we can notice that increasing the parameterβ decreases the number of interfaces due to the fact that in a ferromagnetic scenario all the units tend to be in the same configuration (either folded or unfolded). Furthermore, in this plot, we can easily notice that the force thresholds that are responsible of the synchronized unfolding of the chain units are the same for different values ofβ. This is due to the fact that in the second row we have considered chain units that present the same elastic constant for both wells (α = 1). Consider now the first row (α = 1/3), when the second well is more compliant (softening regime). In this case, the fully unfolded configuration force-extension curve has a lower slope than the homogeneous folded one. The most interesting aspect here, is that the force plateau occurs at different thresholds depending on the value ofβ and, then, on the temperature T . This behavior is further confirmed by the shift of the peak in the average number of interfaces ι , corresponding to a force threshold depending on the value ofβ. In this plot, in fact, we can observe that the force threshold increases withβ, meaning that the force plateau decreases with the temperature. All these consideration holds even in the third row (α = 3), when the second state is stiffer (hardening regime), with the only difference that now the force plateau increases with the temperature.
While these results, representing the main physical effects of the proposed model, can already be qualitatively discussed, as described above, in the limiting cases of N → ∞ andβ ≫ 1 this dependence can be quantitatively described, so that a more detailed analysis is postponed later.

Ising model within the Helmholtz ensemble
Consider now the isometric loading condition, described by the Helmholtz ensemble (see Fig.1, bottom panel). In this case, the total elongation of the chain is fixed by assigning x N . As shown in Ref. [54], one may use an inverse Laplace transform to obtain the canonical partition function in the Helmholtz ensemble, starting from the Gibbs one given in Eq. (6). Using the change of variable f → −iω/β, we can write and we obtain (see Appendix B for the detailed calculation) where This partition function allows us to evaluate the expectation values f of the force conjugated to the assigned displacement, the average value of the number of unfolded units n u , and the expectation value of the number of interfaces ι . We have In Fig.4, we show the behavior of the system within the Helmholtz ensemble. As before, we have different values for α in different rows (namely, α = 1/3, 1 , 3 in the first, second and third row, respectively) and, in each single plot,β is variable (namely,β = 1/2, 1, 3/2 for the blue, yellow and red curves, respectively). The average values f , n u and ι are represented versusx in the first, second and third column of Fig.4. We can observe that the expectation values of the normalized force versus the applied extension always exhibits the typical saw-tooth path associated to a non-synchronized phase transition of the units. In the Helmholtz ensemble, in fact, due to the prescribed total elongation of the chain, we observe a sequential unfolding rather than a synchronized unfolding of the units, as seen previously within the Gibbs ensemble. In the plots of the average number of interfaces (third column), we can notice that at the highestβ value, there is only one single domain wall throughout all the unfolding process. This means that with high values of the Ising parameter (or similarly for lower value of the temperature) it is strongly disadvantaged to create more than one interface in the course of the sequential unfolding process. The important result to underline, as already shown in Ref. [54], is that when a large ferromagnetic Ising coefficient is considered (β ≫ 1), the system favors the generation of a single propagating interface, with a nucleation stress peak corresponding to the sudden transition of a chain fraction in the unfolded conformation. The asymmetry of the curves representing the number of interfaces and the dependence of the force plateau on the temperature (throughβ) are related once again to the presence of different elastic constants for the folded and unfolded wells (for the cases with α = 1). We remark that from a theoretical point of view, the unfolding of the last unit is characterized by a downward force peak that represents the coalescence of the folded phase to the unfolded one. Actually, in real experiments this down peak is typically not attained due to the presence of grips forbidding full propagation of the new phase in the terminal region. Moreover, whenβ is low (meaning that the thermal fluctuations are high compared to the Ising parameter J), we can observe a reduction of the force peaks resulting in a smoothing of the force-extension curve, as shown in the first column of the figure. The behavior of the first peak is largely analyzed in the following Sections. The steps observed in the curves representing the number of unfolded units n u versus the rising normalized extension are also smoothed with a decreasingβ. Similarly, the number of interfaces ι increases with lower values ofβ.
3 Thermodynamic limit with strongly ferromagnetic behavior in the Gibbs ensemble As anticipated, the plateau force observed in the force-extension curves, corresponding to the so called Maxwell stress in the purely mechanical case [48], sensibly depends onβ and, therefore, on temperature in both the cases with applied force (Gibbs) and with total fixed elongation (Helmholtz). Here, to obtain an analytical measure of such an important effect, we consider strongly Ising interactions, i.e. withβ = J kB T sufficiently large, and we study the system in the thermodynamic limit (N → ∞). An analogous limit was considered in Ref. [52] for a purely mechanical system. We remark that this limit, due to the complexity of the system, is performed here only for the Gibbs ensemble. In the following Sections, we will introduce particular hypothesis to get more analytical results also in the Helmholtz ensemble.
Let us consider first the thermodynamic limit for the Gibbs ensemble. Thus, from Eq.
where we exploited the propertyλ 1 >λ 2 (see Eq. (B13) for details). From Eq. (10) we find A direct evaluation of the derivative in Eq. (20) gives where c − and c + also depend onf and are defined in Eq. (B11) and ∆ in Eq. (B14) (see Appendix B). In Observe that here we are considering a fixed value of α = 1 3 with the second phase softer than the first one. On the right panel of Fig.5, instead, we can see that keepingβ constant and increasing α we induce an increase in the force plateau, proving that the Maxwell force strictly depends on the asymmetry of the elastic constants of the two energy wells, as already observed in the purely mechanical case when constant force plateaux depending on α are observed.
Interestingly, it is possible to link the expectation value of the number of unfolded units to the assigned force in the thermodynamic limit. Indeed, by using Eq. (11) for N → ∞, we get where again c − and c + depend onf and are defined in Eq. (B11). Thus we may write Eq.(21) in the form We notice that Eq. (23), derived from Eqs. (21) and (22), gives a direct physical interpretation of the transition process with the expected value of the elongation obtained as a convex combination of the stretch-force relations in the purely folded and purely unfolded phases, depending on their percentage. In particular, when n u = 0, all units are in the folded phase and the system follows the first branch as obtained without temperature and non local interaction effects. On the other hand, in the opposite extreme case with n u = N , all the bistable units are unfolded and the elongation is now given by x = χN + (N/αk)f . These two homogeneous regimes, where the energy is convex and temperature has no effects, are connected by a force plateau whose height and slope depend on both temperature and Ising coefficient J through β.
in Eq. (B14), it is possible to consider the approximation Taking into account the form of c ± (see Eqs. (B7) and (B11) in Appendix B), this condition corresponds tõ where In particular, one verifies that in the limitβ → +∞ (with δ = 0) the ratio in Eq. (25) goes to zero. Using this approximation, Eq. (21) shows the two following different regimes, corresponding to the two homogeneous folded and unfolded fractions, depending on the forcef through the parameters c − = c − (f ) and c + = c + (f ) Thus under the assumption of largeβ, we get a horizontal force plateau, as shown in Fig.6 (panel a), obtained by imposing c − = c + , or, by using Eq. (B7), given by the condition where δ has been defined in Eq. (27). In this way we obtain the main result of the paper, i.e. that for large enough interface energy J and for a large value of N , the system is characterized by a temperature dependent transition force. This is controlled by a new temperature dependent term in Eq. (29), not present in Eq.(1). Moreover, we can obtain an explicit analytic expression of the dimensionless Maxwell forcef M as follows where Due to the different elastic constants (α = 1), the Maxwell force obtained from Eq. (30) depends on the temperature T (throughβ). A simpler dependence can be obtained by expanding Eq. (30) up to the first order in the difference (1 − α)/α, measuring the difference from the case when the two wells are identical α = 1. We obtaiñ which is valid for | (1 − α)/α |≪ 1. When α = 1 we retrieve the expressioñ f = ∆Ẽ/(χ − 1), well known in literature [68,69,78,79], and discussed in the Introduction. On the other hand, if T approaches zero, Eq.(29) simplifies to Eq.(1). The dependence of the Maxwell force onβ is shown in Fig.6 (panel b). We can observe the opposite behaviors exhibited in the hardening case with α > 1 (solid blue curve with α = 1/3) and in the softening case with α < 1 (solid red curve with α = 3). We can also notice that whenβ → ∞ the Maxwell force readsf M ≃k(χ − 1) with a first order correction O(1/β). In particular, when α < 1, ifβ decreases (as the temperature increases), the Maxwell force decreases as well. The corresponding variation in the stress-strain diagram is shown in Fig.6 (panel c). In particular, from Eqs. (30) and (31), we deduce that when the value ofβ satisfies the equality log(1/α)/2∆Ẽ =β, the Maxwell force characterizing the plateau isf M = 0 (see the point B in Fig.6, panel b). Thus, at this valueβ =β 0 , the Maxwell plateau at zero applied force connects the two natural configurations with x /N ≃ 1 and x /N ≃ χ corresponding to the completely folded and completely unfolded states, respectively. In other word, the system undergoes a transition to the second state at zero applied force. It is important to observe that solutions with zero or even negative Maxwell forces have been considered in the case of phase transitions in shape memory alloys [120]. In this case a consistent definition of the reference configuration, assumed to be stable, should be temperature dependent, so that for positive (negative) Maxwell stress the reference configuration is the homogeneous folded (unfolded) phase configuration.
It is important to observe that, for this case with α < 1, the Maxwell stress can decrease until the limit valuef F , as defined in Eq.(3), because after this value the energy of the unfolded state is not defined (see the point A in Fig.6, panel b, and the floor plateau in Fig.6, panel c). We denoteβ * the value ofβ, for which we attain the plateau atf F . Interestingly, for α < 1, the equatioñ f M (β * ) =f F is solved bỹ In the case with α > 1, represented in Fig.6, panel b (red solid curve), we observe that the Maxwell stress always increases with the temperature (decreasingβ). In this situation, the Maxwell force is always positive. However, the plateau stops existing atβ * corresponding to the value of forcef R , defined in Eq.(2), when the energy of the first well ceases to be defined (see the point C in Fig.6, panel b, and the roof plateau in Fig.6, panel d). Interestingly, for α > 1, the equationf M (β * ) =f R is solved bỹ It is important to point out that, due to the use of the spin variables method, the temperature is always limited by the fact that the barrier between the two energy wells must always be sufficiently larger than K B T . In this Section, we obtained an analytical expression (δ = 0) linking the Maxwell force to the temperature in the thermodynamic limit (N → ∞) and with the assumption of strong ferromagnetic interaction (under the isotensional condition). This result is related to the difference of stiffness between the two energy wells (α = 1) and explains the mechanical response observed in several nanosystems [53,[89][90][91][92][93][94][95][96][97][98][99][100][101].
Interestingly enough, if for α > 1 (α < 1) we increase (decrease) the stress starting from the homogeneous folded (unfolded) state, the system keeps this configuration even for values of the stress for which this configuration is characterized by a higher elastic energy than the other homogeneous state, due to entropic effects. This counter-intuitive behavior, tending to stabilize the softer phase, is observable only in the case of different wells and was named entropic stabilization [120].

The zipper model
In this Section, we further extend the study of the non-local non-convex chain of bistable units under the assumption of strong ferromagnetic behavior (β ≫ 1). As previously discussed, this hypothesis leads to the existence of a single domain wall [52,54], i.e a single interface between the folded and unfolded regions that propagates continuously through the chain (see Fig.7). In this framework, we derive simplified analytical expressions for both the Gibbs and Helmholtz ensembles and we prove their equivalence in the thermodynamic limit. Accordingly, we consider a chain composed by N − ξ units in the folded state and the remaining ξ units in the unfolded one, being ξ a discrete variable assuming values in the set {0, 1, 2, . . . , N }, representing the position, in terms of chain units, of the moving domain wall between folded and unfolded regions. This discrete variable can vary depending on both mechanical and thermal effects and this simplified scheme is typically called zipper model [113][114][115][116].

Zipper model within the Gibbs ensemble
Let us evaluate the canonical partition function in the Gibbs ensemble by using the Hamiltonian in Eq. (37). By definition we have that can be evaluated by a Gaussian integration giving (39) Following the same reasoning as in Section 2, one may evaluate the expectation values of the mechanical macroscopic observables of the system, namely x , n u and ι (see Eqs.(10), (11) and (12)), where we use againx to indicatex N . In Fig.8, the behavior of the mechanical quantities x , n u and ι is shown by varying the values ofβ (different colors) and α (different rows). Each row of Fig.8 corresponds to a different constant α (namely, 1/3, 1 and 3 in the first, second and third row, respectively) whileβ = J/K B T varies within each single plot (β = 1, 2, 3 for the blue, yellow and red curves, respectively). The first column shows the average chain length x , the second one the average number of unfolded units n u and the third one the average number of interfaces ι . By looking at the force-extension curve and at the average number of unfolded units, one observes again the typical synchronized behavior of the Gibbs ensemble, where all the units unfold cooperatively. A difference with respect to the exact Ising model is that the number of interfaces is smaller, limited by the presence of the zipper assumption.
In particular, we may study x , evaluated through Eq.(10) and the partition function in Eq. (39). After some straightforward calculations, one obtains where δ is defined in Eq. (27) and (41) Discrete lattices with non-local interactions and non-convex energy We can now use the following sums N ξ=0 which allow us to obtain the explicit expressions Thus, Eq.(40) reads which provides the Gibbs force-extension relation in the approximations of single interface and large ferromagnetic interaction (zipper model). Observe that, thanks to this hypothesis, we are able to obtain for the zipper model exact explicit relations for the material response of the system also in the discrete case with finite arbitrary N . On the other hand, in the thermodynamic limit (N → ∞), we can confirm the results obtained in Section 3 for the complete Ising scheme. In particular, when N → ∞, the value of x /N depends on the sign of δ. If δ < 0 we have lim N →∞ and we get (in the limit N → ∞) Hence, a manipulation of Eq.(46) gives (in the limit N → ∞) representing the unfolded elastic branch. Thus, the two elastic branches for δ < 0 and δ > 0 are linked by a force plateau (see Fig.6, panel a), corresponding to Eq. δ = 0. This is exactly the same condition found in Eq. (29), proving that the results obtained within the zipper approximation in the thermodynamic limit coincide with the ones obtained in Sect.3 under the strong ferromagnetic assumption.

Zipper model within the Helmholtz ensemble
Let us then consider the Helmholtz ensemble under the zipper assumption. In this case, the partition function can be obtained as in Eq.(13), by the application of a Fourier transform. Indeed, using Eq.(39), we get The average macroscopic property of the zipper model within the Helmholtz ensemble can be found through Eqs. (16), (17) and (18). In Fig.9 we show f , n u and ι as function of the prescribed extensionx N for different values of the Ising coefficientβ and different elastic constants ratio α. We recover the same macroscopic behavior of the exact Ising model characterized by the sequential unfolding of the units and represented by the typical saw-tooth path of the force-extension curve. In particular, it is important to observe that the value of the number of interface remains constantly at 1, due to the zipper assumption, and it is 0 only when ξ = 0 or ξ = N , coherently with the model introduced. It is interesting to notice that in the force-extension curve, we can see the first upward force peak, corresponding of the nucleation of a new phase with unfolded units, and the last downward peak, representing the coalescence of the folded region to the unfolded one. Observing the sequential unfolding of the units in the force-extension plots of Fig.9, we can notice that there are three different behaviors depending on the elastic constants ratio: in the case with α = 1/3, the force needed to unfold a single unit decreases while the number of unfolded units increases; in the case with α = 1, except for the nucleation (first peak) and the coalescence (last peak) forces, all the other forces needed to unfold the units are equal; finally, in the case with α = 3, we have the opposite behavior to the case with α = 1/3, i.e. the unfolding forces increase with the number of unfolded units. This result is perfectly coherent with the previously introduced concept of entropic stabilization of the softer configuration. Moreover we remark that by considering higher number of elements (not reported for compactness), it is possible to observe an initial nucleation of a segment with more than one single element passing to the new unfolded state. This aspect is discussed in detail in the next section where we consider the system behavior for large values of N . In this case we can analyze the amplitude of the first peak and of the first nucleated segment analytically leading to a useful result for comparing the model with experiments and numerical simulations.

Stationary phase analysis of the zipper model within the Helmholtz ensemble
Here, we want to study the thermodynamic limit of the Helmholtz zipper model. The aim of this development is to obtain three important results: (i) a simplified expression of the force-extension relation under isometric condition valid for large values of N and forβ ≫ 1, (ii) an explicit expression for the amplitude of the first force peak, representing the nucleation of the unfolded phase, and (iii) a rigorous demonstration of the equivalence of the Gibbs and Helmholtz ensembles for the zipper model. Let us then introduce the average chain stretchλ = xN ℓN , prescribed to the chain under isometric condition. By using Eqs. (16) and (51), we obtain the force-extension relation in the form where we introduced the following functions we defined the rescaled global stiffness of the system (with 0 ≤ t ≤ 1) and we used the compact notations Let us now consider the behavior in the thermodynamical limit N → ∞. By following the approach suggested in Refs. [65,66] to obtain explicit analytical results, not only for the stress plateau but also for the stress peak, we may consider the Euler-MacLaurin (EM) approximation for a given function φ where higher order terms of the EM approximation would lead to more detailed, but analytically cumbersome results. Thus, from Eq.(52) we get To simplify the calculation, we may rewrite the integrals as where we introduced the phase fraction η = ξ/N and the auxiliary functions Now, we use the following general results [121,122]. Let then, we have 1. Suppose F is bounded and continuous on (a, b), and F (a)F (b) = 0. Suppose also that G is strictly monotone and differentiable and that F (a) G ′ (a) and F (b) both exist as finite reals, defined as limits if either endpoint is infinite. Assume also that the integral in Eq. (61) exists for all x > 0. Then, we have that 2. Suppose F is bounded and continuous on (a, b), that G has unique maximum at some c in the open interval (a, b), G is differentiable in some neighborhood of c, G ′′ (c) exists and is G ′′ (c) < 0, and that F (c) = 0. Then, we have that We use the results in Eqs. (62) and (63) to analyze the behavior of the integrals defined in Eqs. (57) and (58). First, by using Eq.(59), we obtain Solving ∂g ∂η = 0, we obtain the solution η 0 that represents the stationary point to be computed. Thus, following the definition in Eq.(60), we also introduce the valuef We observe that the equation ∂g ∂η = 0 can be written in terms off 0 as that coincides with Eq.(29) (i.e. with δ = 0). Thus, we identifyf 0 with the Maxwell forcef M . Now, we can use Eq.(63) (stationary phase theorem) to simplify Eqs. (57) and (58) only if 0 < η 0 < 1 (the stationary point must be within the integration interval). For η 0 = 0, we havef 0 = (λ − 1)k, and for η 0 = 1, we havef 0 = (λ − χ)αk, as given by Eq. (65). Hence, the interval 0 < η 0 < 1 is equivalent to 1 +f 0 k <λ < χ+f 0 αk , which corresponds to the force plateau region between the two elastic branches. Therefore, only in the plateau interval we have a stationary point and we can approximate Eqs. (57) and (58) with Eq. (63). The application of the stationary phase method is further justified by the relation g ′′ (η 0 ) < 0, simply proved by a direct evaluation Note that from now on, we use the notationf M for the quantityf 0 , to be consistent with previous Sections. We obtain from Eq.(56) the expression where This important analytical expression for the force-extension relation contains all the physical features describing the Helmholtz ensemble, as the first nucleation peak and the coalescence one. In particular, it can be used to evaluate these peaks as function of N and the temperature T . Moreover, considering that ∆g(0) < 0 and ∆g(1) < 0 (since g(η 0 ) is the maximum value of g(η) in the interval 0 < η < 1), we may prove that in the thermodynamic limit which is valid for 1 +f M k <λ < χ +f M αk , and meaning that the Maxwell force is the same for both the Helmholtz and Gibbs ensembles, proving their equivalence for the zipper model.
In order to conclude the analysis, we have to simplify Eq.(56) also for the external regions, i.e. forλ < 1 +f M k andλ > χ +f M αk . Hence, we suppose that the critical point η 0 is external to the interval (0, 1). We can have either η 0 < 0 or η 0 > 1. In these cases, the asymptotic behavior of the integrals in Eqs. (57) and (58) is described by Eq. (62). Then, Eq.(56) assumes the form where A = (αk) Now, by using Eq.(66), we determine the quantity g(1) − g(0) that eventually reads Accordingly, in the limit of N → ∞ we obtain a result representing the elastic branches in the external regions, corresponding to the fully folded (left) and fully unfolded (right) phases. This completes the proof of the equivalence of the ensembles in the thermodynamic limit for the zipper model. In Fig.10, we compare the force-extension curves obtained through the exact zipper model with the ones delivered by the stationary phase method. Observe that while this approximation, stated in Eqs. (68) and (73), is finer as N grows, our results show the possibility of a good approximation of just the nucleation peak already for small values of N (we used N = 10 in Fig.10). Also, the approximation correctly reproduces the last downwards peak, corresponding to the coalescence of the remaining folded regions to the fully unfolded phase. The most important result, quite evident in this framework, is the dependence of the Maxwell (propagation) force on the temperature, through the parameterβ. We observe that for the represented case with α < 1 (softening) asβ decreases (e.g. temperature increases), the value of the propagation force f M decreases, also affecting the height and the width of the nucleation and coalescence peaks with an entropic stabilization of the unfolded phase. In Fig.11, we observe that, when the number of units N increases and the thermodynamic limit is approached, the mechanical response of the system in the Helmholtz ensemble (see Eq.(68)) approaches the Gibbs behavior, where there is a force plateau without nucleation and coalescence force peaks. In particular, we notice that by means of the stationary phase method (when N → ∞) the typical saw-tooth path of the isometric condition disappears whereas only the first and the last peaks remain, as expected from a macroscopic point of view. We remark that, in experiments and simulations of phase transformations in nanowires or other nanostructures, it is possible to observe only the first nucleation peak, while the last one is rarely attained because of the breaking of the specimen or because of the boundary effects induced by the grips of the traction device.
Given the importance of the first nucleation force peak, we now look for an explicit expression of its magnitude. Consider Eq. (68) and observe that in the region of the first peak only the term e ∆g0N is relevant, being e ∆g1N negligible. Then, we can write a reduced form of Eq. (68), as follows where we introduced the change of variable y =λ−1−f M k . Now, let us consider as reference for the force its Maxwell value, so as to obtain We search for the stationary point of Eq. (79) by derivation, easily obtaining which is a Lambert equation of the form we w = s, solved by w = W 0 (s), where W 0 is the Lambert function (more precisely, W 0 is the principal branch of the Lambert function) [123,124]. Thus, we may writẽ Eventually, by using Eq. (81) in Eq.(79) we get wheref P is the value of the nucleation force peak. It is worth to remark that Eq. (82) gives the explicit dependence of the force peak on bothβ and the number of units of the chain N .
In the first row of Fig.12 we study the nucleation peak with respect to N (parametrized byβ) and with respect toβ (parametrized by N ). Similarly, in the second row, we show the plots off P −f M . The trend can be summarized by stating thatf P −f M increases for increasing values ofβ and decreases for increasing values of N . Moreover,f P converges tof M in the thermodynamic limit.
The consistency of the approaches presented in this work can be appreciated by observing Fig.13. Here, we compare the force-extension relation evaluated within the Helmholtz ensemble and obtained with (i) the complete Ising scheme discussed in Section 2.2, (ii) the simplified zipper model introduced in Section 4.2, and (iii) the stationary phase approximation given in Eq. (68). Moreover, we also show the amplitude of the first peak calculated by Eq. (82). This plot proves that, for a sufficiently large value of N and with the hypothesis of strong ferromagnetic interactions (large values ofβ), the three frameworks give the same results and hence the hypothesis introduced are consistent with the problem under investigation.

Applications to the tensile behavior of nanowires
While, as pointed out in the Introduction, our model can be applied to different multistable systems with softening or hardening effects (i.e. with different stiffness of the different phases) and non local interactions, we here focus our attention on the specific application to phase transition in metallic nanowires. We can show the ability of the model to predict the main physical phenomena observed during the phase transition at variable temperature. Single-crystalline metal nanowires with nanometric cross-section can exhibit a pseudo-elastic behavior characterized by very large eleongations which can be up to 50% of the original length [53,[90][91][92]. This behavior, which is exceptional compared to all other shape memory alloys, is due to a reversible lattice reorientation process with twin boundary propagation between two differently oriented facecentered cubic (FCC) crystalline structures. This behavior is typical of Copper (Cu) and Nickel (Ni), where we can identify the two following configurations: the original one named < 110 > / {111} with axis < 110 > and surfaces {111}, and the deformed one named < 001 > / {100} with axis < 001 > and surfaces {100}, as shown e.g. in Figs.2 and 3 of Ref. [53].
To compare the behavior with our one dimensional system, let s be the area on the cross-section pertaining to a single longitudinal chain of atoms in the nanowire crystal structure. Let then M the average number of atoms in the cross section so that the total area is given by S = M s. Similarly, we use ℓ to indicate roughly the lattice constant of the crystal structure and hence L = N ℓ is the total nanowire length, being N the number of atoms in the longitudinal direction. Accordingly, SL is the total volume of the system, M N is the total number of atoms, sℓ is the average volume pertaining to one atom, and one can observe also that s ≃ ℓ 2 and sℓ ≃ ℓ 3 .
It is important to remark that both the Maxwell force in Eqs. (29) and the force-extension relation in Eq.(68), have been obtained for a single chain of units (here atoms) and therefore the continuous parameters must be introduced as follows. We define the Young moduli of the folded (f ) and unfolded (u) phases, respectively. Accordingly, k 0 and h 0 are the elastic constants of the two crystals in the harmonic limit. Moreover, the jump energy between the two configurations can be written as ∆E = sℓ∆e where ∆e [J/m 3 ] is the energy density difference between the two crystal states. In addition, the Ising energy can be rewritten as J = Λ s, where Λ [J/m 2 ] is the surface energy density of the twin boundary separating original and deformed crystals. We remark that the parameters ∆e and Λ may depend on the total area S when elastic surface effects are relevant, as in the case of metallic nanowires. Further, the longitudinal (normal) stress is defined by σ = f /s, with f given in Eq. (68) if we consider the stationary phase approximation. The strain is defined asε = x−N ℓ N ℓ , and therefore we havē ε =λ − 1.
We use the introduced continuous quantities to rewrite previous expressions, for the purpose of mimicking the tensile behavior of one-dimensional structures. Indeed, Eqs. (29) and (27) can be rewritten as that is used to obtain the temperature dependent Maxwell stress σ M . We observe that, due to the hypothesis of homogeneity in the cross-section, the wire specific response coincides with that of a longitudinal chains of atoms. Correspondingly, the nucleation stress can be defined, using Eq.(82) as follows where C = CN 2 k e is, by using (71), given by Young modulus Value 122 GPa E f (T = 200K) 114 GPa E f (T = 300K) 105 GPa Eu (T = 100K) 10 GPa Eu (T = 200K) 11 GPa Eu(T = 300K) 12 GPa Table 1 Mechanical parameters used to reproduce the molecular dynamics simulations on a Cu nanowire of Ref. [53].
Moreover, by recalling the definition of β = 1/(K B T ), we introduce also the the parameterizations to obtain the two values of the characteristic energies. Here, we use the reference temperature T 0 = 300 K. In Ref. [53] a series of molecular dynamics simulations have been performed to obtain the stress-strain behavior of Cu nanowires at different temperatures and under quasi-static tensile deformation, which is compatible with our equilibrium statistical mechanics approach. The simulations are based on the embedded-atom-method interatomic potential for Cu [53], and the uniaxial displacement-controlled loading strategy is applied coherently with our Helmholtz ensemble. The stress-strain behaviors of this Cu nanowire (1.96 nm × 1.96 nm section) at temperatures of 100, 200, and 300 K is compared with our model in Fig.14, left panel. The theoretical results are based on the stationary phase approximation stated in Eq. (68). As one can see, the model is able to reproduce the magnitude of both the nucleation stress and the Maxwell stress plateau. In Table 1, one can find the temperature dependent Young moduli we adopted to reproduce the simulations of Ref. [53]. Moreover, we adopt the following parameters: N = 206, ℓ = 0.14 nm, s = ℓ 2 and χ = 4.4. By introducing the parameterizations in Eq.(87), we obtained the adimensional values γ = 7.7 andβ 0 = 6.8. We remark that sinceβ 0 , representing our non dimensional parameterβ for the higher value of the temperature, is sufficiently larger than 1, the use of the zipper model, simplified by the stationary phase method, is justified. Moreover, the propagation of a single twin boundary is confirmed by the simulations in Ref. [53].
We propose also a second example of application of our model to ultrathin semiconductor Germanium (Ge) nanowires that have been recently studied due to their peculiar properties. Here, we compare our results with molecular dynamics simulations of the helix Ge nanowire with radius of 1.8Å(see Fig.1 in Ref. [96]). The stress-strain curves has been obtained for temperatures of 100, 200, 300, and 400 K [96]. The simulations for germanium were based on the Stillinger-Weber potential, coupled with the Nosé-Hoover thermostat to Young modulus Value 561 GPa E f (T = 300K) 552 GPa 83 GPa Table 2 Mechanical parameters used to reproduce the molecular dynamics simulations on a Ge ultra thin nanowire of Ref. [96].
impose the system temperature. In Fig. 14, right panel, we show the comparison of the simulations with our model, and we observe that the nucleation and Maxwell stresses are in fairly good agreement. Unfortunately, the numerical stress plateaus are rather noisy because they correspond to a single simulation (without averaging) and the amorphous phase is strongly fluctuating. Anyway, the temperature dependent behavior of such plateaus is clearly visible. The helical structure is perfectly ordered initially (see Fig.1 in Ref. [96]), and shows a transition to an amorphous structure (see Fig.5 in Ref. [96]) following a transition process that begins with a stress peak, namely the nucleation stress. After this nucleation, the unfolded region evolves 'smoothly' into the one-atom chain structure. In our theoretical formulation, we used the parameters listed in Table 2 corresponding to the data in Ref. [96]. Moreover, we adopted the parameters: N = 5562, ℓ = 0.208 nm, s = ℓ 2 and χ = 3.36. In particular, by introducing the parameterizations in Eq.(87), we obtained the adimensional values γ = 2.3 andβ 0 = 25.3 (orβ = 19 in the case of highest temperature T = 400K) so that the assumption of the zipper model is widely justified. To sum up our results, we considered here systems where two possible configurations exist. In the first application they correspond to two different crystal structures whereas in the second example to a regular helical structure and an irregular amorphous structure. It is worth noticing that the model based on statistical mechanics here developed is able to reproduce the complex behavior of these systems. It is important to highlight that in our model there is no crystallographic or morphological information of the three-dimensional structures, whereas the informations of the microstructure are embedded in the elastic and energetic properties of the system. Thus the model captures the main physical behavior in term of nucleation and propagation stress and their dependence on temperature, without exploiting the microstructural details. Nevertheless, the generality of the model allows its application to different systems ranging from material science to cooperative biological structures. Of course, if it is necessary to take into account crystallographic details, structural anisotropies, twin-boundary geometries and other morphological features, then it is necessary to turn to another class of models that are specifically adapted to the problem at hand, e.g. the model in Ref. [53] for the pseudo-elasticity of metallic nanowires. On the other hand, the possibility of deducing explicit equations of the temperature dependence of the nucleation and propagation stress, depending on parameters with a clear physical interpretation, represents in our opinion the main advantage of the proposed model both in the perspectives of interpreting the transition behavior of biological materials and phase transition alloys and in the field of new material design.

Conclusions
In this work, we elaborated some models, with different level of complexity, to describe the temperature-dependent behavior of one-dimensional non-local non-convex systems. The paradigmatic system under investigation is composed of a sequence of units, exhibiting a bistable behavior described by a two-state potential energy, and which are in interaction with each other. On the one hand, this scheme is able to represent the specific feature of several macromolecules of biological origin (mainly proteins). In this case, the non-convexity (bistability) describes the possible switching of each domain of the macromolecule between its folded and unfolded states. In addition, the non-locality, introduced by means of the Ising scheme, depicts the cooperativity observed in most of the biological structures. From the other hand, the proposed model is equally able to represent the phase transformations in materials. The two states of the bistable behavior represent in this case two microstructures of the solid material, corresponding e.g. to two crystal structures or to austenitic and martensitic phases. This enables the study of pseudo-elasticity and shape memory effects in solid systems such as whiskers, nanowires or nanocomposites. For these systems, the Ising interaction scheme reproduces the energy cost of creating an interface between the two different states of matter due to the absence of kinematic compatibility. Our model then describes in a prototypical way the origin of the nucleation, propagation and, possibly, coalescence phenomena and their possible dependence on temperature when low dimensions are considered as in nanowires so that entropic contributions can compete with bulk elastic energy terms.
Indeed, in the proposed approach, the thermal effects are carefully considered by introducing the statistical mechanics analysis of the problem, based on the method of the spin variables. This allows us to develop a complete thermodynamic framework of the transformation processes in systems with non-locality and non-convexity. The first proposed approach is the most refined and implements the complexity of the system in its entire generality, without simplifying assumptions. The system is constituted by a chain of bistable elements interacting through an Ising scheme. We can study both the Gibbs (isotensional) and the Helmholtz (isometric) ensembles of the statistical mechanics by evaluating the pertinent partition functions in closed form. These results are valid for both the ferromagnetic (positive Ising interaction) and antiferromagnetic cases (negative Ising interaction). However, in this work, we mainly focus on the ferromagnetic behavior. In this model, each unit can freely assume one of the two admitted states, depending on the system temperature and on the mechanical actions applied to the system. The number of interfaces between folded and unfolded regions is therefore free to vary and represents a measure of the entropic effects on the system. Hence, we can observe an increase in the number of interfaces with increasing temperature for a ferromagnetic system. The important point is that for strongly ferromagnetic systems (with high energetic cost of interface generation), which are the most common in practical applications, only one interface is observed between the two different states and it propagates along the chain when the traction or elongation are applied to the system. In addition to the interface propagation phenomenon, in the isometric Helmholtz case, we can also observe an initial upward peak of force, representing the nucleation of the new unfolded phase, and a final downward peak of force, representing the coalescence of the folded phase into the unfolded one. An important feature of this system is that the force plateau, describing the interface propagation within the chain is in general temperature dependent for both isotensional and isometric boundary conditions. The origin of this dependence lies in the difference between the elastic constants of the unfolded and folded phases. In fact, when these elastic constants are equal, the force plateau becomes independent of temperature and known results from the literature are retrieved.
Since the thermodynamic limit analysis is extremely complex for the exact Ising model (in particular, for the isometric case), and since the important case for applications is the strongly ferromagnetic one, we have introduced the zipper model where we consider only one interface propagating along the system. This assumption makes the calculation of the partition function much simpler in both the Gibbs and the Helmholtz ensembles. The thermodynamic limit for the Gibbs ensemble can be analyzed straightforwardly and the results are perfectly coherent with the exact Ising model studied previously. In particular, the temperature dependent behavior of the force plateau is confirmed also by means of this zipper approach.
The investigation of the thermodynamic limit concerning the zipper model under isometric condition, i.e. within the Helmholtz ensemble, represents the third approach of this work. This analysis is based on the stationary phase method and it is useful to obtain three important results: (i) a simplified expression of the force-extension relation under isometric condition valid for very long and strongly ferromagnetic chains (which perfectly describes the nucleation, propagation and coalescence phenomena), (ii) an explicit expression for the amplitude of the first force peak (based on the Lambert function), representing the nucleation of unfolded units, and (iii) a rigorous demonstration of the equivalence of the Gibbs and Helmholtz ensembles for the zipper model (nucleation and coalescence peaks disappear as length increases). These results, based on the stationary phase approximation, properly describe the transformation processes in one-dimensional objects (macromolecules or solid nanowires) in terms of thermal fluctuations and mechanical actions applied to the system. Some examples of quantitative comparisons are shown for the microstructural evolution in metallic and semiconductor nanowires. • All authors certify that they have no affiliations with or involvement in any organization or entity with any financial interest or non-financial interest in the subject matter or materials discussed in this manuscript. • All data generated or analysed during this study are included in this published article • All authors contributed to the study conception and design. Numerical procedures were performed by AC and LB. Theoretical developments were coordinated by GF, GP and SG. All authors read and approved the final manuscript.
Appendix A Non-local behavior: relation between the next to nearest neighbor (NNN) interaction strategy and the Ising scheme The non-local interactions in discrete elastic chains are typically introduced through next to nearest neighbor (NNN) elastic elements [51,54]. Under specific hypotheses, we prove here that the scheme with the NNN elements can be reconducted to a typical Ising model, widely adopted in previous studies on similar topics [63,64]. To begin with, we consider a chain of bistable units with additional NNN linear springs to model non-local effects [54]. The Hamiltonian reads where R is the elastic constant of the NNN springs, and the other quantities are defined in Section 2. In particular, we name here H 0 the energy corresponding to the bistable nearest neighbor (NN) elements and H I the energy of the NNN linear springs, such that H = H 0 + H I . We assume that the equilibrium length of the NNN elements is fixed at L(−1)+L(+1) = (1+χ)ℓ, and we introduce two different behaviors: an antiferromagnetic one when R > 0, and a ferromagnetic one when R < 0. In the former case (R > 0), two adjacent units entail a lower energy if they are in two different states (folded and unfolded) while, in the latter one (R < 0), two adjacent units result in a lower energy if they are in the same state (either both folded or both unfolded), as shown in Fig.A1. Supposing that the elastic constants k and αk are sufficiently large, the Hamiltonian H I can be simplified by assuming that the lengths of the units can be approximated with the equilibrium lengths of the explored wells. Then, by using the relations λ i+1 ≃ λ 0 (S i+1 ) and λ i ≃ λ 0 (S i ), we get Thus, the interaction Hamiltonian becomes where This proves that we can approximate the behavior of the chain with the NNN elements by means of a classical Ising chain. We observe that when J > 0 (R < 0) we are in a ferromagnetic case, and when J < 0 (R > 0) in the antiferromagnetic one. For this reason, in this work we adopted the overall Hamiltonian given by Eq.(4), where we neglected the irrelevant constant in H I .

Appendix B Helmholtz and Gibbs ensembles in the Ising Model
In this appendix we show the details of the calculation of the Gibbs and Helmholtz partition functions.

B.1 Gibbs ensemble
For the Gibbs ensemble, using Eq.(5), the partition function can be evaluated as wheref is the dimensionless force applied to the last unit of the chain. Each sum on S i (i = 1, . . . , N ) must be interpreted as a sum over the values +1 and −1. By a Gaussian integration we obtain so that we obtain It order to explicitly evaluate the summation we can use the transfer matrix method [81]. We obtain where we have defined the transfer matrix T and the vector w (taking care of the boundary conditions) as follows with, see Eq.(B7), c + c(+1), c − c(−1). (B11) By using the standard matrix functions theory [125,126], we have where I is the 2 × 2 identity matrix andλ 1,2 are the eigenvalues of T , namelŷ In Eq.(B13),λ 1 (λ 2 ) corresponds to the + (−) sign and we have also defined By substitutingλ 1,2 into Eq.(B9) and Eq.(B12), we get the partition function given in Eq. (8).

B.2 Helmholtz ensemble
We consider here the case with fixed x N (isometric condition), described by the Helmholtz ensemble, and evaluate the canonical partition function Z H (x N ). The partition functions in the Gibbs and Helmholtz ensembles are linked by a Laplace transform [82] Z Thus, one can write Z H (x N ), using the change of variable f → −iω/β, as an inverse Fourier transform in the complex plane To simplify the notation and perform the calculation we define δ 1,2 = e −βλ 1,2 = Thus, the Gibbs partition function can be written as By using the Newton binomial rule to expand the powers δ N 1 and δ N 2 , we find Here, we can separate the even and odds terms as follows where the square brackets in the sums stand for the floor function defined as [x] = max {n ∈ Z|n ≤ x}. Then, we get We can further develop the powers ∆ k and (c + + c − ) N −2j through the Newton binomial rule, obtaining Finally, Eq. (B16) can be evaluated, yielding the canonical partition function within the Helmholtz ensemble given in Eqs. (14) and (15).