Statistical mechanics of rate-independent stick-slip on a corrugated surface composed of parabolic wells

The stick-slip phenomenon, at the basis of friction, is crucial for several applications ranging from nanotechnology and biophysics to mechanics and geology. Deep understanding of friction mechanisms and, in particular, the methodologies for its reduction must be sought in its nanoscopic nature, where atomic interactions and stick-slip processes play a crucial role. At this scale, thermal fluctuations clearly have a major effect on the physics of the problem. Hence, we develop here a theory for rate-independent stick-slip, based on equilibrium statistical mechanics. In particular, we introduce suitably modified Prandtl–Tomlinson and Frenkel–Kontorova models in order to study the system with one particle and the chain with N particles, respectively. The adopted corrugated substrate is composed of a sequence of quadratic wells. Interestingly, the calculation of corresponding partition functions shows a conceptual link with the theory of Jacobi and Riemann theta functions, allowing an efficient determination of the average static frictional force and other relevant quantities. We show some applications including the study of structural lubricity and thermolubricity.

From the historical point of view, the empirical laws governing the frictional force occurring between a slider and a substrate have been investigated and established by Leonardo da Vinci (1452-1519), Guillaume Amontons (1663-1705), and Charles Augustin de Coulomb (1736-1806) [29]: (i) frictional force is independent of the apparent contact area between slider and substrate, (ii) the frictional force is proportional to the loading normal force, and (iii) kinetic friction is independent of the sliding velocity and is typically smaller than static friction. These macroscopic observations are purely empirical and ultimately derive from the complex physical processes at the interface between slider and substrate, involving the attachment and Communicated by Andreas Öchsner.
S. Giordano (B) Univ. Lille, CNRS, Centrale Lille, Univ. Polytechnique Hauts-de-France, UMR 8520 -IEMN -Institut d'Electronique de Microélectronique et de Nanotechnologie, Lille 59000, France E-mail: stefano.giordano@univ-lille.fr detachment of asperities and the resulting variation in the real contact area [30,31]. The nature of friction at the nanoscale has begun to be studied more systematically with the advent of the atomic force microscope [32], the surface force apparatus [33], and the development of efficient molecular dynamics simulations and multiscale models [34][35][36][37]. These approaches to nanotribology have enabled the experimental and numerical observation of the characteristic stick-slip motion of the slider atoms, which interact with the periodically corrugated potential energy representing the atoms of the substrate. This result allowed the validation of an older conceptual model, which is usually referred to as the Prandtl-Tomlinson model [38,39]. This is the most efficient model currently known for describing stick-slip motion or nanoscale friction and, for this reason, has been largely investigated and compared with experimental and numerical results [40,41]. It consists in a single point mass moving on a one-dimensional periodic potential (typically sinusoidal) pulled by a linear spring characterized by its elastic constant. The analysis of this model revealed different regimes of temperature and velocity, and these theoretical achievements explained several atomic force microscope experiments [42][43][44][45][46][47][48][49][50][51][52]. The most natural generalization of the Prandtl-Tomlinson model considers a one-dimensional elastic chain of interacting atoms (or particles) moving over the periodic potential mimicking the substrate. This scheme is known as Frenkel-Kontorova model and perfectly describes the mutual sliding of two crystalline interfaces. It has been first introduced by Dehlinger [53] and then adopted by Kontorova and Frenkel to study the dislocations motion and explain the plastic deformation in crystals [54][55][56]. The Frenkel-Kontorova model is able to describe a large range of nonlinear phenomena, such as solitons, chaos, friction, dislocations, incommensurate phases, and glass-like systems. Because of its simplicity and universality, it has been widely studied, eventually becoming one of the most popular models used in low-dimensional nonlinear physics [57][58][59]. We remark that the Frenkel-Kontorova model is sometimes named discrete sine-Gordon model, since the continuous limit of the Frenkel-Kontorova model exactly yields the sine-Gordon differential equation. Interestingly, this equation was originally obtained by Edmond Bour in 1862 while conducting researches on surfaces of constant negative curvature [60] and only rediscovered by Frenkel and Kontorova in 1939 in their study on dislocations [54][55][56]. The tribological importance of the Frenkel-Kontorova model and of the sine-Gordon equation is that their dynamics is characterized by topological solitons named kinks (in case of the chain compression) and anti-kinks (in the case of the chain extension), which move in the chain easily than particles, controlling the overall chain motion. Indeed, the shift of the whole chain on the left or the the right by one period of the corrugated substrate corresponds to the propagation of a kink (or anti-kink) along the chain. The energy barrier to be crossed to generate the kink movement is referred to as Peierls-Nabarro barrier (introduced in the study of dislocations motion) and is smaller than the amplitude of energetic oscillations of the corrugated substrate potential [61][62][63]. In this context, Aubry proved an important transition for an infinite incommensurate chain between a pinned state and an unpinned (frictionless) state as its elastic constant increases [64][65][66] (see Sect. 3 for details). This transition has been much studied in the past years, generalized to chains of finite length and observed experimentally [67][68][69][70][71][72][73][74][75].
Although the absence of friction can be dangerous in some cases, facilitating for example slipping on ice [76] or a banana peel [77], a great scientific challenge consists in reducing its value in order to limit energy consumption in technological process [16]. In this regard, the superlubricity phenomenon, corresponding to a state of very low friction between two surfaces, has been observed and thoroughly investigated [78,79]. This very special condition between surfaces can be achieved by various techniques, the most important of which are the structural superlubricity and the thermolubricity. In a structural superlubric state, the two interacting surfaces have different or incommensurate lattices whose mismatch leads to a nearly vanishing frictional force [80][81][82]. This effect has been experimentally observed by scanning tunneling microscopy for W(011) tip on Si(001) surfaces [83], by frictional force microscopy for metal tips over graphite [84,85] and by noncontact atomic force microscopy for graphene nanoribbons on gold surfaces [86], just to name a few. These results have been complemented by computational studies [87,88]. The effect of temperature on friction and the thermolubricity phenomenon are well described by the Prandtl-Tomlinson minimalistic model introduced above [38,39]. Actually, in the early pioneering article by Prandtl [38], one already finds that friction decreases with temperature since thermal energy fosters the crossing of energy barriers and eventually enables slip [42][43][44][45][46][47][48][49][50][51][52]. Once again, experimental evidence for this principle has been obtained for various materials using the atomic force microscope [89,90].
In order to have a comprehensive view of the effect of temperature on friction, we develop here a theory for the rate-independent stick-slip phenomenon, based on equilibrium statistical mechanics. To do this, we implement the so-called spin variable approach useful to deal with nonconvex potential energies [91]. This technique has already been successfully applied to different problems such as the physics of muscles [92,93], the folding of macromolecules [94][95][96][97], the adhesion processes [98,99], and the phase transformations in solids [100,101]. We remark that this approach is complementary to the methodologies commonly used to study the behavior of physical systems with multiple stable and metastable states [102][103][104][105]. Here, we modify the Prandtl-Tomlinson and Frenkel-Kontorova models by substituting the sinusoidal energy profile of the substrate with a sequence of quadratic potentials, each of which represents a substrate well. It is important to remark that the sinusoidal profile and the sequence of quadratic wells have a qualitatively similar behavior due to their periodicity but different features from the quantitative point of view because of the quite different shape. In this sense, the quadratic substrate should not be seen as an approximation of the sinusoidal one, but rather as an independent model with application, e.g., to frictional motion of defects in solids or sliding of macromolecules on biological structures [12,28]. Our approach allows a simpler calculation of the partition functions for both (modified) Prandtl-Tomlinson and Frenkel-Kontorova systems. The spin variables approach requires the introduction of a discrete variable (spin) for each particle of the system, useful to identify the substrate well locally explored by the particle itself. The application of this methodology to the (modified) Prandtl-Tomlinson and Frenkel-Kontorova models leads to explicit partition functions written in terms of Jacobi or Riemann theta functions. This allows an efficient calculation of all thermodynamic variables, including the average frictional force (when the energy barrier between the wells is sufficiently larger than K B T ). Eventually, this force depends on temperature and other geometrical and physical variables, thus explaining possible structural and thermal lubricity phenomena. In particular, we proved the existence for finite chains of specific commensurate values of the lattice constant ratio between slider and substrate, which are able to sensibly reduce the friction. We finally studied the effect of the temperature on this friction reduction. The models introduced, although rather simple, offer the future possibility of taking into account more complex phenomena not accessible to current theories. We refer for example to periodic but complex energy profiles and corrugated but deformable substrates where the slider affects the shape of the substrate itself (this is relevant in soft matter applications).
The paper is structured as follows. In Sect. 2, we develop the theory for the modified Prandtl-Tomlinson model, and we discuss the results in Sect. 3. Then, in Sect. 4, we elaborate the statistical mechanics for the modified Frenkel-Kontorova model, and we show some applications in Sect. 5. The conclusions (Sect. 6) and a mathematical "Appendix" close the paper.

One-dimensional stick-slip model
We introduce a simple model concerning the one-dimensional stick-slip behavior of a single particle in contact with a substrate described by a periodic potential energy. In particular, we propose **here a modification of the classical Prandtl-Tomlinson model [38][39][40][41], able to consider the effect of thermal fluctuations on the static system behavior. Indeed, these fluctuations can significantly modify the friction experienced by the particle deposited on the substrate. While the interplay between temperature and velocity in the Prandtl-Tomlinson model has been largely studied by means of several different approaches [42][43][44][45][46][47][48][49][50][51], the static stick-slip behavior is directly investigated here through statistical mechanics. In our approach, we assume to be close to the thermodynamic equilibrium, and therefore, we consider rate-independent processes observed with the slider at rest or at sufficiently low speed. The characteristic sinusoidal profile of the Prandtl-Tomlinson substrate potential is substituted here by a periodic sequence of quadratic wells with elastic constant k 0 and spacing d (see Fig. 1). The energy barrier between adjacent wells is given by W = k 0 d 2 /8 (it must be much larger than K B T , as discussed below). Of course, we do not consider the natural transition rates between adjacent wells (Kramers rate theory [106]) and the interaction of this dynamical processes with the (low) speed of the traction device. These effects can be taken into account only by means of the out of equilibrium statistical mechanics, and this analysis has not been performed here. In order to identify the well occupied by the particle, we have to introduce a new discrete variable n ∈ Z within the phase space of the system, which represents the number corresponding to the concerned substrate site. It is important to underline that the energy barrier of the corrugated substrate (quantified by W or k 0 ) has been directly related to the normal load applied between slider and substrate [52]. This gives still more importance to the Prandtl-Tomlinson model, whose features can be easily put in relation to the experimental parameters.
In our model, the manipulation of the particle is performed by a first linear spring (constant k b ) tethering the particle at the origin of the x-axis and a second linear spring (constant k) pulling the particle at position (see Fig. 1). These springs represent the device used to investigate the stick-slip behavior. Of course, if we want to simulate a non-constrained particle interacting with the substrate, we can let k b = 0, thus eliminating the bonding at the origin. Moreover, must be considered slowly variable to allow the particle to explore quasi-statically the substrate. These premises enable us to write the potential energy of the system as with x ∈ R and n ∈ Z. We suppose that the system is in contact with a thermal reservoir at temperature T and we denote the Boltzmann constant by K B . Hence, to introduce the partition function of the canonical distribution, we have to sum over the discrete variable n and to integrate over the continuous one x. Importantly, since we integrate the variable x for any potential well from −∞ to +∞, we have to consider an energy barrier W sufficiently larger than the thermal energy K B T . So doing, the effects generated by the overlapping of quadratic potentials are negligible. For example, if we consider a silicon cantilever sliding on a NaCl single crystal, we can observe that a normal load in the range between 0 and 6 nN corresponds to an energy barrier W between 0 and 0.6 eV [52]. Then, it is not difficult to prove that the relation W K B T (we can say W = 10K B T ) is verified when W > 0.3 eV, which means a normal load larger than 3 nN. This proves the applicability of our theory to real cases. In addition, we work within the Helmholtz ensemble of statistical mechanics since the distance is prescribed [107,108]. These remarks lead to the expression where we use Eq. (1). The integral can be now performed by means of the classical gaussian integration which is able to deliver the following result. To write it in compact form, we define the parameters where k t = k 0 + k b + k. We therefore observe that the partition function can be written as follows where we have neglected the noninfluential multiplicative constant and being ϑ 3 (z, τ ) the third Jacobi theta function defined below Here z = −iξ d ∈ C is the argument, and τ = iηd 2 /π ∈ C is the so-called lattice parameter satisfying the condition m(τ ) > 0 [109][110][111][112][113]. Sometimes, also the nome q = e iπτ = e −ηd 2 is introduced with the assumption |q| < 1, assuring the series convergence. The introduced function satisfies the two relations stating that the function is completely determined in the entire complex plane by the values it assumes in the parallelogram identified by the four points z 0 , z 0 +πτ, z 0 +π +πτ and z 0 +π ∀z 0 ∈ C (the fundamental domain) [109][110][111][112][113]. There is a wide variety of notations for the theta functions. Here, rather than the original Jacobi notation, we prefer the more modern notation introduced in Refs. [109][110][111][112][113]. Interestingly, theta functions, historically used for the study of elliptic functions and various problems in number theory, are also able to completely describe the statistical mechanics of rate-independent friction. From the thermodynamics of the system, we know that the average value of the force necessary to maintain the position of the device is given by By using the following relation giving the logarithmic derivative of the third theta function [109][110][111][112][113] where ϑ 3 = ∂ϑ 3 ∂z , q = e iπτ , and which is valid for | m(z)| < π 2 m(τ ), we easily obtain with validity in the range |ξ | < η/2. For the simplest case in which the particle is not bound to the origin of the reference system (k b = 0), we obtain with η = 1 2K B T k 0 k k 0 +k and | | < d/2. Concerning the average value of the discrete variable n, directly from Eq. (5), we can write n = +∞ n=−∞ ne −ηd 2 n 2 e 2ξ d n +∞ n=−∞ e −ηd 2 n 2 e 2ξ d n , (13) or, equivalently n = 1 2i Using again Eq. (10), we get which is valid if |ξ | < η/2. If k b = 0, this result can be easily simplified and gives Although correct and rather interesting, Eqs. (11)- (12) and (15)- (16) are not very useful to give a physical interpretation of the phenomenon under study. To get a new form for these solutions, we can introduce the following Jacobi functional identity for the theta function [112,113] where the square root is to be interpreted as the principal value (see "Appendix" for details). If we apply this identity to our partition function given in Eq. (5), we eventually obtain which is a new form completely equivalent to the previous one. From Eq. (9), by using Eq. (10) we can obtain a new expression for the static friction force as follows After the use of the Jacobi functional identity, this expression, contrary to Eq. (11), is always convergent and is composed of a linear term in added to a Fourier series in the same variable . The convergence everywhere of the new form in Eq. (19) depends on the fact that we are now dealing with a Fourier series, whereas before in Eq. (11) we had an arbitrary function series with limited convergence. The first term represents the effect of the bonding spring k b , whereas the second one comes from the periodic behavior of the corrugated substrate. The physical interpretation of Eq. (19) is therefore evident and transparent. When k b = 0, the linear term vanishes and we obtain which is the periodic rate-independent friction force. This expression allows a simple analysis of the low temperature regime. Indeed, by using the limit lim where we can recognize the Fourier series of a sawtooth wave thus obtaining lim T →0 F s = k 0 k /(k 0 + k), for any such that −1/2 < /d < 1/2. This simple result exactly corresponds to the purely mechanical force (without temperature effects) necessary to maintain the particle at position when it is attracted by the well centered at x = 0 (elastic constant k 0 ) and pulled by the spring with constant k. The theory is therefore perfectly consistent with the pure mechanics at T = 0. We can now use Eq. (14) combined with Eq. (18) to find the new expression for the average value n of the discrete or spin variable n. The result is which is again always convergent and composed of a linear term and a periodic one. When k b = 0, we finally obtain where the linear term now remains present because it represents the increase in n during the sliding, i.e. the sequential transitions between different wells of the substrate potential. In conclusion, the latest results obtained through the Jacobi functional relation are represented by Fourier series that always converge, converge quickly, and have an immediate physical interpretation.

Results for the modified Prandtl-Tomlinson model
All obtained expressions can be numerically implemented and some examples are found in Fig. 2. In panels (a) and (b) of Fig. 2 one can find the behavior of the static force F s /(k 0 d) and the spin variable n versus /d for a system with k b /k 0 = 0 (without bonding spring). It means that the device slowly spans some periods of the corrugated substrate and the main physical quantities are calculated for different temperatures. By observing the energy profile of the substrate, we see that the transitions of the particle between different wells occur at /d = 1/2, 3/2, ..., and so on. Then, in correspondence to these positions, we have peaks in the force [see panel (a) of Fig. 2] and switches between two adjacent integers in the average spin variable [see panel (b) of Fig. 2]. This is coherent with the Fourier series in Eqs. (20) and (24). Indeed, since k b /k 0 = 0, the curve representing F s /(k 0 d) in panel (a) is periodic in /d with period 1. In panels (c) and (d) of Fig. 2 we can find the behavior of the same quantities F s /(k 0 d) and n versus /d for a system with k b /k 0 = 0.2 (with a spring tethering the particle at the axis origin). In this case, the periodicity of the friction force is evidently lost, since a larger device force will be necessary moving away from the origin of the axis, as described by Eq. (19), with k b > 0. It is important to remark that, in both panels (a) and (c) of Fig. 2, we observe a reduction in the friction force with an increasing temperature, clearly explaining the thermolubricity phenomenon [78,79]. In practice, the increase in temperature promotes thermal fluctuations and thus the wider exploration of substrate energy basins. This then allows energy barriers between wells to be more easily overcome and ultimately reduces friction. We also note that in both panels (b) and (d) the average value n exhibits a sharp behavior with net transitions between the adjacent wells for low values of the temperature and a more blurred response for high temperatures. This is consistent with Eqs. (23) and (24).
To conclude, in panel (e), we show the behavior of the friction in terms of the temperature and of the elastic constants ratio k/k 0 . To do this, we quantify the static friction force by means of the maximum value max { F s /(k 0 d), 0 < /d ≤ 1} and we represent it as function of K B T /(k 0 d 2 ) and k/k 0 . It is important to recognize that this maximum force value corresponds exactly to the force that the device must be able to apply to allow the particle to advance from one period to another of the substrate corrugated potential. This concept is in close analogy with the Peierls-Nabarro stress introduced for the study of the onset of movement of dislocations [61][62][63]. We observe that the maximum value of F s /(k 0 d) in the period is decreasing with an increasing adimensional temperature K B T /(k 0 d 2 ) and is increasing with an increasing ratio k/k 0 . Importantly, we remark that the validity of the results is in the range with W sufficiently larger than K B T . In Fig. 2, we explored the behavior of the system up to the limit W K B T in order to better show the trend of the obtained results. The plot in panel (e) is useful to study the behavior of the system when k is variable and k 0 is constant. But we can show that the dual behavior can be deduced directly from the previous one. Indeed, if we define the functions we can rewrite Eqs. (12) and (20) as which represents the plot shown in panel (e) of Fig. 2. Similarly, from Eqs. (12) and (20) we can also easily obtain the dual result which can be used when we consider k 0 variable and k constant. We remark that the function G of Eq. (25) has been obtained with the original partition function in Eq. (5), while the function F of Eq. (26) with the partition function in Eq. (18), modified by means of the Jacobi functional relation. Since the behavior of Eqs. (27) and (28) is exactly the same, we deduce that in the panel (e) of Fig. 2 we can substitute F s /(k 0 d), k/k 0 and K B T /(k 0 d 2 ) with F s /(kd), k 0 /k and K B T /(kd 2 ) and the curves remain unchanged. Therefore, the maximum value of F s /(kd) is decreasing with an increasing adimensional temperature K B T /(kd 2 ) and is increasing with an increasing ratio k 0 /k. This result is typically observed in several experimental measurements where the friction is an increasing function of the relative corrugation k 0 /k [44,52].

Thermally modified Frenkel-Kontorova model
We consider a one-dimensional chain composed of N particles linked by linear springs with elastic constant k m and equilibrium length s (see Fig. 3). This chain interacts with a substrate lattice composed by a sequence of pinning sites with stiffness k 0 and spaced by a uniform distance d (see Fig. 3). This arrangement corresponds to the Frenkel-Kontorova model, where a competition between the intrinsic chain spacing s and the substrate lattice spacing d is introduced [53][54][55][56][57][58][59]. At the same time, the elastic potential energy of the moving chain competes with the pinning potential energy of the substrate. At the extremity of the chain, we also introduce two more linear springs, with elastic constant k 1 and k N , which are mimicking the devices adopted to manipulate the chain at the nanoscale. They try to keep the two ends in the positions 1 and N , respectively, to the left and right of the chain (see Fig. 3). Actually, in order to study the stick-slip phenomenon at the nanoscale, it is sufficient to use only one device, let's say at the right of the chain, but we keep in the formalism the two lateral springs for symmetry reasons. In the original Frenkel-Kontorova model, the lattice substrate is represented by a periodic potential energy [53][54][55][56][57][58][59].
where W is the energy barrier between adjacent sites and x j are the particle positions. As before, we define here the substrate potential energy by introducing a sequence of quadratic wells (elastic constant k 0 ). Then, we introduce thermal fluctuations in the system and we work at thermodynamic equilibrium, corresponding to rate-independent or static friction. In order to identify the pinning site concerned for each particle, we add N discrete or spin variables n 1 , ..., n N to the phase space of the system. It means that the value of n j corresponds to the number of the site occupied by the j-th particle. To drawn a comparison with the original Frenkel-Kontorova model, we observe that the energy barrier in our modified approach is given by W = k 0 d 2 /8, being defined by the intersection of two adjacent quadratic wells. As previously discussed, W must be sufficiently larger than the thermal energy K B T . These assumptions lead to total energy where x = (x 1 , ..., x N ) ∈ R N is the particle positions vector and n = (n 1 , ..., n N ) ∈ Z N is the spin variables vector. We stress again that both those vectors belong to the phase space of the system. Then, the Helmholtz partition function can be calculated as where we integrated over x and summed over n. To simplify the calculations, it is not difficult to rewrite the total energy U ( x, n) in the following compact form Here, we introduced the tridiagonal matrix A defined below with all the subdiagonal and superdiagonal elements equal to -1 and the diagonal elements given by Moreover, we defined the vector v as and the scalar quantity (independent of x) These premises allow us to rewrite the partition function in the form that can be simplified using the gaussian integration holding for any symmetric and positive definite matrix N . This expression represents the multidimensional generalization of Eq. (3). The integration leads to the result where I ∈ M N ,N (R) is the identity matrix. The partition function can be finally written in terms of the so-called Riemann theta function [112,113] ( z| ) = which is defined when the matrix is in the following set: H N = ∈ M N ,N (C) : = T , m( ) > 0 . This is the set of symmetric square matrices whose imaginary part is positive definite. The set H N is called the Siegel upper half-space and is the multi-dimensional analog of the upper half-plane in C. Here, z ∈ C N is an arbitrary complex vector. It is known that the Riemann theta function converges absolutely and uniformly on compact subsets of C N × H N [112,113]. It can be seen as the multidimensional generalization of the function ϑ 3 (z, τ ) introduced in Eq. (6) to study the single particle friction. Through the identifications we obtain the partition function as Since it is possible to prove that the real matrix I − k 0 k m A −1 is positive definite, defined in Eq. (41) belongs to H N , and then ( z| ) in Eq. (43) is always convergent. The Riemann theta function is able therefore to fully describe the complexity of the rate-independent stick-slip phenomena modeled by our thermally modified Frenkel-Kontorova approach. In particular, it perfectly handles the geometric and energetic compromise between chain and substrate.
We can give a physical interpretation of the exponential term summed in Eq. (39), i.e., summed within the Riemann theta function. If we look at the total energy of the system given in Eq. (32), we can minimize it with respect of the particle positions x by stating that This condition delivers the following relation between the minimizing vector x 0 and n Now, we can substitute this value x 0 in U ( x, n), by obtaining where we have omitted all terms not depending on n. We can note that the n-dependent part of U ( x 0 , n) in Eq. (46) represents the exponent of the terms summed in Eq. (39) (or within the Riemann theta function), except for the Boltzmann factor −1/(K B T ). This means that everything goes as if the Boltzmann statistic were applied only to the discrete spin variables after minimizing the energy with respect to the continuous coordinate variables. This is possible only because all interactions are based on quadratic energies and therefore on a linear response with respect to the applied forces. The nonlinearity of the system emerges only thanks to the jumps of the spin variables, which make the overall potential energy non-convex.
Similarly to the one-particle case, also for the N -particle chain we can obtain a second mathematical form for Z H . To this aim, it is possible to introduce the following reciprocal formula or modular transformation (sometimes also called Jacobi identity) for the Riemann theta function [112][113][114] which is valid for any ∈ H N (see "Appendix" for details). Here, the square root assumes its principal value. This property can be applied to Eq. (43) and allows to obtain a second form for the partition function of our system, as follows where B = A − k 0 k m I. The possibility of having two different forms for the partition function is useful from a numerical point of view. In fact, if we look at the structure of the theta function, we observe that its convergence is faster if all the eigenvalues of 2π m( ) are much larger than one [113]. As a matter of fact, if this condition is verified, the Gaussian shape summed within the theta function is very narrow, and therefore, only a few values of n are needed to find the whole sum accurately. So, when 2π m( ), defined in Eq. (41), has large eigenvalues (larger than one) it is better to use Eq. (43), and when 2π m( ) has small eigenvalues (smaller than one) it is better to use Eq. (48). We have observed numerically that frequently (but not always) the matrix 2π m( ) is well-conditioned in the sense that the eigenvalues are typically either all larger than one or all between zero and one. We can also affirm that Eq. (43) is a partition function expression rapidly convergent in the low temperature regime, whereas Eq. (48) is more adapted in the high temperature regime (always with K B T W ). Anyway, numerical algorithms for the calculation of the Riemann theta function can be found in the literature [115].
The knowledge of the partition function in both forms given in Eqs. (43) and (48) allows the determination of the force applied by the device to hold the chain in a given position. In the most general case, we fix the left device at position 1 and the right device at position 2 . Consequently, the two applied forces can be determined as N ) being the Helmholtz free energy of the system. Two particular cases are of interest: (i) if the first extremity is left free, which means that k 1 = 0, then the partition function depends only on N and the right force reads where we have considered N = to simplify the notation; (ii) if the two extremities are rigidly connected, we can write 1 = − (N − 1)s and N = [which means that the distance N − 1 is maintained equal to the equilibrium value (N − 1)s] and the total force is representing the sum of the two left and right forces. Examples of application of Eqs. (51) and (52) are given below.
To determine the average value n of the spin variables vector, we can simply use the expression of the partition function given in Eq. (39), and we eventually get This expression can be written in terms of the theta function and its derivatives as follows. First of all, we can write the gradient of the Riemann theta function with respect of the vector variable z in the form Then, we can rewrite Eq. (53) as which is the multi-dimensional generalization of Eq. (14). Unfortunately, there are no explicit formulas for the logarithmic derivatives of the Riemann theta function, and therefore, we cannot obtain more explicit expressions for n or other quantities as we had done for a single particle system. By deriving Eq. (38) with respect to b, we obtain the integral expression which can be used to determine the average value x of the particles position vector. Indeed, by definition we have and by identifying N with k m A/(K B T ) and b with (k 0 d n −sk m v)/(K B T ), we easily obtain the simple relation Remarkably, it means that the knowledge of n is sufficient to determine the vector x . Moreover, it is also interesting to remark that Eq. (58) has exactly the same mathematical form of Eq. (45). It means that the minimizing of the total energy with respect to x delivers correct results even if the thermal fluctuations act on the system. In order to evaluate the effect of the substrate on the sliding chain, we can define this sort of adimensional order parameter which assumes the value zero if the chain is at equilibrium and larger values for an heterogeneous nonequilibrium spacing. Of course, depends on during the sliding of the chain, and therefore, we can also introduce its mean value over a period of the corrugated substrate which is useful to have an overall quantification of the trade-off between chain and substrate energies.

Results for the modified Frenkel-Kontorova model
The theoretical results previously shown have been numerically implemented in order to explore the behavior of the thermally modified Frenkel-Kontorova model. All examples have been studied under the assumption W K B T to ensure that we stay within the validity of the theory. We explore here the following issues: the Frenkel-Kontorova stick-slip with sliding boundary devices, the effect of thermal fluctuations, the effect of the elastic contrast between chain and substrate, the effect of the boundary devices intrinsic elasticity, and finally, the effect of the geometric contrast between chain and substrate.

Frenkel-Kontorova stick-slip with sliding boundary devices
A first example is found in Fig. 4, where we show the evolution of the static friction force, the particles position, and the order parameter as function of the position of the device for the two cases (i) and (ii) defined above. It means that we considered the case (i) with the device acting only on the right of the chain and the case (ii) where the device acts on both extremities, keeping their distance fixed at s(N − 1). In both cases, the device spans one period of the corrugated substrate and the previously mentioned physical quantities are recorded for two chains of length N = 4 and N = 5. In this first example, all parameters are fixed as declared in the caption of Fig. 4. We observe that in both cases (i) and (ii), to keep the device at /d = 0, we have to apply a given force since this configuration is not at equilibrium. Of course, this issue depends in a complex way on the ratios s/d and k m /k 0 . Indeed, we can observe that there are some values of /d for which the force is zero, and then, they correspond to the equilibrium configurations. Anyway, in the curves of F s /(k 0 d), we can see some jumps corresponding to the transitions of the particles between adjacent substrate wells. We adopted a low temperature and a low ratio k m /k 0 to amplify and then easily identify this effect. This behavior is consistent with the evolution of x /d, where we can observe the transitions of each single particle (stick-slip processes). It is interesting to see that a transition of one particle between two adjacent substrate wells may induce a step behavior (smoothed by the temperature) in the position of other particles not transiting between the wells. This is due to the strong energetic coupling of all variables within the Frenkel-Kontorova model. Furthermore, we can observe that at transitions there is a significant decrease in the order parameter . In fact, it represents a measure of the equilibrium state of the chain and, therefore, must decrease as a result of the transitions that occur due to the geometric and energetic compromise between chain and substrate. In other words, transitions occur because they are energetically favored and thus bring the system into a configuration closer to equilibrium, with a lower . From Fig. 4, we can also deduce that the friction force patterns are more  Fig. 4. We considered the chain lengths N = 4 (first column) and N = 5 (second column) and the two boundary conditions (i) (first row) and (ii) (second row). All parameters are indicated in the caption of Fig. 4 complicated in the case of the boundary condition (ii) since this configuration is more constraining for the system, thus generating a more complex sequence of transitions.
We are studying the geometric and energetic competition between a mass-spring chain and a corrugated substrate. In this context, Aubry proved an interesting transition originated by this competition in the Frenkel-Kontorova model with an infinite (N → ∞) and incommensurate (s/d irrational) chain without temperature effects [64][65][66]. For values of the ratio k 0 /k m below a certain threshold, chain elasticity prevails over the energy barrier of the substrate and we observe an unpinned sliding phase. In this case, the static friction vanishes and we say that the sliding is superlubric (an infinitesimal force can move the chain). We remark that even if the static friction is zero, the kinetic friction is nonzero since the positive speed induces the excitation of phonons in the chain, requiring a certain amount of energy. On the other hand, above the critical value of k 0 /k m the particles are pinned by the corrugated substrate since the substrate energy is predominant over the chain elastic energy. In this case, we have a nonzero static friction to move the chain since it must be able to overcome the Peierls-Nabarro barrier (and in this case the static friction force is analogous to the Peierls-Nabarro stress) [61][62][63]. The threshold or critical value of k 0 /k m that discriminates unpinned and pinned states depends in a very intricate way on the incommensurability ratio s/d [67]. In particular, it reaches its highest value when s/d is equal to the golden ratio ( √ 5 − 1)/2 [59]. The shift between the unpinned sliding phase and the pinned condition with stick-slip motion generates the so-called analyticity breaking of the particle positions. It means that the curves generated by plotting the actual particle position x j versus the unperturbed position x j,0 (without substrate) exhibit discontinuities or gaps, induced by the jump of the particle between substrate wells. In general, the relation between x /d and x 0 /d (with variable) is called hull function and shows discontinuities in the pinned Aubry state. Then, in this pinned state the hull function analyticity breaking, the crossing of the Peierls-Nabarro barrier and the stick-slip motion are three different aspects of the same phenomenology. Although the Aubry transition exists only for infinite and incommensurable chains at T = 0, it is interesting to observe the shape of the hull functions for our system as well. In our case, the reference positions are given by x j,0 = − (N − j)s ∀ j = 1, ..., N where is the variable position of the sliding device. Of course, these references positions can be used for both boundary conditions (i) and (ii), previously introduced. In Fig. 5, one can find the hull functions for the chains considered in Fig. 4. We considered the chain lengths N = 4 (first column) and N = 5 (second column) and the two boundary conditions (i) (first row) and (ii) (second row). We can easily recognize the jumps in the hull functions, but they are not real discontinuities because of the effect of the thermal fluctuations that are able to smooth out their shape. We observe that our finite-size finite-temperature Frenkel-Kontorova model is always statically pinned, even for an irrational value of s/d since we always have a positive static friction corresponding to a positive Peierls-Nabarro barrier. It is also interesting to remember that an Aubry-like transition, based a symmetry-breaking phenomenon, has been identified for finite chain at zero temperature [68][69][70]. However, in our finite chain at finite temperature all transitions are smooth as one can see in Fig. 5, and the response is controlled by the compromise between pinning intensity and thermal fluctuations that can facilitate the chain sliding. Of course, the slope of the jumps in our hull functions is larger for lower temperatures. Besides, from Fig. 5 we deduce that the jumps of the transitions are larger for the boundary condition (i) since the second condition (ii) is more constraining, limiting the positions discontinuity. However, under condition (ii), the height of the jumps is smaller but the pattern of transitions is more intricate, as shown in Fig. 4. Interestingly, Aubry-type transitions have been experimentally observed in different finite systems, as discussed in the recent literature for one-dimensional [71][72][73][74] and two-dimensional [75] structures.

Effect of thermal fluctuations
For the moment, we have considered all the fixed parameters of the system to show some initial examples. We can now study the effects induced by varying the temperature T . In Fig. 6, we plotted the adimensional static friction force F s /(k 0 d), the average particles positions x /d, and the order parameter as function of the sliding ratio /d for a chain with length N = 4, and for both boundary conditions (i) and (ii). The main important effect of the temperature is found in Fig. 6 where the static friction quantified by F s /(k 0 d) is clearly decreasing with increasing temperature. At the same time, the transitions in the curves of F s /(k 0 d) versus /d are smoothed with increasing temperature. The friction reduction with increasing temperature has been already pointed out for the system with one particle and can be also related to the decreasing of the Peierls-Nabarro stress with increasing temperature classically observed in materials science [116][117][118][119]. In the plot of x /d versus /d in Fig. 6, we observe that particle transitions between different substrate wells are smoother when the temperature is elevated. This means that the temperature fosters these transitions allowing a more complete exploration of the phase space, even where the total energy is higher. This behavior is directly reflected in the graph of , which decreases significantly as the temperature increases. In particular, the minimum points of the function are in correspondence of the transitions, occurring as soon as they are energetically favored and, therefore, bringing the system closer to equilibrium. As the chain is already closer to equilibrium as a result of these transitions, temperature has a lesser effect in reducing the value of at these points. This can be seen in the graphs of versus /d in Fig. 6, where at the transitions the curves at different temperatures are very closed together as opposed to the regions without transitions where the curves are more spaced.

Effect of the elastic contrast between chain and substrate
The effects induced by varying the chain elastic constant k m are shown in Fig. 7. First of all, from the plots of F s /(k 0 d) versus /d, we deduce that the static friction force rises as the elastic constant k m increases. This behavior is similar to the one obtained for a single particle in Fig. 2. However, in the present case with N = 4 particles, we also observe that the transitions occur at positions that depend significantly on the elastic constant k m [for both boundary conditions (i) and (ii)]. This is easily explained by the fact that k m controls the forces in the system and thus modifies the energy competition between chain and substrate. The transitions shift with k m can be clearly seen in both plots of x /d versus /d and versus /d. The latter was represented on a semi-logarithmic scale because the effects of varying k m on are very pronounced. It is in fact evident that a strong increase in k m brings the system very close to the equilibrium configuration with the particles equispaced at distance s.
The effect of varying the temperature T and the elastic constant k m can be summarized by means of the plots in Fig. 8. Here, one can find the influence of k m and T on the behavior of the maximum static friction  interval (0, d). This allowed us to determine the maximum friction force and the average value of the order parameter, as defined in Eq. (60). The final result is that the maximum of friction force F s /(k 0 d) increases with k m /k 0 and decreases with K B T /(k 0 d 2 ) (for each value of N and for both boundary conditions (i) and (i)). Similarly, the average order parameter decreases with both k m /k 0 and K B T /(k 0 d 2 ). This analysis is useful when k 0 is constant and k m variable. The dual situation with k 0 variable and k m constant is not represented in Fig. 8, but, as discussed for the system with one particle, its overall behavior can be summarized as follows. The maximum value of F s /(k m d) is decreasing with an increasing temperature K B T /(kd 2 ) and is increasing with an increasing ratio k 0 /k m . At the same time, decreases with K B T /(k m d 2 ) and increases with k 0 /k m .

Effect of the boundary devices intrinsic elasticity
To further show the complexity of the system under investigation, we show in Fig. 9 the effect of the boundary devices on the transitions pattern. We considered a chain with length N = 4 under boundary condition (i), i.e., with a single device pulling the chain from the right extremity. This device is here modified by changing its elastic constant k N . It is interesting to observe that for the first two values adopted k N /k 0 = 0.5, 2 we have no transitions for /d < 0.5 while for the other values k N /k 0 = 8, 32 we observe the emergence of one transition around /d 0.3. This transition origination can be easily identified in all panels of Fig. 9 (see dashed ellipses). We give a direct interpretation of this result as follows. The initial positions of the ideal chain particles for /d = 0 would be 4 and x 4 /d 0, without interaction with the corrugated substrate. Hence, at the beginning, the first two particles are within the well with n = −1 (−3/2 < x/d < −1/2) and the others in the well with n = 0 (−1/2 < x/d < 1/2). So, in an ideal situation, we should expect a first transition for /d 0.3, when x 2 /d −0.5, i.e., when the second particle reaches the transition threshold between the two neighboring wells. However, this transition is not observed for low values of k N since in this case the device, instead of generating chain sliding, generates its own elongation due to its weak elasticity. In this case, the chain pinning is stronger than the device pulling. Instead, for higher values of k N , the device deforms much less and the chain has to move, generating the transition of the second particle when /d 0.3. This example explains that there is not only a geometric and energetic trade-off between chain and substrate but also a competition between device elasticity and chain elasticity. This brings additional complexity to the friction phenomenon we are studying. Up to this point, we have always considered the s/d ratio to be constant and therefore assumed that we were working with the same geometric competition between chain and substrate. Now it is interesting to explore the effects of the ratio s/d on the stick-slip process and on the static friction. The results of this analysis are found in Figs. 10,11,12, and 13, where we considered chains of length N = 2, N = 3, N = 4, and N = 5, respectively. In each of these figures, we plotted the static friction surface F s /(k 0 d) as function of /d and  The results shown in Figs. 10, 11, 12, and 13 concern only one value of temperature. We performed the same analysis by changing the temperature over the interval K B T /(k 0 d 2 ) = 0.01 j, ∀ j = 1, .., 10, and we collected the outcomes in order to show the maximum static friction force max {| F s /(k 0 d)|, 0 < /d ≤ 1} as function of s/d and parametrized by the temperature. The result of this operation is found in Fig. 14, where we considered the chains of length N = 2, N = 3, N = 4 and N = 5 in the first, second, third and fourth row, respectively. Moreover, we used the chain to substrate elastic constants ratio k m /k 0 = 1 in first column, k m /k 0 = 10 in the second column, and k m /k 0 = 50 in the third column. The collection of these graphs shows well the effect of geometric compromise between chain and substrate as the temperature increases. In fact, the first thing we can observe is that the maximum static friction is always decreasing with temperature regardless of the structural parameter s/d. As far as the geometric ratio s/d is concerned, we note that there are particular values of this parameter that significantly reduce the maximum friction force when k m is sufficiently large. For We can explain this behavior as follows. We suppose to work with k m and k N sufficiently large, and therefore, we consider here a sliding finite chain at equilibrium, corresponding to the positions x j = − s(N − j) ∀ j = 1, ..., N . It means that we introduce the geometrical competition between chain and substrate, but we neglect the corresponding energetic trade-off. Similarly, we also neglect the thermal fluctuations. We also consider in this calculation the original Frenkel-Kontorova potential stated in Eq. (29). So, we can write where we defined η = 2π /d and ξ = 2πs/d. Now, we can easily develop this expression as follows where we used the sum N j=1 . We can now multiply numerator and denominator of the internal fraction by 1 − e −iξ and we get Here, we can use a prosthaphaeresis formula obtaining the simplified result To conclude, we can use another prosthaphaeresis formula to get the final result This expression represents the Frenkel-Kontorova potential energy of a non-deformable chain interacting with the corrugated substrate.
The force necessary to keep the chain at /d = η/(2/π) is given by F s = ∂U F K /∂ = (2π/d)∂U F K /∂η, or more explicitly by The maximum value of this force over the period of the corrugated substrate is then simply obtained as max η∈(0,2π) This result perfectly explains the behavior of the friction force observed in Fig. 14. Indeed, the numerator of Eq. (67) vanishes when N ξ/2 = jπ, ∀ j = 0, ..., N , which means s/d = ξ/(2π) = j/N , ∀ j = 0, ..., N . However, the two extremes values s/d = 0 and s/d = 1 must not be considered because they are zeros of the denominator and they in fact correspond to the maximum value of Eq. (67). This scenario exactly corresponds to the properties (I) and (II) stated above, as one can also see in Fig. 15. Finally, the particular values of s/d that significantly reduce the maximum friction force are explained by the pure geometrical competition between chain and substrate. It is interesting to remark that for an infinite chain the Aubry sliding phase is observed with incommensurate structures while for a finite chain the maximum reduction in friction can be obtained for commensurate structures. This is not a contradiction since the incommensurate Aubry chain is elastically deformable, whereas we considered here a rigid chain to explain the friction reduction. Moreover, in our analysis of Fig. 14 we have considered the complete behavior of the system including elasticity and thermal fluctuations.

Conclusions
In this work, we have considered the equilibrium statistical mechanics of stick-slip processes, with application to the nanoscale friction. In other words, we can say that we studied the effect of thermal fluctuations on the rate-independent nanoscale interactions between a slider and a corrugated substrate. We have conveniently modified the one-particle Prandtl-Tomlinson model and the N -particle Frenkel-Kontorova model so that we can easily calculate the partition function of both systems. In particular, instead of considering the classical sinusoidal energy profile of the corrugated substrate, we introduced a sequence of quadratic potentials, each of which represents a substrate well. It is important to remark that the sinusoidal and the quadratic profiles have a qualitatively similar but quantitatively different behavior. Indeed, both profiles are periodic, but their shapes are quite different. Therefore, the scheme with the sequence of parabolic wells should be considered as an independent model with specific applications. First of all, it is interesting because of the possibility of a comprehensive analytical development, and moreover, it is well suited to the study of the propagation of dislocations and fractures in solids or the sliding of biological structures [12,28]. Of course, the quadratic potentials are integrable and then allow a rather simple analytic treatment of the problem. Following the classical spin variables approach, we had to add to the phase space a discrete variable for each particle useful to identify the well explored by the particle itself. In this context, it is important to remark that the calculation of the partition function based on the spin approach assumes that for all substrate well, all possible particle positions are considered (also outside the well extension). This effectively corresponds to a multivalued energy function (see the superpositions of the quadratic curves in Figs. 1 and 3). As shown in Refs. [91,96], the effect of this approximation is statistically negligible since the potential energy in the regions outside the actual well extension is sensibly higher than the energy in the internal region. Of course, this is true if the barrier W between two adjacent wells is sufficiently larger than K B T . On the one hand, the application of this technique to the (modified) one-particle Prandtl-Tomlinson model leads to a partition function which can be written in terms of the third Jacobi theta function. Then, the application of the classical Jacobi functional identity allows us to obtain the closed form Fourier series describing the static friction force experienced by the single particle sliding on the corrugated periodic substrate through a stick-slip motion. Importantly, this force depends on the temperature, and therefore, the obtained results explain clearly the thermolubricity phenomenon. On the other hand the spin variable method applied to the (modified) N -particle Frenkel-Kontorova model leads to a partition function written in terms of the N -dimensional Riemann theta function. It is interesting to underline that the Riemann theta function is able to fully describe the complexity of the rate-independent stick-slip phenomenon controlled by the geometric and energetic compromise between chain and substrate. In fact, in this case, we have the geometric compromise induced by the difference of the characteristic lattice lengths of chain and substrate and the energetic compromise induced by the difference between the elastic constants of chain and substrate. We discussed the effect of the elastic constants and of the temperature of the stick-slip and friction responses. Moreover, we plotted the hull functions for the particles of the chain in order to compare the transitions between the wells in our system and in the infinite incommensurate Aubry configuration. We have shown that in our case with finite chains and nonzero temperature the transitions are smooth, in contrast with those observed by Aubry for infinite incommensurable chains at zero temperature, where one finds the so-called analyticity breaking of the particle positions. To conclude, we explored the effects of the geometric trade-off induced by the ratio s/d. In particular, we proved the existence for finite chains of specific commensurate values of the ratio s/d, which are able to strongly reduce the friction experienced by the slider. We investigated the effect of the temperature on this friction reduction. To conclude, we can say that our proposed method to take into account the temperature effect in friction phenomena can be adopted for several problems, including the application to the motion of defects in solid materials or the sliding of biological structures [12,28]. The originality of this approach will allow the study of a number of problems that are for the moment only partially or not yet theoretically modeled. These perspectives concern the generalization of this approach to the two-dimensional case [120], to more complex corrugated substrates having for example two or more energy levels for the intercalated wells [121], to corrugated soft substrates with elastically moving wells [122], to the interplay between adhesion and friction with application to bio-systems [123], to the kinetic or dynamic case by means of the Langevin and/or Fokker-Planck methodologies [124,125], and so forth. In particular, specific techniques for studying the dynamics of systems with multiwell energy landscapes will be exploited [106,126,127].
Funding No funding was received to assist with the preparation of this manuscript.

Declarations
Financial interest The author certifies that he has 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.
Published article All data generated or analyzed during this study are included in this published article.

Appendix A: reciprocal relation for Jacobi and Riemann theta functions
In Appendix, we briefly discuss the conceptual relationship between direct and reciprocal Bravais lattices, Dirac comb, Poisson summation formula, and reciprocal relation for the Riemann theta function and third Jacobi theta function. We consider a Bravais lattice composed of points r = n 1 a 1 + ... + n N a N ∈ R N where n j ∈ Z and a j ∈ R N are the primitive vectors. We introduce a periodic sampling function f ( x) on this Bravais lattice constructed by delta functions (Dirac comb) n 1 a 1 ... − n N a N ) . (A1) To simplify the notation, we define a matrix C ∈ M N ,N (R) where the columns represent the vectors a j . It means that C = [ a 1 |...| a N ] and we can write with n = (n 1 , ..., n N ) T . We assume that C is not singular. Then, if we define η such that x = C η, we have Now, the function f ( η) is multi-periodic with period one along all directions η 1 , ..., η N . Hence, it can be developed in Fourier series where the coefficients are given by as we can easily prove by applying the substitution y = C η. Coming back to the variable x, we finally obtain which is the Fourier series of the Dirac comb defined in Eq. (A1) or (A2). We can also determine the Fourier transform of the same function Since δ( q) = R N e i q· x d x/(2π) N , we easily get It means that if r = C n ∀ n ∈ Z N is the direct Bravais lattice, then k = 2πC −T m ∀ m ∈ Z N is the reciprocal Bravais lattice. In other words, the Dirac comb on the direct lattice is Fourier transformed into a Dirac comb on the reciprocal lattice. This property can be applied to find the so-called N -dimensional Poisson summation formula. We consider an arbitrary function φ( x) for which the Fourier transform ( k) exists. We define the replication or periodic summation where f is defined in Eq. (A2) and * means convolution. Since the Fourier transform of the convolution is the product of the two Fourier transforms (convolution theorem), we have Now, we can remember that the Fourier transform of e i x· v is given by (2π) N δ( k − v), and then from Eq. (A11) we come back to the original function g( x), eventually obtaining which is the Fourier series of the periodic summation. For x = 0, this result delivers the Poisson summation formula We take now into account the following particular function φ( x) with its Fourier transform ( k) φ( x) = e 2πi 1 2 x·T x+ x· σ , (A14) where T ∈ M N ,N (C) with T = T T , m(T ) > 0 and σ ∈ C N . To conclude, we substitute Eqs. (A14) and (A15) in the Poisson summation formula given in Eq. (A13) and we get, after straightforward calculations, the N -dimensional Jacobi reciprocal relation for the Riemann theta function where we have identified = C T T C and z = C T σ . This proves Eq. (47) of the main text. Other more refined properties of ( z| ) can be found in the literature [112,113]. The result given in Eq. (A16) can be specialized to the case with N = 1 by obtaining the original Jacobi identity where we have conveniently compared the definitions of ϑ 3 (z, τ ) and ( z| ) with N = 1 [112,113]. This finally proves Eq. (17) of the main text as well.