<figure><img src="image_1.png" /></figure>

HAL
open science

# One step replica symmetry breaking and overlaps between two temperatures

Bernard Derrida, Peter Mottishaw

► To cite this version:

Bernard Derrida, Peter Mottishaw. One step replica symmetry breaking and overlaps between two temperatures. Journal of Physics A: Mathematical and Theoretical, 2021, 54 (4), pp.045002. 10.1088/1751-8121/abd4ad. hal-03146644

HAL Id: hal-03146644

https://hal.sorbonne-universite.fr/hal-03146644

Submitted on 19 Feb 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
PAPER • OPEN ACCESS

# One step replica symmetry breaking and overlaps between two temperatures

To cite this article: Bernard Derrida and Peter Mottishaw 2021 *J. Phys. A: Math. Theor.* **54** 045002

View the [article online](#) for updates and enhancements.

<figure><img src="image_2.png" /></figure>

IOP ebooks™

Bringing together innovative digital publishing with
leading authors from the global scientific community.

Start exploring the collection-download the
first chapter of every title for free.
# One step replica symmetry breaking and overlaps between two temperatures

Bernard Derrida<fnref n="1" />,<fnref n="2" /> ID and Peter Mottishaw<fnref n="3" />,*, ID

E-mail: bernard.derrida@college-de-france.fr and peter.mottishaw@ed.ac.uk

Received 8 September 2020, revised 10 December 2020
Accepted for publication 17 December 2020
Published 7 January 2021

<figure><img src="image_3.png" /></figure>

CrossMark

## Abstract

We obtain an exact analytic expression for the average distribution, in the thermodynamic limit, of overlaps between two copies of the same random energy model (REM) at different temperatures. We quantify the non-self averaging effects and provide an exact approach to the computation of the fluctuations in the distribution of overlaps in the thermodynamic limit. We show that the overlap probabilities satisfy recurrence relations that generalise Ghirlanda–Guerra identities to two temperatures. We also analyse the two temperature REM using the replica method. The replica expressions for the overlap probabilities satisfy the same recurrence relations as the exact form. We show how a generalisation of Parisi’s replica symmetry breaking ansatz is consistent with our replica expressions. A crucial aspect to this generalisation is that we must allow for fluctuations in the replica block sizes even in the thermodynamic limit. This contrasts with the single temperature case where the extremal condition leads to a fixed block size in the thermodynamic limit. Finally, we analyse the fluctuations of the block sizes in our generalised Parisi ansatz and show that in general they may have a negative variance.

Keywords: spin glass, replica, temperature chaos, random energy model

(Some figures may appear in colour only in the online journal)

<figure><img src="image_4.png" /></figure>

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

<footnote n="1">¹ Collège de France, 11 Place Marcelin Berthelot, 75005 Paris, France</footnote>

<footnote n="2">² Laboratoire de Physique de l'Ecole Normale Supérieure, ENS, Université PSL CNRS, Sorbonne Université, Université de Paris, F-75005 Paris, France</footnote>

<footnote n="3">³ SUPA, School of Physics and Astronomy, University of Edinburgh, Peter Guthrie Tait Road, Edinburgh EH9 3FD, United Kingdom</footnote>

<footnote>* Author to whom any correspondence should be addressed.</footnote>
# 1. Introduction

Since replica symmetry breaking (RSB) was invented by Parisi, 40 years ago [1], it has been used in many different contexts and the subtle physical meaning of the scheme he used has been elucidated [2–4] (for reviews see [5] or [6]). Here we would like to provide a simple example to explore how the Parisi scheme could be extended to calculate correlations between different temperatures.

In the replica approach, a central role is played by the overlaps which represent the correlations between pure states. For a system of N Ising spins with the interactions sampled from some disorder distribution (as in the Sherrington–Kirkpatrick model [7], for example), the overlap between a configuration C and a configuration C' is defined by

$$q(C, C') = \frac{1}{N} \sum_{i=1}^{N} \sigma_{i}^{C} \sigma_{i}^{C'}, \qquad (1)$$

where $\sigma_i^C = \pm 1$ is the value of the spin at site i in configuration C. The distribution $P(q)$ of this overlap at a single inverse temperature $\beta$ for a particular sample is then given by [2]

$$P(q) = \sum_{C, C'} \frac{e^{-\beta E(C)}}{Z(\beta)} \frac{e^{-\beta E(C')}}{Z(\beta')} \delta(q - q(C, C')), \qquad (2)$$

where $Z(\beta) = \sum_C e^{-\beta E(C)}$ is the partition function at inverse temperature $\beta$ and $E(C)$ is the energy of configuration C for the particular sample. In a disordered system the energies are quenched random variables and $P(q)$ is itself a random quantity, sample dependent in the sense that it depends on the energies $E(C)$. One of the achievements of Parisi’s theory of spin glasses was to predict that $P(q)$ remains sample dependent even in the thermodynamic limit, and to allow the calculation of various averages and moments which characterize its sample to sample fluctuations [3, 4, 8, 9].

The notion of overlap distribution can be generalized when the two configurations are at different temperatures

$$P_{\beta, \beta'}(q) = \sum_{C, C'} \frac{e^{-\beta E(C)}}{Z(\beta)} \frac{e^{-\beta' E(C')}}{Z(\beta')} \delta(q - q(C, C')). \qquad (3)$$

This clearly reduces to (2) when $\beta = \beta'$. These multiple temperature overlaps have mostly been studied in the context of temperature chaos, in order to see how the random free energy landscapes are correlated at different temperatures (see for example [10]). Several spin glass models exhibit temperature chaos, meaning that the overlap between different temperatures vanishes in the thermodynamic limit (in which case, the question of the fluctuations of $P_{\beta, \beta'}(q)$ becomes superfluous). One way to predict temperature chaos is to show that $P_{\beta, \beta'}(q)$ vanishes exponentially with the system size when $\beta \neq \beta'$ and $q > 0$ [11]. There are however models for which these multiple temperature overlaps do not vanish and the question of how the Parisi theory has to be modified in these cases is, to our knowledge, not fully understood (see for example the multi-p-spin models analysed in [12]). Here we attack this question in the simplest model which exhibits RSB, the random energy model (REM, see [13, 14]) which has the advantage of being open to both exact and replica analysis. This will allow us to propose a way to adapt Parisi’s scheme for the two temperature case, in order to be compatible with our exact
results of section 2. The absence of chaos in the REM has been discussed in [15] as well as its dynamical effects such as rejuvenation in [16].

One of our motivations is the hope that the insight developed here into the use of the replica method in two temperature situations will be applicable to spin glass problems where there is no alternative to the replica method such as the multi-p-spin spherical models discussed in [12].

In the REM, the energies $E(C)$ are $2^N$ independent random variables distributed according to a Gaussian distribution of width proportional to $N$. The overlap can then only take two values

$$q(C, C') = \delta_{C,C'}.$$

Therefore $P(q)$ consists of two delta function peaks [5, 17],

$$P(q) = (1 - Y_2) \delta(q) + Y_2 \delta(q - 1), \quad (4)$$

where $Y_2$ is the probability, at equilibrium, of finding two copies of the same sample in the same configuration.

$$Y_2 = \sum_C \left( \frac{e^{-\beta E(C)}}{\sum_C e^{-\beta E(C)}} \right)^2 = \frac{Z(2\beta)}{Z(\beta)^2}. \quad (5)$$

In the large $N$ limit, $Y_2$ vanishes in the high temperature phase ($\beta < \beta_c$), while in the low temperature phase ($\beta > \beta_c$) it takes non zero values with sample to sample fluctuations.

A direct calculation [8, 9] as well as a replica calculation [4] lead to

$$\langle Y_2 \rangle = 1 - \mu; \quad \langle Y_2^2 \rangle - \langle Y_2 \rangle^2 = \frac{\mu - \mu^2}{3} = \frac{\langle Y_2 \rangle - \langle Y_2 \rangle^2}{3}, \quad (6)$$

where $\langle \cdot \rangle$ denotes the disorder average i.e. the average over the random energies $E(C)$ and

$$\mu = \frac{\beta_c}{\beta}. \quad (7)$$

The quantity $Y_2$ can be generalized to the probabilities $Y_k$ of finding $k$ copies of the same sample in the same configuration

$$Y_k = \sum_C \left( \frac{e^{-\beta E(C)}}{\sum_C e^{-\beta E(C)}} \right)^k = \frac{Z(k\beta)}{Z(\beta)^k}. \quad (8)$$

As for $Y_2$, the large $N$ limits of the disorder averages of these overlaps are known [3, 4, 18, 19]

$$\langle Y_k \rangle = \frac{\Gamma(k - \mu)}{\Gamma(1 - \mu) \Gamma(k)}. \quad (9)$$

Since $\mu = 1 - \langle Y_2 \rangle$, equation (9) implies

$$(k\beta - \beta_c) \langle Y_k \rangle = k\beta \langle Y_{k+1} \rangle, \quad (10)$$

which can be seen as simple cases of the Ghirlanda–Guerra identities [20–23].

At two different temperatures, the overlap distribution (3) for the REM is still a sum of two delta functions

$$P_{\beta, \beta'}(q) = (1 - Y_{1,1}) \delta(q) + Y_{1,1} \delta(q - 1), \quad (11)$$
where the random variable

$$Y_{k,k'} = \sum_{c} \left( \frac{e^{-\beta E(C)}}{Z(\beta)} \right)^k \left( \frac{e^{-\beta' E(C)}}{Z(\beta')} \right)^{k'} = \frac{Z(k\beta + k'\beta')}{Z(\beta)^k Z(\beta')^{k'}} \quad (12)$$

is the probability that $k$ copies of a given sample at temperature $\beta$ and $k'$ at temperature $\beta'$ are all in the same configuration.

In section 2 and in the appendix A we will derive the following exact expressions of the sample averages of these generalized overlaps

$$\langle Y_{k,k'} \rangle = \frac{\beta}{\beta_c} \frac{1}{\Gamma(k)} \frac{1}{\Gamma(k')} \int_0^\infty dv \, v^{k'-1} \frac{\Psi\left(k + k' \frac{\beta'}{\beta}; v\right)}{(-\psi(v))}, \quad (13)$$

where the functions $\psi(v)$ and $\Psi(v)$ are given by

$$\psi(v) = \int_0^\infty du \, (e^{-u-vu^\frac{\beta'}{\beta}} - 1) u^{-1-\frac{\beta}{\beta}} \quad (14)$$

and

$$\Psi(z; v) = \int_0^\infty du \, e^{-u-vu^\frac{\beta'}{\beta}} u^{z-1-\frac{\beta}{\beta}}. \quad (15)$$

These two functions are generalisations of the Gamma function. Like the Gamma function which satisfies $\Gamma(z+1) = z\Gamma(z)$, they obey some recursion relations (which can be obtained from (15) via integrations by parts) leading to the following relations

$$(k\beta + k'\beta' - \beta_c) \langle Y_{k,k'} \rangle = k\beta \langle Y_{k+1,k'} \rangle + k'\beta' \langle Y_{k,k'+1} \rangle \quad (16)$$

which generalize (10).

We will also show that

$$\langle (Y_{1,1})^2 \rangle = \frac{\beta}{\beta_c} \int_0^\infty dv \, v \left[ \frac{\Psi\left(2 + 2\frac{\beta'}{\beta}; v\right)}{(-\psi(v))} + \frac{\Psi\left(1 + \frac{\beta'}{\beta}; v\right)^2}{(-\psi(v))^2} \right]. \quad (17)$$

By varying $\beta$ one can draw, using the exact expressions (13) and (17), the variance $\langle (Y_{1,1})^2 \rangle - \langle Y_{1,1} \rangle^2$ versus the average $\langle Y_{1,1} \rangle$ as in figure 1. Clearly the relation (6) (which is a direct consequence of Parisi's ansatz) is no longer satisfied when the two temperatures are different ($\beta \neq \beta'$).

In the rest of the paper we will show in section 2 how the expressions (13)–(15) and (17) can be derived directly. We will also give the generalisation of these expressions when the overlaps are weighted by the partition functions to some power (in the replica language when the number of replicas is non-zero). Then in section 3 we will show what needs to be done in order to calculate these overlaps in a replica approach and in section 4 we will propose a scheme which generalises Parisi's ansatz and is compatible with our exact results. Finally, in section 5 we will explore the nature of the fluctuations in block size that we observe in this generalisation of Parisi's ansatz to two temperatures.
<figure><img src="image_5.png" /><figcaption>Figure 1. The variance of $Y_{1,1}$ versus its average when $\beta = \beta'$ and $\beta = 2\beta'$. The curves are obtained by varying $\beta'$ between $\beta_c$ and $\infty$. The lines represent the expression (6) in the case $\beta = \beta'$ and (13) when $\beta = 2\beta'$. The points are the results of Monte Carlo simulations in these two cases.</figcaption></figure>

## 2. The direct calculation of the overlaps

In this section, after recalling the definition of the REM, we explain how the expressions of the overlaps such as (13)–(15) can be derived. In the REM, a sample is determined by the choice of $2^N$ random energies $E(C)$ chosen independently from a Gaussian distribution

$$P(E) = \frac{1}{\sqrt{N\pi J^2}} \exp\left[-\frac{E^2}{NJ^2}\right]. \qquad (18)$$

It is known that, in the thermodynamic limit ($N \to \infty$), there is a phase transition [13, 14] at an inverse temperature

$$\beta_c = \frac{2\sqrt{\log 2}}{J} \qquad (19)$$

and that in the frozen phase $\beta > \beta_c$ (which is the only phase with non-zero overlaps) the partition function is dominated by energies close to the ground state which itself has fluctuations of order 1 around a characteristic energy $E_0 = -JN\sqrt{\log 2} + J \log N/(4\sqrt{\log 2})$. As only energies at a distance of order 1 (i.e. $\ll N$) contribute to the partition function one can replace the REM by a Poisson REM [24] which has, in the frozen phase and for large system sizes, the same properties as the original REM [13]. In this Poisson REM (PREM), the values of the energies, for a given sample, are the points generated by a Poisson process [18] on the real line with intensity

$$\rho(E) = C \exp[\beta_c(E - E_0)] \quad \text{with} \quad C = \frac{1}{J\sqrt{\pi}}. \qquad (20)$$

One way to think of it is to slice the real axis into infinitesimal energy intervals ($E, E + dE$) indexed by $\nu$, and to say that there is an energy $E_\nu$ in the interval $\nu$ with probability
$p_ν = ρ(E)dE$. In other words the partition function is given by

$$Z(\beta) = \sum_{\nu=-\infty}^{\infty} y_{\nu} \exp[-\beta E_{\nu}], \qquad (21)$$

where the $y_ν$ are independent binary random variables such that $y_ν = 1$ with probability $p_ν$ and $y_ν = 0$ with probability $1 - p_ν$ (because the intervals are infinitesimal, there is no interval $ν$ occupied by more than one energy).

The details of the calculation leading to the expressions (13)–(15) are given in appendix A.

One can generalize (13) to obtain (see (76) of appendix A) the average over disorder of $Y_{k_1,k_2;k'_1,k'_2}$ defined by

$$Y_{k_1,k_2;k'_1,k'_2} = \frac{\sum_{\nu \neq \nu'} y_{\nu} y_{\nu'} e^{-(\beta k_1 + \beta' k'_1)E_{\nu} - (\beta k_2 + \beta' k'_2)E_{\nu'}}}{Z(\beta)^{k_1+k_2} Z(\beta')^{k'_1+k'_2}}.$$

In particular this allows one to obtain (17) from (76) as one has from the definition (12)

$$(Y_{1,1})^2 = Y_{2,2} + Y_{1,1;1,1}.$$

In the replica approach, as we will see, it is often convenient to deal with weighted overlaps defined as

$$\langle Y_{k,k'} \rangle_{n,n'} = \frac{\langle Z(k\beta + k'\beta') Z(\beta)^{n-k} Z(\beta')^{n'-k'} \rangle}{\langle Z(\beta)^n Z(\beta')^{n'} \rangle} \qquad (22)$$

(see (12)). These averages can be performed (see appendix A) to get expressions (82) and (83) (valid for $n < 0$ and $n' < 0$) which generalize (13) and (16)

$$\langle Y_{k,k'} \rangle_{n,n'} = \frac{(-r) \Gamma(-n) \Gamma(-n')}{\Gamma(k-n) \Gamma(k'-n')} \times \frac{\int_0^\infty dv \, v^{k'-n'-1} \Psi(k+k'\frac{\beta'}{\beta}; v) (-\psi(v))^{r-1}}{\int_0^\infty dv \, v^{-n'-1} (-\psi(v))^r}, \qquad (23)$$

where

$$r = \frac{n\beta + n'\beta'}{\beta_c} \qquad (24)$$

and

$$(k\beta + k'\beta' - \beta_c) \langle Y_{k,k'} \rangle_{n,n'} = (k-n)\beta \langle Y_{k+1,k'} \rangle_{n,n'} + (k'-n')\beta' \langle Y_{k,k'+1} \rangle_{n,n'}. \qquad (25)$$

**Remark 2.1.** The $n, n' \to 0$ limit

As explained in appendix A (see (85)) the expression (23) reduces to (13) in the limit $n \to 0^-$ and $n' \to 0^-$. Very much like the integral representation of the Gamma function the expression (23) would take a different form for $n$ and/or $n' > 0$, and so we would need to use these alternative expressions to verify that the limits $n \to 0^+$ and $n' \to 0^+$ lead also to (13).

**Remark 2.2.** The $\beta' = \beta$ case
In the particular case where $\beta = \beta'$, one can perform the integrals in (14), (15) as in (79), (80) and obtain more explicit expressions of the overlaps as in (81). In particular one gets (see (81)) that

$$ \langle Y_{k,k'} \rangle = \frac{1}{\Gamma(k + k')} \frac{\Gamma(k + k' - \frac{\beta_c}{\beta})}{\Gamma(1 - \frac{\beta_c}{\beta})} \qquad (26) $$

which agrees with (9) as $Y_{k,k'} = Y_{k+k'}$ when $\beta = \beta'$ (see the definitions (9) and (12)). Similarly (23) becomes when $\beta = \beta'$

$$ \langle Y_{k,k'} \rangle_{n,n'} = \frac{\Gamma(1 - n - n') \Gamma(k + k' - \frac{\beta_c}{\beta})}{\Gamma(k + k' - n - n') \Gamma(1 - \frac{\beta_c}{\beta})} \qquad (27) $$

As in this single temperature case $Y_{k,k'} = Y_{k+k'}$, (16) reduces to equation (10). Therefore once $\langle Y_2 \rangle$ is known, all the other $\langle Y_k \rangle$ can be determined by the relations (10). Clearly for $\beta' \neq \beta$ this is not the case. However, (16) or (25) would allow to determine all the $\langle Y_{k,k'} \rangle$ from the knowledge of all the $\langle Y_{k,1} \rangle$.

In the rest of the paper we will see how expressions (13) or (23) can be interpreted in terms of replica symmetry breaking.

### 3. The replica method

In this section we apply the replica method to the REM. To illustrate the approach we first recall the computation of the free energy and overlap probability $\langle Y_k \rangle_n$ for a single temperature. It is well known that a single step in Parisi's replica symmetry breaking scheme [1] gives the correct low temperature solution [14, 17]. Here we will start with a slightly more general approach that allows for fluctuations in the block sizes. In the single temperature case the need to allow fluctuations in the block sizes has been discussed in [25, 26] and used in [24, 27] to compute finite size corrections in the REM.

In the two temperature case block size fluctuations have also been discussed in the context of temperature chaos in spin glasses (see [28, 29] appendix G for a detailed discussion), but computing the full overlap distribution between two temperatures by averaging over these block fluctuations has proved challenging. In this section we outline a replica symmetry breaking scheme for the two temperature case and show that it satisfies the same recursion (16) as the exact solution.

#### 3.1. The REM at a single temperature

To implement the replica method, the first step is to calculate the integer moments of the partition function, $\langle Z(\beta)^n \rangle$, then one assumes that the expression is valid for non-integer $n$ and finally one makes use of

$$ \langle \log Z \rangle = \lim_{n \to 0} \frac{\log \langle Z(\beta)^n \rangle}{n} \qquad (28) $$

to obtain the disorder average of the free energy.
### 3.1.1. Integer moments of the partition function.
As shown in appendix B the integer moments of the partition function for the REM are given by

$$ \langle Z(\beta)^n \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\mu_1 \ge 1} \cdots \sum_{\mu_r \ge 1} C_{n,r}(\{\mu_i\}) \langle Z(\beta\mu_1) \rangle \langle Z(\beta\mu_2) \rangle \cdots \langle Z(\beta\mu_r) \rangle, \quad (29) $$

(see (91) in appendix B) where

$$ C_{n,r}(\{\mu_i\}) = \frac{n!}{\mu_1!\mu_2!\cdots\mu_r!} \delta\left[\sum_{i=1}^{r} \mu_i = n\right] \quad (30) $$

and the Kronecker delta $\delta[\sum_{i=1}^{r} \mu_i = n]$ ensures that the $\{\mu_i\}$ sum to $n$. As

$$ \langle Z(\beta) \rangle = \int_{-\infty}^{\infty} \rho(E)e^{-\beta E} dE = e^{N f(\beta)} \quad (31) $$

with

$$ f(\beta) = \log 2 + \frac{(\beta J)^2}{4}. \quad (32) $$

One can rewrite (29) as

$$ \langle Z(\beta)^n \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\mu_1 \ge 1} \cdots \sum_{\mu_r \ge 1} C_{n,r}(\{\mu_i\}) \ e^{NA(r, \{\mu_i\})}, \quad (33) $$

where

$$ A(r, \{\mu_i\}) = \sum_{i=1}^{r} f(\mu_i, \beta) = r \log 2 + \frac{(\beta J)^2}{4} \sum_{i=1}^{r} \mu_i^2. \quad (34) $$

One can interpret a given term of the sum (29) or (33) as $n$ replicas distributed over $r$ distinct configurations with $\mu_i$ replicas in configuration $i$. Before using the replica method to compute the free energy from these expressions we compute the overlap probability for integer $n$.

An expression for the overlap probability $Y_k$, defined in (8), can be obtained from the ratio of integer moments

$$ \langle Y_k \rangle_n = \frac{\langle Y_k Z(\beta)^n \rangle}{\langle Z(\beta)^n \rangle} = \frac{\langle Z(k\beta) Z(\beta)^{n-k} \rangle}{\langle Z(\beta)^n \rangle}. \quad (35) $$

The denominator in the rightmost expression is the single temperature moment in (33). If $n$ and $k$ are positive integers with $n > k$ then the numerator is a two temperature moment of the partition function that is computed (see (96) and (97) in appendix B) and one gets

$$ \langle Y_k \rangle_n = \left\langle r \frac{\mu_1(\mu_1 - 1) \cdots (\mu_1 - k + 1)}{n(n-1) \cdots (n-k+1)} \right\rangle_{\{\mu_i\}}, \quad (36) $$

where the average $\langle \cdot \rangle_{\{\mu_i\}}$ means that for any function $F(\{\mu_i\})$,

$$ \langle F(\{\mu_i\}) \rangle_{\{\mu_i\}} = \frac{\sum_{r \ge 1} \sum_{\{\mu_i\} \ge 1} F(\{\mu_i\}) W_r(\{\mu_i\})}{\sum_{r \ge 1} \sum_{\{\mu_i\} \ge 1} W_r(\{\mu_i\})}. \quad (37) $$
with (see (30) and (33))

$$W_r(\{\mu_i\}) = \frac{C_{n,r}(\{\mu_i\})}{r!} e^{NA(r, \{\mu_i\})}.$$

### 3.1.2. The thermodynamic limit and the extremal condition.
In the thermodynamic limit $\langle Z(\beta)^n \rangle$ in equation (33) should be dominated by terms which maximize $A(r, \{\mu_i\})$ in equation (34). At high temperatures the maximum corresponds to all $n$ replicas being in different configurations. Thus $r = n$, $\mu_i = 1$ for all $i$. Then

$$\langle Z(\beta)^n \rangle \approx e^{N n f(\beta)} \qquad (38)$$

which gives (see (28))

$$\langle \log Z \rangle = N f(\beta) = N \left[ \log 2 + \frac{(\beta J)^2}{4} \right]. \qquad (39)$$

The entropy of this solution is

$$\langle S \rangle = N [f(\beta) - \beta f'(\beta)] = N \left[ \log 2 - \frac{(\beta J)^2}{4} \right]. \qquad (40)$$

There is a critical inverse temperature $\beta_c$ where this entropy vanishes. It is the solution of

$$f(\beta_c) - \beta_c f'(\beta_c) = 0 \qquad (41)$$

and therefore given by (19). When $\beta < \beta_c$ the entropy in (40) is positive and (39) is indeed the right free energy [13]. On the other hand at low temperatures, when $\beta > \beta_c$, the entropy is negative and one must look for a different solution.

To do so we proceed as Parisi did in his original papers [1, 2, 30, 31] on replica symmetry breaking. To identify the terms that dominate the sum in (33) in the thermodynamic limit ($N \to \infty$) in the low temperature phase we make the following three assumptions:

(a) We expect all the dominant terms to have a large $N$ behaviour of the form $\exp NA(r, \{\mu_i\})$ with the same value of $A$ and the same value of $r$.

(b) The dominant terms in the $n \to 0$ limit correspond to the minimum of $A(r, \{\mu_i\})$ and not the maximum. This seems an unreasonable assumption, but gives the correct result when replica symmetry is broken. One argument to support this assumption is that when the number of independent parameters we are maximising over is negative the maximum becomes a minimum [5]. In (34) there are $r-1$ independent parameters $\mu_i$ (due to the constraint $\sum_{i=1}^r \mu_i = n$) and, as we will see in (44) below, $r-1$ is negative when $n < \frac{\beta_c}{\beta}$.

(c) We allow $n, r, \mu_i$ to become real parameters when we compute the minimum of $A(r, \{\mu_i\})$.

As for Parisi’s original ansatz, it is clear that these assumptions, as such, have no rigorous justification. However, the free energy obtained using these assumptions has been verified for a number of spin glass models by a rigorous mathematical analysis (for reviews see [32, 33]). They also lead to the correct free energy of the REM in the low temperature phase [14]. In addition, there is a representation of the free energy as a contour integral in the complex plane [24, 27], where the minimum of $A(r, \{\mu_i\})$ and the non-integer values of $r$ and the $\mu_i$ appear naturally as a consequence of taking the saddle point along the contour.
The minimum of $A(r, \{\mu_i\})$ in (34) with respect to the $\{\mu_i\}$ subject to the constraint $\sum_{i=1}^r \mu_i = n$ can be found using a Lagrange multiplier and it corresponds to all $\mu_i$ taking the same value. The constraint then gives immediately

$$ \mu_i = \frac{n}{r}. \qquad (42) $$

for all $i$. Then (34) gives $A(r, \{\mu_i\}) = rf(\frac{n\beta}{r})$ and taking the extremal value with respect to $r$ one gets

$$ f\left(\frac{\beta n}{r}\right) - \frac{\beta n}{r} f'\left(\frac{\beta n}{r}\right) = 0. \qquad (43) $$

Comparison with (41) gives

$$ r = n \frac{\beta}{\beta_c}, \quad \text{so that} \quad \mu = \frac{\beta_c}{\beta} \qquad (44) $$

so that for large $N$ we can approximate (29) as

$$ \langle Z(\beta)^n \rangle \sim \exp \left\{ Nn \left[ \frac{1}{\mu} \log 2 + \frac{\beta^2}{4\mu} \right] \right\} = \exp \left\{ Nn \frac{\beta}{2\beta_c} \right\}, \qquad (45) $$

where we have defined $\mu = \frac{\beta_c}{\beta}$. The free energy (28) in the frozen phase is therefore $\langle \log Z \rangle = N \frac{\beta}{2\beta_c}$ which is known to be the correct expression [13]

The extremal condition (42) tells us that $\mu_i = \mu$ and $r = \frac{n}{\mu}$ for the dominant terms. So the $\mu_i$ do not fluctuate and we can immediately write

$$ \langle Y_k \rangle_n = r \frac{\mu(\mu-1)\cdots(\mu-k+1)}{n(n-1)\cdots(n-k+1)} = \frac{\Gamma(k-\mu)\Gamma(1-n)}{\Gamma(1-\mu)\Gamma(k-n)} \qquad (46) $$

which in the $n \to 0$ limit gives the well known result (9) and confirms that $\langle Y_2 \rangle = 1 - \mu$. We are now going to see that in the two temperature case this simple last step is not possible because fluctuations of $\mu_i$ remain even in the thermodynamic limit.

## 3.2. The REM at two temperatures

### 3.2.1. Integer moments of the partition function.
In the two temperature case our starting point is the following expression for the moments (see appendix B)

$$ \langle Z(\beta)^n Z(\beta')^{n'} \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\{\mu_i\}} \sum_{\{\mu_i'\} \ge 0} \delta[\mu_i + \mu_i' \ge 1] \\ \times C_{n,r}(\{\mu_i\}) C_{n',r}(\{\mu_i'\}) e^{NA(r, \{\mu_i, \mu_i'\})}, \qquad (47) $$

where the sum on $\{\mu_i\}, \{\mu_i'\}$ is over all non-negative integers, the $C_{n,r}(\{\mu_i\})$ and $C_{n',r}(\{\mu_i'\})$ are combinatorial factors defined in (30) and

$$ A(r, \{\mu_i, \mu_i'\}) = \sum_{i=1}^{r} f(\mu_i \beta + \mu_i' \beta') = r \log 2 + \frac{J^2}{4} \sum_{i=1}^{r} (\mu_i \beta + \mu_i' \beta')^2. \qquad (48) $$

As in (33) each term in the sum (47) corresponds to a different grouping of the $n + n'$ replicas: in configuration $i$ there are $\mu_i$ replicas at inverse temperature $\beta$ and $\mu_i'$ replicas at inverse
temperature $\beta'$. There is an additional constraint associated with each configuration $i$ that $\mu_i + \mu'_i \ge 1$; in other words we need at least one replica, which can be from either $n$ or $n'$, in each configuration.

We are interested in how the single temperature overlap calculation leading to (46) generalises to the two temperature case. We start with the weighted form of $Y_{k,k'}$ defined in (22)

$$ \langle Y_{k,k'} \rangle_{n,n'} = \frac{\langle Y_{k,k'} Z(\beta)^n Z(\beta')^{n'} \rangle}{\langle Z(\beta)^n Z(\beta')^{n'} \rangle} = \frac{\langle Z(\beta k + \beta' k') Z(\beta)^{n-k} Z(\beta')^{n'-k'} \rangle}{\langle Z(\beta)^n Z(\beta')^{n'} \rangle}. \quad (49) $$

The denominator in the rightmost expression is the two temperature moment in (47). The numerator is a three temperature moment of the partition function (95) computed in appendix B. By a direct generalisation of the derivation (96) and (97) of (36) one gets

$$ \langle Y_{k,k'} \rangle_{n,n'} = \left\langle r \frac{\mu_1(\mu_1 - 1) \cdots (\mu_1 - k + 1)}{n(n-1) \cdots (n-k+1)} \frac{\mu'_1(\mu'_1 - 1) \cdots (\mu'_1 - k + 1)}{n'(n'-1) \cdots (n'-k'+1)} \right\rangle_{\{\mu_i, \mu'_i\}}, \quad (50) $$

where the average $\langle \cdot \rangle_{\{\mu_i, \mu'_i\}}$ means that for any function $F(\{\mu_i, \mu'_i\})$,

$$ \langle F(\{\mu_i, \mu'_i\}) \rangle_{\{\mu_i, \mu'_i\}} = \frac{\sum_{r \ge 1} \sum_{\{\mu_i \ge 0\}} \sum_{\{\mu'_i \ge 0\}} F(\{\mu_i, \mu'_i\}) W_r(\{\mu_i, \mu'_i\})}{\sum_{r \ge 1} \sum_{\{\mu_i \ge 0\}} \sum_{\{\mu'_i \ge 0\}} W_r(\{\mu_i, \mu'_i\})} \quad (51) $$

with (see (30) and (48))

$$ W_r(\{\mu_i, \mu'_i\}) = \frac{C_{n,r}(\{\mu_i\}) C_{n',r}(\{\mu'_i\})}{r!} e^{NA(r, \{\mu_i, \mu'_i\})} \prod_{i=1}^{r} \theta[\mu_i + \mu'_i \ge 1]. $$

Here $\theta[\mu_i + \mu'_i \ge 1]$ is one if $\mu_i + \mu'_i \ge 1$ and zero otherwise.

### 3.2.2. The thermodynamic limit and the extremal condition.
Let us focus on the low temperature phase and take $\beta > \beta_c$ and $\beta' > \beta_c$ when replica symmetry is broken. We proceed as we did in the single temperature case, by making a similar set of three assumptions on how to take the thermodynamic limit as in the single temperature case. As before we look for the minimum of $A(r, \{\mu_i, \mu'_i\})$ in (48) and the only difference is that we now have the additional parameters $n'$ and $\{\mu'_i\}$. Using Lagrange multipliers we find that the minimum corresponds to $\beta\mu_i + \beta'\mu'_i$ being independent of $i$. Summing on $i$ and using the constraints $\sum_{i=1}^r \mu_i = n$ and $\sum_{i=1}^r \mu'_i = n'$ we find that

$$ \beta\mu_i + \beta'\mu'_i = \frac{\beta n + \beta'n'}{r} \quad \text{for all } i. \quad (52) $$

The value of $r$ that gives the minimum of $A(r, \{\mu_i, \mu'_i\}) = rf(\frac{n\beta+n'\beta'}{r})$ is then given by

$$ f\left(\frac{\beta n + \beta' n'}{r}\right) - \frac{\beta n + \beta' n'}{r} f'\left(\frac{\beta n + \beta' n'}{r}\right) = 0. \quad (53) $$

so that from (41)

$$ r = \frac{\beta n + \beta' n'}{\beta_c}. \quad (54) $$
Together with equation (52) this gives

$$ \beta\mu_i + \beta'\mu'_i = \beta_c, \qquad (55) $$

which constrains the fluctuations of $\mu_i$ and $\mu'_i$, but, unlike the single temperature case, (55) does not eliminate them completely.

One can however, without any further assumption, recover (25) from (50). Using the replica form (50) one can see that

$$ \begin{aligned} & (k-n)\beta\langle Y_{k+1,k'} \rangle_{n,n'} + (k'-n')\beta'\langle Y_{k,k'+1} \rangle_{n,n'} \\ & = \left\langle (\beta k + \beta' k' - \beta \mu_1 - \beta' \mu'_1) \right. \\ & \qquad \left. \times r \frac{\mu_1(\mu_1-1)\cdots(\mu_1-k+1)}{n(n-1)\cdots(n-k+1)} \frac{\mu'_1(\mu'_1-1)\cdots(\mu'_1-k'+1)}{n'(n'-1)\cdots(n'-k'+1)} \right\rangle_{\mu_i, \mu'_i}. \end{aligned} \qquad (56) $$

If we take the large $N$ limit of the right-hand side we expect (see (55)) that the extremal condition $\beta\mu_1 + \beta'\mu'_1 = \beta_c$ should apply. Then (56) simplifies to give

$$ (k-n)\beta\langle Y_{k+1,k'} \rangle_{n,n'} + (k'-n')\beta'\langle Y_{k,k'+1} \rangle_{n,n'} = (k\beta + k'\beta' - \beta_c)\langle Y_{k,k'} \rangle_{n,n'} \qquad (57) $$

the same recursion relation (25) as in the direct calculation. This gives at least some confidence in the assumptions that we have made in developing the replica approach so far for the two temperature problem.

Ideally we would like to go a step further and recover the exact solution (23) directly from (50) using the replica approach. The challenge is to find a way to compute the average over $\{\mu_i, \mu'_i\}$ in (50) subject to the constraints (54) and (55) that is valid when $n, n'$ are no-longer integers. We have not found an approach that is sufficiently convincing to merit inclusion here. The problem essentially has to do with the $r$ and the $\{\mu_i, \mu'_i\}$ becoming non-integer. We will see however in section 5 that by matching with the exact expressions of section 2 one can obtain the generating function of the $\{\mu_i, \mu'_i\}$.

## 4. Parisi overlap matrices

For the REM the replica method in section 3 could be implemented without explicitly using the replica overlap matrices. However, for more complex problems such as the Sherringon-Kirkpatrick model [7] the saddle point equations are expressed in terms of replica overlap matrices and so it is useful to see what they look like in the case of the REM. In this section we describe the structure of these matrices in the single and two temperature case of the REM. The Parisi ansatz is used in the single temperature case and we show how it can be generalised to two temperatures in the case of the REM.

We could have approached this by applying the replica method to the large $p$ limit of the p-spin models introduced in [14] as was done in [17] for the single temperature case. However, the corresponding two temperature calculation is rather long and is not essential to understanding how to generalise the Parisi ansatz.

### 4.1. Single temperature case

In the single temperature case $\langle Z(\beta)^n \rangle$ is expressed in equation (33) as a sum over the parameters $r$ and $\mu_1, \mu_2, \dots, \mu_r$ with the constraint $\sum_{i=1}^r \mu_i = n$. We can use these parameters to
define an $n \times n$ replica overlap matrix $\mathbf{Q}(\{\mu_i\})$. We divide the $n$ replicas into $r$ groups of sizes $\mu_1, \mu_2, \dots, \mu_r$. The overlap matrix is then defined as

$$ \mathcal{Q}_{a,b}(\{\mu_i\}) = \begin{cases} 1, & \text{if replicas } a, b \text{ are in the same group,} \\ 0, & \text{otherwise.} \end{cases} \qquad (58) $$

This means that, up to a permutation of the replica indices, the $n \times n$ matrix $\mathbf{Q}(\{\mu_i\})$ consists of $r$ blocks of size $\mu_1 \times \mu_1, \mu_2 \times \mu_2, \dots, \mu_r \times \mu_r$, along the diagonal where the matrix elements take value unity and they are zero elsewhere. For example if $n = 6, r = 3, \mu_1 = 2, \mu_2 = 3, \mu_3 = 1$

$$ \mathbf{Q}(\{\mu_i\}) = \begin{pmatrix} \mathbf{1} & \mathbf{1} & 0 & 0 & 0 & 0 \\ \mathbf{1} & \mathbf{1} & 0 & 0 & 0 & 0 \\ 0 & 0 & \mathbf{1} & \mathbf{1} & \mathbf{1} & 0 \\ 0 & 0 & \mathbf{1} & \mathbf{1} & \mathbf{1} & 0 \\ 0 & 0 & \mathbf{1} & \mathbf{1} & \mathbf{1} & 0 \\ 0 & 0 & 0 & 0 & 0 & \mathbf{1} \end{pmatrix}. \qquad (59) $$

(Here we have taken $Q_{a,a} = 1$ for simplicity). In terms of this overlap matrix, one can rewrite (34) as

$$ A(r, \{\mu_i\}) = r \log 2 + \frac{(\beta J)^2}{4} \sum_{a=1}^{n} \sum_{b=1}^{n} Q_{a,b}(\{\mu_i\}), \qquad (60) $$

where we have used that $\sum_{a=1}^{n} \sum_{b=1}^{n} Q_{a,b}(\{\mu_i\}) = \sum_{i=1}^{r} \mu_i^2$.

The thermodynamic limit gives us the extremal condition (42) which indicates that the dominant form of $\mathbf{Q}(\{\mu_i\})$ has all the $\mu_i$ equal. This fixed block structure gives the one step RSB form of the overlap matrices introduced by Parisi [1] to solve mean field spin glass models such as the Sherrington-Kirkpatrick model [7]. It should be noted that any permutation of the replica indices will also give a $\mathbf{Q}(\{\mu_i\})$ that satisfies the extremal condition and we should sum over all these saddle points when computing physical properties such as $P(q)$ (see [34]).

## 4.2. Two temperature case

In the two temperature case $\langle Z(\beta)^n Z(\beta')^{n'} \rangle$ is expressed in equation (47) as a sum over the parameters $r, \mu_1, \mu_2, \dots, \mu_r$ and $\mu'_1, \mu'_2, \dots, \mu'_r$ with the constraints $\sum_{i=1}^r \mu_i = n$ and $\sum_{i=1}^r \mu'_i = n'$. We can use these parameters to define three different replica overlap matrices. We divide the $n$ replicas into $r$ groups of size $\mu_1, \mu_2, \dots, \mu_r$ and the $n'$ replicas into $r$ groups of size $\mu'_1, \mu'_2, \dots, \mu'_r$. We then have the single temperature $n \times n$ replica overlap matrix $\mathbf{Q}(\{\mu_i\})$ defined in equation (58) and the equivalent $n' \times n'$ matrix $\mathbf{Q}'(\{\mu'_i\})$ at inverse temperature $\beta'$. We can also define an $n \times n'$ overlap matrix $\mathbf{R}(\{\mu_i, \mu'_i\})$ between the inverse temperature $\beta$ and the inverse temperature $\beta'$ by

$$ R_{a,b'}(\{\mu_i, \mu'_i\}) = \begin{cases} 1, & \text{if replicas } a, b' \text{ are in the same group;} \\ 0, & \text{otherwise.} \end{cases} \qquad (61) $$

So all the matrix elements of the rectangular matrix $\mathbf{R}(\{\mu_i, \mu'_i\})$ are zero except $r$ blocks of sizes $\mu_1 \times \mu'_1, \mu_2 \times \mu'_2, \dots, \mu_r \times \mu'_r$, along the diagonal where they take value unity. It has
## 5. The fluctuations in block sizes $\mu_i, \mu'_i$

In this section we analyse the fluctuations of the block sizes $\mu_i, \mu'_i$ in the thermodynamic limit. We first obtain the mean and the variance of the $\mu_i$ and $\mu'_i$. We will then obtain the moment generating function for the distribution of these block sizes from the exact expression (23) for $\langle Y_{k,k'} \rangle_{n,n'}$. One outcome of our results is that the $\mu_i$ and the $\mu'_i$ do fluctuate even in the equal temperature case. However extracting the distribution of $P(\mu_i, \mu'_i)$ from this generating function is not an easy task and can be interpreted as a signed measure, i.e. a measure with negative probabilities. Also not all properties of the distribution of these block sizes have physical implications: for example, we will see that in the limit $n \to 0^-$ and $n' \to 0^-$ these distributions depend on the ratio $n'/n$ although all physical properties have a limit independent of this ratio.

In this section $n, n'$ are negative real numbers because our analysis is based on the exact expression (23) which is only valid for this range of values.

### 5.1. The first moments of $\mu_i$

As $Y_{1,0} = 1$ (see (12)) and as $r$ does not fluctuate (see (54)) one can show from (50) that

$$ \langle \mu_i \rangle_{\{\mu_i, \mu'_i\}} = \frac{n\beta_c}{n\beta + n'\beta'} \\ \langle \mu_i^2 \rangle_{\{\mu_i, \mu'_i\}} = \frac{n\beta_c(1 - (1-n)\langle Y_{2,0} \rangle_{n,n'})}{n\beta + n'\beta'} = \frac{n\beta_c(\beta_c - n'\beta'\langle Y_{1,1} \rangle_{n,n'})}{\beta(n\beta + n'\beta')} \quad (64) $$

where we have used the relation between $\langle Y_{2,0} \rangle_{n,n'}$ and $\langle Y_{1,1} \rangle_{n,n'}$

$$ \langle Y_{2,0} \rangle_{n,n'} = \frac{\beta - \beta_c + n'\beta' \langle Y_{1,1} \rangle_{n,n'}}{(1-n)\beta} $$

which follows from (25) and the fact that $Y_{1,0} = 1$.

From these expressions (64) one can notice first that the limit of the first moment $\langle \mu_i \rangle$, when $n \to 0^-$ and $n' \to 0^-$, depends on the ratio $n'/n$. This means that not all the properties of the $\mu_i$ have a physical meaning, since one expects all physical properties to be independent of this ratio when $n$ and $n'$ vanish.

One can also notice that the variance of $\mu_i$ is in general non-zero. Depending on $n, n', \beta, \beta'$, this variance may change its sign, implying that the distribution of $\mu_i$ is not really a probability distribution. For example when $\beta = \beta'$, one has (27)

$$ \langle Y_{1,1} \rangle_{n,n'} = \langle Y_{2,0} \rangle_{n,n'} = \frac{1 - \frac{\beta_c}{\beta}}{1 - n - n'} $$

and

$$ \langle \mu_i^2 \rangle_{\{\mu_i, \mu'_i\}} - \langle \mu_i \rangle_{\{\mu_i, \mu'_i\}}^2 = \frac{nn'\beta_c(\beta(n+n') - \beta_c)}{\beta^2(n+n')^2(n+n'-1)} $$

which in general does not vanish and can be of either sign. We already observed such negative variances of the block sizes when we tried to reproduce, using replicas, finite size corrections of the REM at a single temperature [24].

Expressions similar to (64) can be obtained for $\mu'_i$ by using either the symmetry $n, \beta \leftrightarrow n'\beta'$ or the fact that the sum $\beta\mu_i + \beta'\mu'_i = \beta_c$ does not fluctuate (55).
## 5.2. The generating function of $\mu_i$ and $\mu_i'$

We are now going to obtain the exact expression of the generating function $\langle x^{\mu_1} y^{\mu_1'} \rangle_{\{\mu_i, \mu_i'\}}$

where the average $\langle \cdot \rangle_{\{\mu_i, \mu_i'\}}$ is defined in equation (51). Taking the Taylor expansions of $x^{\mu_1}$, $y^{\mu_1'}$ about $x = 1, y = 1$ the generating function can be written as

$$ \langle x^{\mu_1} y^{\mu_1'} \rangle_{\{\mu_i, \mu_i'\}} = \sum_{k \ge 0} \sum_{k' \ge 0} \frac{(x-1)^k (y-1)^k}{k! k'!} \times \langle \mu_1(\mu_1-1) \cdots (\mu_1-k+1) \mu_1'(\mu_1'-1) \cdots (\mu_1'-k+1) \rangle_{\{\mu_i, \mu_i'\}} \quad (65) $$

In the thermodynamic limit we can express the average on $\{\mu_i, \mu_i'\}$ on the right-hand side in terms of $\langle Y_{k,k'} \rangle_{n,n'}$ using equation (50). This gives

$$ \langle \mu_1(\mu_1 - 1) \cdots (\mu_1 - k + 1) \quad \mu_1'(\mu_1' - 1) \cdots (\mu_1' - k + 1) \rangle_{\{\mu_i, \mu_i'\}} \\ = \frac{1}{r} n'(n' - 1) \cdots (n' - k + 1) n(n - 1) \cdots (n - k + 1) \langle Y_{k,k'} \rangle_{n,n'}, \quad (66) $$

where we have used the fact that $r$ does not fluctuate in the thermodynamic limit (see (24)). Using the exact expression (23) for $\langle Y_{k,k'} \rangle_{n,n'}$ we obtain

$$ \langle \mu_1(\mu_1 - 1) \cdots (\mu_1 - k + 1) \quad \mu_1'(\mu_1' - 1) \cdots (\mu_1' - k + 1) \rangle_{\{\mu_i, \mu_i'\}} \\ = (-1)^{k+k'} \frac{\int_0^\infty dv \, v^{k'-n'-1} \Psi(k+k'\frac{\beta}{\beta}; v) (-\psi(v))^{r-1}}{\int_0^\infty dv \, v^{-n'-1} (-\psi(v))^r}. \quad (67) $$

Finally, substituting into equation (65) and summing on $k, k'$ we obtain

$$ \langle x^{\mu_1} y^{\mu_1'} \rangle_{\{\mu_i, \mu_i'\}} = \frac{x^{\frac{\beta_c}{\beta}} \int_0^\infty dv \, v^{-n'-1} \left( -\psi(v) \frac{v}{y} x^{-\frac{\beta'}{\beta}} \right) (-\psi(v))^{r-1}}{\int_0^\infty dv \, v^{-n'-1} (-\psi(v))^r}, \quad (68) $$

where $r$ is given by (24). (Note that, as mentioned at the beginning of this section the above expression (65) is valid for $n < 0$ and $n' < 0$.)

### Remark 5.1. One recovers (55)

Making the substitution $x = z^\beta, y = z^{\beta'}$ in (68) gives

$$ \langle z^{\beta\mu_1 + \beta'\mu_1'} \rangle_{\{\mu_i, \mu_i'\}} = z^{\beta_c}. \quad (69) $$

This confirms the fact that the sum $\beta\mu_1 + \beta'\mu_1'$ does not fluctuate and takes the value $\beta_c$, as expected from (55).

### Remark 5.2. The $n, n' \to 0^-$ limit
One can show using the asymptotics (85) that, in the limit $n \to 0^-$ and $n' \to 0^-$, the generating function (68) becomes

$$ \langle x^{\mu_1} y^{\mu'_1} \rangle_{\{\mu_i, \mu'_i\}} = \frac{\beta n x^{\frac{\beta c}{\beta}} + \beta' n' y^{\frac{\beta c}{\beta}}}{\beta n + \beta' n'} + \frac{\beta \beta' n n' (y^{\frac{\beta c}{\beta}} - x^{\frac{\beta c}{\beta}})}{\beta_c (\beta n + \beta' n')} \log \left( \frac{\Gamma\left(1 - \frac{\beta c}{\beta}\right)}{\Gamma\left(1 - \frac{\beta c}{\beta}\right)} \right) \\ + \frac{\beta n n'}{\beta n + \beta' n'} x^{\frac{\beta c}{\beta}} \int_0^{\infty} \log \frac{dv}{v} \frac{\psi(yx^{-\frac{\beta'}{\beta}} v)}{\psi(v)} dv + o(n, n'). \quad (70) $$

To leading order one finds that the distribution of $\mu_1$ and $\mu'_1$ consists of two delta functions

$$ P(\mu_1, \mu'_1) = \left[ \frac{n\beta}{n\beta + n'\beta'} \delta\left(\mu_1 - \frac{\beta_c}{\beta}\right) + \frac{n'\beta'}{n\beta + n'\beta'} \delta(\mu_1) \right] \beta' \delta(\beta\mu_1 + \beta'\mu'_1 - \beta_c). $$

Clearly this expression does not contain any information on the overlaps $\langle Y_{k,k'} \rangle$. In fact the generating function of these overlaps only appears in the first order term in (70). Note also that as for the variance of $\mu_1$, the $n \to 0^-$, $n' \to 0^-$ limit depends on the ratio $n'/n$.

### 5.3. Trying to describe $P(\mu_1, \mu'_1)$ the probability distribution

As the variance of the $\mu_i$ can become negative, it is clear from the very beginning that it is not possible to find a meaningful distribution of the block sizes compatible with the generating function (68). We made a number of attempts which became rather complicated and we do not think it is of much interest to mention them here. Let us however discuss briefly one case for which we could get a rather simple picture, the equal temperature case, i.e. when $\beta = \beta'$.

In this case we have an explicit expression (79) of the function $\psi(v)$. Then (68) becomes

$$ \langle x^{\mu_1} y^{\mu'_1} \rangle_{\{\mu_i, \mu'_i\}} = \frac{\int_0^{\infty} v^{-n'-1} (1+v)^{n+n'-\frac{\beta c}{\beta}} (x+vy)^{\frac{\beta c}{\beta}} dv}{\int_0^{\infty} v^{-n'-1} (1+v)^{n+n'} dv}. $$

After a simple change of variable $v = (1 - t)/t$ this becomes

$$ \langle x^{\mu_1} y^{\mu'_1} \rangle_{\{\mu_i, \mu'_i\}} = \frac{\Gamma(-n-n')}{\Gamma(-n)\Gamma(-n')} \int_0^1 dt \, t^{-1-n} (1-t)^{-1-n'} (tx + (1-t)y)^{\frac{\beta c}{\beta}}. \quad (71) $$

In order to give an interpretation to (71) let us consider a random variable $s$, sum of $m$ i.i.d. random variables $\tau_i$ which take the value $\tau_i = 1$ with probability $t$ and $\tau_i = 0$ with probability $1-t$. The distribution of $s$ is a binomial distribution and one has

$$ \langle z^s \rangle = (zt + 1 - t)^m. \quad (72) $$

Let us further consider that the parameter $t$ is itself randomly distributed according to some distribution $\rho(t)$ so that the distribution of $s$ becomes a superposition of binomial distributions. Then the generating function of $s$ becomes

$$ \langle z^s \rangle = \int_0^1 \rho(t)(zt + 1 - t)^m dt. $$
This is exactly the form we have in (71) (by taking $x = z$ and $y = 1$) if one chooses for $\rho(t)$

$$ \rho(t) = \frac{\Gamma(-n-n')}{\Gamma(-n)\Gamma(-n')} t^{-1-n} (1-t)^{-1-n'} $$

(remember that here $n$ and $n'$ are negative).

Therefore the distribution of $\mu_1$ can be thought as a superposition of binomial distributions. The only odd aspect is that $s$ is a sum of $m = \frac{\beta}{\beta'}$ binary variables, that is $s$ is a sum of a non-integer number of random variables!

**Remark 5.3.** A signed measure

If one takes a non-integer $m$ in (72) one gets by expanding in powers of $z$

$$ \langle z^s \rangle = \sum_{p=0} (1-t)^{m-p} t^p \frac{m(m-1)\cdots(m-p+1)}{p!} z^p $$

which one can interpret, for $t < \frac{1}{2}$, as the probability $P(s)$ of $s$ being a signed measure concentrated on positive integers. Expanding in powers of $1/z$ leads, for $\frac{1}{2} < t < 1$, to a different signed measure. Combining these two representations by cutting the integral (71) into two parts ($t < \frac{1}{2}$ and $t > \frac{1}{2}$) leads to a signed measure concentrated on the points $(\mu_i = p, \mu'_i = \frac{\beta}{\beta'} - p)$ and $(\mu_i = \frac{\beta}{\beta'} - p, \mu'_i = p)$ for all positive integers $p \ge 0$.

## 6. Conclusion

In this paper we have analysed the distribution of overlaps (11) between two copies of the same REM at two temperatures. A direct calculation was used to obtain exact expressions (23) for the two temperature overlaps (12) in the thermodynamic limit. Generalising this approach allows us to quantify (17) the non-self-averaging effects illustrated in figure 1.

An alternative approach using the replica method enables us to obtain expressions for the two temperature overlaps in terms of replicas (50). In the thermodynamic limit the exact and replica expressions satisfy the same Ghirlanda–Guerra type recurrence relation, (25) and (57), giving confidence that the replica expressions are valid. We also proposed a way to generalise the Parisi ansatz (63), in the one step RSB form, to the two temperature case which is consistent with the replica expressions for the overlaps. In contrast to the single temperature case we find that the block sizes at the two different temperatures fluctuate even in the thermodynamic limit subject the constraint (55). We characterised these fluctuations in terms of a moment generating function for the block sizes (68). It is well known that in the single temperature case the strange properties of Parisi’s RSB ansatz (non-integer number of blocks and block sizes) lead to a perfectly good physical description in terms of overlaps [3]. In the two temperature case the distribution of $\mu_i$ and $\mu'_i$ analysed in section 5 also has strange properties, but it remains an open question as to which of these properties have a clear physical interpretation.

It would be interesting to extend both the exact and replica approaches to the generalised random energy model [35], directed polymer in a random medium [36] and other models where exact methods are likely to be tractable. In contrast to the single temperature case, the multi-temperature overlaps should be different in the REM and in the directed polymer problem on a tree because the lowest energies of the directed polymer can be thought of as a decorated Poisson process and it has been proved that the multi-temperature overlaps depend on the decoration [37]. One could also look at spin models where one step RSB occurs to see if the two temperature ansatz with fluctuating block sizes is applicable. An obvious starting point
would be the p-spin spherical model proposed in [38]. In order to address these spin problems, where exact expressions for the two temperature overlaps are not currently available, it would be essential to develop a systematic approach to summing over the fluctuations in block sizes in the replica expressions.

## Appendix A. Direct calculation of the overlaps (13) and (23)

To begin with, it is easier to think that the energies can take only a discrete set of values $E_ν$ indexed by $ν$ and that the partition function at inverse temperature $β$ is given by

$$Z(β) = \sum_{ν} y_{ν} e^{-βE_{ν}},$$

where

$$y_{\nu} = \begin{cases} 1 & \text{with probability } p_{\nu} \\ 0 & \text{with probability } 1 - p_{\nu} \end{cases}.$$

So a given sample is specified by the value of all these binary random variables $y_ν$. Then the probability of finding $k$ copies at temperature $β$ and $k'$ copies at temperature $β'$ in the same configuration is given by

$$Y_{k,k'} = \frac{\sum_{\nu} y_{\nu} e^{-(\beta k + \beta' k') E_{\nu}}}{Z(\beta)^k Z(\beta')^{k'}}.$$

These $Y_{k,k'}$ are random quantities as they depend on the realization of the $y_ν$'s. Using the identity $Z^{-k} = \Gamma(k)^{-1} \int_0^\infty dt e^{-tZ} t^{k-1}$ one gets

$$Y_{k,k'} = \sum_{\nu} y_{\nu} e^{-(\beta k + \beta' k') E_{\nu}} \int_{0}^{\infty} \frac{t^{k-1} dt}{\Gamma(k)} \times \int_{0}^{\infty} \frac{t^{k'-1} dt'}{\Gamma(k')} \exp \left[ -\sum_{\nu'} y_{\nu'} \left( t e^{-\beta E_{\nu'}} + t' e^{-\beta' E_{\nu'}} \right) \right].$$

Averaging over the $y_ν$'s leads to

$$\langle Y_{k,k'} \rangle = \int_0^\infty \frac{t^{k-1} dt}{\Gamma(k)} \int_0^\infty \frac{t^{k'-1} dt'}{\Gamma(k')} \sum_\nu p_\nu e^{-(\beta k + \beta' k') E_\nu} \times \exp \left[ -t e^{-\beta E_\nu} - t' e^{-\beta' E_\nu} \right] \prod_{\nu' \neq \nu} \left( 1 - p_{\nu'} + p_{\nu'} \exp \left[ -t e^{-\beta E_{\nu'}} - t' e^{-\beta' E_{\nu'}} \right] \right).$$

Now if we go to the continuum limit, by saying that each energy interval $(E, E + dE)$ is either occupied by an energy level or empty and if we choose as in (20)

$$p_{\nu} = C' e^{\beta_c E} dE \quad \text{with} \quad C' = C e^{-\beta_c E_0}$$

one gets

$$\langle Y_{k,k'} \rangle = \int_0^\infty \frac{t^{k-1} dt}{\Gamma(k)} \int_0^\infty \frac{t^{k'-1} dt'}{\Gamma(k')} W(k, k'; t, t') e^{w(t,t')},$$
where

$$w(t, t') = C' \int e^{\beta_c E} dE \left( \exp \left[ -te^{-\beta E} - t'e^{-\beta' E} \right] - 1 \right)$$

and

$$W(k, k'; t, t') = C' \int e^{\beta_c E} dE e^{-(\beta k + \beta' k')} \exp \left[ -te^{-\beta E} - t'e^{-\beta' E} \right].$$

Then these expressions can be simplified by noticing that

$$w(t, t') = \frac{C'}{\beta} t^{\frac{\beta_c}{\beta}} \psi\left(\frac{t'}{t^{\frac{\beta'}{\beta}}}\right)$$

and

$$W(k, k'; t, t') = \frac{C'}{\beta} t^{\frac{\beta_c}{\beta} - k - \frac{\beta'}{\beta} k'} \Psi\left(k + k' \frac{\beta'}{\beta}; \frac{t'}{t^{\frac{\beta'}{\beta}}}\right),$$

where

$$\psi(v) = \int_{0}^{\infty} du \left( e^{-u-vu^{\frac{\beta'}{\beta}}} - 1 \right) u^{1-\frac{\beta_c}{\beta}} \quad (73)$$

and

$$\Psi(z; v) = \int_{0}^{\infty} du \, e^{-u-vu^{\frac{\beta'}{\beta}}} u^{z-1-\frac{\beta_c}{\beta}}. \quad (74)$$

This leads to

$$\langle Y_{k,k'} \rangle = \frac{\beta}{\beta_c} \frac{1}{\Gamma(k)} \frac{1}{\Gamma(k')} \int_0^\infty dv \, v^{k'-1} \frac{\Psi\left(k + k' \frac{\beta'}{\beta}; v\right)}{(-\psi(v))} \quad (75)$$

**Remark A.1.** To generalize (75) one can define the probability of finding $k_1$ copies at inverse temperature $\beta$ and $k'_1$ copies of the same system at inverse temperature $\beta'$ in the same configuration, and similarly $k_2$ and $k'_2$ in a different configuration and so on i.e.

$$Y_{k_1 k'_1 \cdots k_p k'_p} = \frac{\sum_{\nu_1 \cdots \nu_p} e^{-(\beta k_1 + \beta' k'_1)E_{\nu_1} - \cdots - (\beta k_p + \beta' k'_p)E_{\nu_p}}}{Z(\beta)^{k_1 + \cdots + k_p} Z(\beta')^{k'_1 + \cdots + k'_p}},$$

where, in the sum, the configurations $\nu_1 \neq \nu_2 \neq \cdots \neq \nu_p$ are all different. By a straightforward extension of the above calculation one gets

$$\langle Y_{k_1 k'_1 \cdots k_p k'_p} \rangle = \frac{\beta}{\beta_c} \Gamma(p) \frac{1}{\Gamma(k_1 + \cdots + k_p)} \frac{1}{\Gamma(k'_1 + \cdots + k'_p)} \times \int_0^\infty dv \, v^{k'_1 + \cdots + k'_p - 1} \frac{\Psi\left(k_1 + k'_1 \frac{\beta'}{\beta}; v\right) \cdots \Psi\left(k_p + k'_p \frac{\beta'}{\beta}; v\right)}{(-\psi(v))^p} \quad (76)$$

**Remark A.2.** In the way the above formulae are written, $\beta$ and $\beta'$ seem to play asymmetric roles. One can however check from the definitions (73), and (74) of $\psi$ and $\Psi$ that

$$\psi_{\beta, \beta'}(v) = \frac{\beta}{\beta'} v^{\frac{\beta_c}{\beta'}} \psi_{\beta', \beta}(v^{-\frac{\beta}{\beta'}}) \quad (77)$$
$$ \Psi_{\beta, \beta'}(z; v) = \frac{\beta}{\beta'} v^{\frac{\beta_c}{\beta} - \frac{\beta}{\beta'}} z \Psi_{\beta', \beta} \left( \frac{\beta}{\beta'} z; v^{-\frac{\beta}{\beta'}} \right) \quad (78) $$

and using these relations one can easily prove that the expressions (75) and (76) are left unchanged by the symmetry

$$ (\beta, \beta', \{k_1, \dots, k_p\}, \{k'_1, \dots, k'_p\}) \longleftrightarrow (\beta', \beta, \{k'_1, \dots, k'_p\}, \{k_1, \dots, k_p\}). $$

**Remark A.3.** When $\beta = \beta'$, the expressions (73) and (74) become

$$ \psi(v) = \Gamma\left(-\frac{\beta_c}{\beta}\right) (1 + v)^{\frac{\beta_c}{\beta}} \quad (79) $$

$$ \Psi(v) = \Gamma\left(z - \frac{\beta_c}{\beta}\right) (1 + v)^{\frac{\beta_c}{\beta} - z}. \quad (80) $$

The integrals in (76) can then be performed and one gets

$$ \langle Y_{k_1, k'_1; \dots, k_p, k'_p} \rangle = \frac{\beta}{\beta_c} \frac{\Gamma(p)}{\Gamma(k_1 + k'_1 + \dots + k_p + k'_p)} \prod_{i=1}^{p} \left( \frac{\Gamma(k_i + k'_i - \frac{\beta_c}{\beta})}{-\Gamma(-\frac{\beta_c}{\beta})} \right) \quad (81) $$

which was already known (see for example [24]).

**Remark A.4.** It is easy to show, using an integration by parts in (74), that (for $z > \frac{\beta_c}{\beta}$)

$$ \left(z - \frac{\beta_c}{\beta}\right) \Psi(z; v) = \Psi(z + 1; v) + \frac{\beta'}{\beta} v \Psi\left(z + \frac{\beta'}{\beta}; v\right) $$

and this leads (see (76)) to relationships between the $\langle Y_{k_1, k'_1; \dots, k_p, k'_p} \rangle$

$$ (k_1 \beta + k'_1 \beta' - \beta_c) \langle Y_{k_1, k'_1; \dots, k_p, k'_p} \rangle = (k_1 + \dots + k_p) \beta \langle Y_{k_1+1, k'_1; \dots, k_p, k'_p} \rangle \\ + (k'_1 + \dots + k'_p) \beta' \langle Y_{k_1, k'_1+1; \dots, k_p, k'_p} \rangle \quad (82) $$

and similar identities for the pairs $k_2, k'_2, \dots, k_p, k'_p$.

**Remark A.5.** In the replica approach, one is usually interested in the limit where the number of replicas $n \to 0$. It is however often easier to first think in terms of a non-zero number of replicas and to take the $n \to 0$ limit afterwards. In this spirit, it is possible to generalize the above formulae (75) and (76) following a very similar calculation.

Defining the weighted overlaps for non-zero numbers $n$ and $n'$ of replicas ($n$ and $n'$ are a priori arbitrary real numbers) by

$$ \langle Y_{k,k'} \rangle_{n,n'} = \frac{\left\langle \sum_{\nu} y_{\nu} e^{-(\beta k + \beta' k') E_{\nu}} Z(\beta)^{n-k} Z(\beta')^{n-k'} \right\rangle}{\langle Z(\beta)^n \rangle \langle Z(\beta')^{n'} \rangle} $$
It turns out that the expressions have somewhat simpler forms when the numbers $n$ and $n'$ take negative values and one gets

$$ \langle Y_{k,k'} \rangle_{n,n'} = \frac{(-r) \Gamma(-n) \Gamma(-n')}{\Gamma(k-n) \Gamma(k'-n')} \times \frac{\int_0^\infty dv \, v^{k'-n'-1} \Psi(k+k' \frac{\beta'}{\beta}; v) (-\psi(v))^{r-1}}{\int_0^\infty dv \, v^{-n'-1} (-\psi(v))^r}, \quad (83) $$

where

$$ r = \frac{n\beta + n'\beta'}{\beta_c}. \quad (84) $$

**Remark A.6.** Using the following asymptotics of $\psi(v)$ which can be derived from expression (77)

$$ \psi(v) \approx \begin{cases} \Gamma\left(-\frac{\beta_c}{\beta}\right) & \text{as } v \to 0 \\ \frac{\beta}{\beta'} \Gamma\left(-\frac{\beta_c}{\beta'}\right) v^{\frac{\beta_c}{\beta'}} & \text{as } v \to \infty \end{cases} \quad (85) $$

one can show that (83) reduces to (75) in the limit $n \to 0^-$ and $n' \to 0^-$.

## Appendix B. Integer moments of the partition function at multiple temperatures

The replica method starts usually with the calculation of integer moments of the partition function. In a two or a multiple temperature case, these are of the form $\langle Z(\beta_1)^{n_1} Z(\beta_2)^{n_2} Z(\beta_3)^{n_3} \dots \rangle$ where $n_1, n_2, n_3, \dots$ are positive integers and $\beta_1, \beta_2, \beta_3, \dots$ are inverse temperatures. In this appendix we obtain the expressions (29) and (47) using a generating function defined for $p$ temperatures as

$$ G(t_1, t_2, \dots, t_p) = \left\langle \exp \left( -\sum_{i=1}^{p} t_i Z(\beta_i) \right) \right\rangle. \quad (86) $$

For the REM (see section 2) the partition function is given by

$$ Z(\beta) = \sum_{C=1}^{2^N} e^{-\beta E(C)}, $$

where the $2^N$ energies $E(C)$ take random values distributed according to $P(E)$ given in (18). Then, because the $E(C)$ are independent,

$$ G(t_1, \dots, t_p) = \left[ \int P(E) dE \exp \left( -\sum_i t_i e^{-\beta_i E} \right) \right]^{2^N} \\ = \exp \left\{ 2^N \log \left( \int P(E) dE \exp \left( -\sum_i t_i e^{-\beta_i E} \right) \right) \right\} $$

which for large $N$ becomes
$$G(t_1, \dots, t_p) \approx \exp \left\{ \int_{-\infty}^{\infty} \rho(E) \left[ \exp \left( -\sum_{i=1}^{p} t_i e^{-\beta_i E} \right) - 1 \right] dE \right\}. \quad (87)$$

(By this approximation, we in fact replace the REM by a Poisson REM of density (see (18))

$$\rho(E) = 2^N P(E) = 2^N \frac{1}{\sqrt{N\pi J^2}} \exp \left[ -\frac{E^2}{NJ^2} \right]. \quad (88)$$

Doing so the error is exponentially small in the system size N as shown in the appendix of [24]). The exponentials on the right-hand side of (87) can be expanded to obtain

$$G(t_1, \dots, t_p) = \sum_{r=0}^{\infty} \frac{1}{r!} \left[ \sum_{\substack{\mu_1=0 \\ (\mu_1+\dots+\mu_p \ge 1)}}^{\infty} \dots \sum_{\mu_p=0}^{\infty} \frac{(-t_1)^{\mu_1}}{\mu_1!} \dots \frac{(-t_p)^{\mu_p}}{\mu_p!} \right. \\ \left. \times \langle Z(\beta_1 \mu_1 + \dots + \beta_p \mu_p) \rangle \right]^r, \quad (89)$$

where we use the fact that for the REM (as well as for the Poisson REM)

$$\langle Z(\beta) \rangle = \int_{-\infty}^{\infty} \rho(E) e^{-\beta E} dE. \quad (90)$$

The general expression for integer moments at p temperatures is obtained by equating powers of $t_i$ in the expansion of the right-hand side of equation (86) with the right-hand side of equation (89). Here we give the three moments that are used in the main text.

The single temperature moments are then given by

$$\langle Z(\beta)^n \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\{\mu_i\} \ge 1} C_{n,r} (\{\mu_i\}) \langle Z(\beta\mu_1) \rangle \langle Z(\beta\mu_2) \rangle \dots \langle Z(\beta\mu_r) \rangle, \quad (91)$$

where we have defined

$$\sum_{\{\mu_i\} \ge 1} = \sum_{\mu_1 \ge 1} \sum_{\mu_2 \ge 1} \dots \sum_{\mu_r \ge 1} \quad (92)$$

and

$$C_{n,r}(\{\mu_i\}) = \frac{n!}{\mu_1! \mu_2! \dots \mu_r!} \delta \left[ \sum_{i=1}^{r} \mu_i = n \right]. \quad (93)$$

The Kronecker delta $\delta[\sum_{i=1}^{r} \mu_i = n]$ ensures that the $\mu_i$ always sum to n. Similarly the two temperature moments are given by

$$\langle Z(\beta)^n Z(\beta')^{n'} \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\{\mu_i\} \ge 0} \sum_{\{\mu_i'\} \ge 0} \theta[\mu_i + \mu_i' \ge 1] C_{n,r}(\{\mu_i\}) C_{n',r}(\{\mu_i'\}) \\ \times \langle Z(\beta\mu_1 + \beta'\mu_1') \rangle \langle Z(\beta\mu_2 + \beta'\mu_2') \rangle \dots \langle Z(\beta\mu_r + \beta'\mu_r') \rangle, \quad (94)$$
where $\theta[\mu_i + \mu_i' \ge 1]$ is unity if the inequality is satisfied for every $i = 1, 2, \dots, r$ and zero otherwise. The three temperature moments are given by

$$ \langle Z(\beta)^n Z(\beta')^{n'} Z(\beta'')^{n''} \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\{\mu_i \ge 0\}} \sum_{\mu_i' \ge 0} \sum_{\mu_i'' \ge 0} \theta[\mu_i + \mu_i' + \mu_i'' \ge 1] \\ \times C_{n,r}(\{\mu_i\}) C_{n',r}(\{\mu_i'\}) C_{n'',r}(\{\mu_i''\}) \langle Z(\beta\mu_1 + \beta'\mu_1' + \beta''\mu_1'') \rangle \\ \times \langle Z(\beta\mu_2 + \beta'\mu_2' + \beta''\mu_2'') \rangle \cdots \langle Z(\beta\mu_r + \beta'\mu_r' + \beta''\mu_r'') \rangle. \quad (95) $$

As a special case of (94) one has

$$ \langle Z(\beta)^{n-k} Z(k\beta) \rangle = \sum_{r \ge 1} \frac{1}{r!} \sum_{\{\mu_i \ge 0\}} \sum_{\mu_i' \ge 0} \theta[\mu_i + \mu_i' \ge 1] C_{n-k,r}(\{\mu_i\}) C_{1,r}(\{\mu_i'\}) \\ \times \langle Z(\beta\mu_1 + k\beta\mu_1') \rangle \langle Z(\beta\mu_2 + k\beta\mu_2') \rangle \cdots \langle Z(\beta\mu_r + k\beta\mu_r') \rangle. \quad (96) $$

In this case $n' = 1$, therefore there is a single $\mu_i' = 1$ all the others being 0. Because of the symmetry between the indices $i$ in the previous formula, one can choose $\mu_1' = 1$ and one gets

$$ \langle Z(\beta)^{n-k} Z(k\beta) \rangle = \sum_{r \ge 1} \frac{1}{(r-1)!} \sum_{\{\mu_i \ge 1\}} C_{n,r}(\{\mu_i\}) \frac{(n-k)!}{n!} \frac{\mu_1!}{(\mu_1 - k)!} \\ \times \langle Z(\beta\mu_1) \rangle \langle Z(\beta\mu_2) \rangle \cdots \langle Z(\beta\mu_r) \rangle, \quad (97) $$

where we take $\frac{1}{(\mu_1 - k)!} = 0$ when $\mu_1 < k$.

## ORCID iDs

Bernard Derrida https://orcid.org/0000-0001-6994-0226

Peter Mottishaw https://orcid.org/0000-0002-0091-4094

## References

[1] Parisi G 1979 Infinite number of order parameters for spin-glasses Phys. Rev. Lett. 43 1754–6

[2] Parisi G 1983 Order parameter for spin-glasses Phys. Rev. Lett. 50 1946–8

[3] Mézard M, Parisi G, Sourlas N, Toulouse G and Virasoro M 1984 Nature of the spin-glass phase Phys. Rev. Lett. 52 1156–9

[4] Mézard M, Parisi G, Sourlas N, Toulouse G and Virasoro M 1984 Replica symmetry breaking and the nature of the spin glass phase J. Phys. 145 843–54

[5] Mézard M, Parisi G and Virasoro M 1987 Spin Glass Theory and beyond: An Introduction to the Replica Method and its Applications vol 9 (Singapore: World Scientific)

[6] Mézard M and Montanari A 2009 Information, Physics, and Computation (Oxford: Oxford University Press)

[7] Sherrington D and Kirkpatrick S 1975 Solvable model of a spin-glass Phys. Rev. Lett. 35 1792–6

[8] Mézard M, Parisi G and Virasoro M A 1985 Random free energies in spin glasses J. Phys. Lett. 46 217–22

[9] Derrida B and Toulouse G 1985 Sample to sample fluctuations in the random energy model J. Phys. Lett. 46 223–8

[10] Rizzo T 2009 Chaos in mean-field spin-glass models Spin Glasses: Statics and Dynamics ed A B de Monvel and A Bovier (Basel: Birkhäuser) pp 143–57
[11] Franz S and Parisi G 1995 Recipes for metastable states in spin glasses *J. Phys. I* **5** 1401–15

[12] Rizzo T 2002 Ultrametricity between states at different temperatures in spin-glasses *Eur. Phys. J.* **B** **29** 425–35

[13] Derrida B 1980 Random-energy model: limit of a family of disordered models *Phys. Rev. Lett.* **45** 79–82

[14] Derrida B 1981 Random-energy model: an exactly solvable model of disordered systems *Phys. Rev.* **B** **24** 2613–26

[15] Franz S and Ney-Nifle M 1995 On chaos in mean-field spin glasses *J. Phys. A: Math. Gen.* **28** 2499–513

[16] Sales M and Bouchaud J-P 2001 Rejuvenation in the random energy model *Europhys. Lett.* **56** 181

[17] Gross D J and Mézard M 1984 The simplest spin glass *Nucl. Phys.* **B** **240** 431–52

[18] Ruelle D 1987 A mathematical reformulation of Derrida's REM and GREM *Commun. Math. Phys.* **108** 225–39

[19] Derrida B 1997 From random walks to spin glasses *Physica D* **107** 186–98

[20] Ghirlanda S and Guerra F 1998 General properties of overlap probability distributions in disordered spin systems. Towards Parisi ultrametricity *J. Phys. A: Math. Gen.* **31** 9149–55

[21] Bovier A and Kurkova I 2003 Rigorous results on some simple spin glass models *Markov Proc. Related Fields* **9** 209–42

[22] Bovier A 2006 *Statistical Mechanics of Disordered Systems* (Cambridge: Cambridge University Press)

[23] Panchenko D 2016 Chaos in temperature in generic 2p-spin models *Commun. Math. Phys.* **346** 703–39

[24] Derrida B and Mottishaw P 2015 Finite size corrections in the random energy model and the replica approach *J. Stat. Mech.* **P01021**

[25] Nieuwenhuizen T M 1996 A puzzle on fluctuations of weights in spin glasses *J. Phys. I* **6** 109–17

[26] Ferrero M E, Parisi G and Ranieri P 1996 Fluctuations in a spin-glass model with one replica symmetry breaking *J. Phys. A: Math. Gen.* **29** L569–74

[27] Campellone M, Parisi G and Virasoro M A 2009 Replica method and finite volume corrections *J. Stat. Phys.* **138** 29–39

[28] Rizzo T 2001 Against chaos in temperature in mean-field spin-glass models *J. Phys. A: Math. Gen.* **34** 5531

[29] Yoshino H and Rizzo T 2008 Stepwise responses in mesoscopic glassy systems: a mean-field approach *Phys. Rev. B* **77** 104429

[30] Parisi G 1980 A sequence of approximated solutions to the S-K model for spin glasses *J. Phys. A: Math. Gen.* **13** L115

[31] Parisi G 1980 The order parameter for spin glasses: a function on the interval 0-1 *J. Phys. A: Math. Gen.* **13** 1101–12

[32] Panchenko D 2013 *The Sherrington-Kirkpatrick Model* (New York: Springer)

[33] Guerra F 2014 Interpolation and comparison methods in the mean field spin glass model *Trends in Contemporary Mathematics* (Berlin: Springer) pp 1–12

[34] Dominicis C D and Young A P 1983 Weighted averages and order parameters for the infinite range Ising spin glass *J. Phys. A: Math. Gen.* **16** 2063

[35] Derrida B 1985 A generalization of the random energy model which includes correlations between energies *J. Phys. Lett.* **46** 401–7

[36] Derrida B and Spohn H 1988 Polymers on disordered trees, spin glasses, and traveling waves *J. Stat. Phys.* **51** 817–40

[37] Pain M and Zindy O 2018 Two-temperatures overlap distribution for the 2D discrete Gaussian free field (arXiv:1811.08313 [math-ph])

[38] Crisanti A and Sommers H-J 1992 The spherical p-spin interaction spin glass model: the statics *Z. Phys. B Condens. Matter* **87** 341–54