# The discrete random energy model and one step replica symmetry breaking

Bernard Derrida, Peter Mottishaw

► To cite this version:

Bernard Derrida, Peter Mottishaw. The discrete random energy model and one step replica symmetry breaking. Journal of Physics A: Mathematical and Theoretical, 2022, 55 (26), pp.265002. 10.1088/1751-8121/ac71e9. hal-03859136

HAL Id: hal-03859136

https://hal.science/hal-03859136v1

Submitted on 13 Jun 2025

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.

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

Distributed under a Creative Commons Attribution 4.0 International License
PAPER • OPEN ACCESS

# The discrete random energy model and one step replica symmetry breaking

To cite this article: Bernard Derrida and Peter Mottishaw 2022 *J. Phys. A: Math. Theor.* **55** 265002

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

You may also like

- *Ultrametric identities in glassy models of natural evolution*
Elena Agliari, Francesco Alemanno, Miriam Aquaro et al.

- *Finite size corrections in the random energy model and the replica approach*
Bernard Derrida and Peter Mottishaw

- *Hessian spectrum at the global minimum of high-dimensional random landscapes*
Yan V Fyodorov and Pierre Le Doussal
# The discrete random energy model and one step replica symmetry breaking

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 28 February 2022, revised 6 May 2022

Accepted for publication 20 May 2022

Published 7 June 2022

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

CrossMark

## Abstract

We solve the random energy model when the energies of the configurations take only integer values. In the thermodynamic limit, the average overlaps remain size dependent and oscillate as the system size increases. While the extensive part of the free energy can still be obtained by a standard replica calculation with one step replica symmetry breaking, it is no longer possible to recover the overlaps in this way. A possible way to adapt the replica approach is to allow the sizes of the blocks in the Parisi matrix to fluctuate and to take complex values.

Keywords: spin glass, replica method, random energy model

## 1. Introduction

Replica symmetry breaking (RSB) invented by Parisi in 1979 to solve spin glass models in their mean field version [1–5] has been used since in a large number of contexts including the theory of neural networks, problems of optimization, directed polymers and glasses [6, 7].

One of the first achievements of Parisi’s theory was to predict the exact expression of the extensive part of the free energy of the mean field model of spin glasses introduced by Sherrington and Kirkpatrick (the SK model) [8]. After 20 years of effort, this exact expression for the free energy was confirmed by rigorous mathematical studies [9–13].

The other striking predictions associated with RSB concern the overlaps, which represent a way of quantifying the fluctuations of the free energy landscapes [14, 15]. Several of the

<figure><img src="image_3.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>
predictions of Parisi's theory, such as ultrametricity or the Ghirlanda-Guerra relations were also mathematically confirmed [16–18] but some of these proofs are dependent on the Gaussian character of the interactions.

In the present note, we will report on a simple case, the discrete random energy model (DREM), for which the standard one step RSB does predict the correct (extensive part) of the free energy but fails to give the correct expression for the overlaps. In contrast to the original random energy model, the energies in the DREM are no longer Gaussian but instead take only integer values. The DREM has already been considered to locate the complex zeroes of the partition function [19] and to investigate phase transitions of the nth moments of the partition function as n varies [20] (see also [21] for a rigorous derivation of the extensive part of the free energy).

For the original REM as well as for the DREM, the system consists of an exponential number $2^N$ of configurations $C$ and the energies $E(C)$ of these configurations are independent and identically distributed (i.i.d.) variables with a distribution which scales with the system size $N$. The goal then is to determine the average free energy, $\langle F \rangle = -\frac{1}{\beta}\langle \log Z \rangle$, where

$$ \langle \log Z(\beta) \rangle = \left\langle \log \left( \sum_{C=1}^{2^N} e^{-\beta E(C)} \right) \right\rangle \qquad (1) $$

(where $\langle \cdot, \cdot \rangle$ denotes an average over disorder, i.e. over the $2^N$ energies $E(C)$) and the average overlaps $\langle Y_k \rangle$ which are defined by

$$ \langle Y_k(\beta) \rangle = \left\langle \frac{\sum_C e^{-\beta k E(C)}}{\left(\sum_C e^{-\beta E(C)}\right)^k} \right\rangle = \left\langle \frac{Z(k\beta)}{Z(\beta)^k} \right\rangle. \qquad (2) $$

In the case of random energy models, $Y_k$ is simply the probability of finding $k$ copies of the system in the same configuration [15, 22, 23], (since the probability of occupying configuration $C$ is

$$ p(C) = \frac{e^{-\beta E(C)}}{\sum_C e^{-\beta E(C)}} $$

the probability of finding $k$ copies of the system in the same configuration is $\sum_C p(C)^k$ which is precisely $Y_k$.

In the original REM [7, 24–27]

$$ P(E) = \frac{1}{\sqrt{\pi N}} e^{-\frac{E^2}{N}}. \qquad (3) $$

For this choice it is known that, in the low temperature phase (the frozen phase) $\beta > \beta_c = 2\sqrt{\log 2}$ the average free energy (1) and the overlaps (2) are given by, in the large $N$ limit,

$$ \langle \log Z(\beta) \rangle = N\beta\sqrt{\log 2} - \frac{\beta}{2\beta_c}\log(4N\pi\log 2) \\ + \frac{\beta}{\beta_c}\log\left(\Gamma\left(1 - \frac{\beta_c}{\beta}\right)\right) + \left(1 - \frac{\beta}{\beta_c}\right)\Gamma'(1) \qquad (4) $$

$$ \langle Y_k(\beta) \rangle = \frac{\Gamma\left(k - \frac{\beta_c}{\beta}\right)}{\Gamma\left(1 - \frac{\beta_c}{\beta}\right) (k-1)!} \qquad (5) $$
up to vanishingly small corrections in the large N limit [22, 23, 28–30].

Our goal in the present work is to see how these expressions are modified when the energies $E(C)$ are discrete, for example when

$$P(E = \nu) = \binom{N}{\nu} (1-p)^{\nu} p^{N-\nu}, \quad \nu = 0, 1 \dots N. \qquad (6)$$

In (6) each energy can be viewed as the sum of N independent binary random variables which take value 0 with probability $p$ and 1 with probability $1-p$ (while in (3) each energy is a sum of N independent Gaussians).

Our main result will be to show that, although the replica approach still predicts the correct extensive part of the free energy, the expression (5) of the overlaps is no longer valid. In section 3, we will show, by an exact calculation, that the main difference between the Gaussian (3) and the discrete energy (6) cases is that the $O(1)$ finite size corrections of free energy as well as the overlaps are no longer constant as in (4) and (5) but oscillate as the system size increases, even in the large N limit. As discussed in section 4, this raises the question of how to adapt the RSB calculation to reproduce these periodic dependences.

## 2. How to calculate the free energy and the overlaps

In this section we establish the following three formulas, valid for both the original REM and the DREM (in their Poisson process versions), which give the average free energy, the negative moments of the partition function and the average overlaps in terms of the generating function $\langle e^{-tZ} \rangle$ of the partition function:

$$\langle \log Z(\beta) \rangle = \int_0^\infty dt \frac{e^{-t} - \langle e^{-tZ} \rangle}{t} = \int_0^\infty dt \log t \left( e^{-t} + \frac{d e^{\phi(t)}}{dt} \right) \qquad (7)$$

$$\langle Z(\beta)^n \rangle = \frac{1}{\Gamma(-n)} \int_0^\infty dt \, t^{-n-1} e^{\phi(t)} = \frac{1}{\Gamma(1-n)} \int_0^\infty dt \, t^{-n} \frac{d e^{\phi(t)}}{dt} \quad \text{for } n < 0 \qquad (8)$$

$$\langle Y_k(\beta) \rangle = \frac{1}{\Gamma(k)} \int_0^\infty dt \, t^{k-1} e^{\phi(t)} \left[ (-)^k \frac{d^k \phi(t)}{dt^k} \right], \qquad (9)$$

where $\phi(t)$ is defined by

$$e^{\phi(t)} = \langle e^{-tZ} \rangle. \qquad (10)$$

These expressions will be used in section 3. (For $n > 0$ alternative expressions could easily be obtained from (8) by successive integrations by parts. We will not use them below. Also for positive $n$, as $n$ varies, the moments $\langle Z^n \rangle$ exhibit phase transitions [20, 31, 32]. Note that the replica method works just as well taking $n \to 0$ from above or below.)

The first two identities (7) and (8) are obviously valid for any positive random variable $Z$. Let us now establish the third one (9) for the DREM, i.e. in the case of where the energies $E$ can only take integer values as in (6). The derivation can easily be extended to the case of continuous energies.

It is well known that, in the low temperature phase (and in the large $N$ limit), only configurations with energies close enough to the ground state (i.e. such that $E(C) - E_{\text{ground state}} \ll N$)
matter [30]. This implies that, up to vanishingly small corrections [33], one can replace a REM of $2^N$ energy levels distributed according to a distribution $P(E)$ by a Poisson REM with energies generated by a Poisson process of density

$$\rho(E) = 2^N P(E).$$

In the case of the DREM, where the energies take discrete values $\nu$, several configurations $C$ may have the same energy $E(C)$ and the partition function is of the form

$$Z = \sum_{\nu} m_{\nu} e^{-\beta \nu}, \qquad (11)$$

where $m_\nu$ is the number of configurations $C$ such that $E(C) = \nu$. Then since the density $\rho(E)$ is of the form

$$\rho(E) = \sum_{\nu=0}^{N} r_{\nu} \delta(E - \nu) \quad \text{with} \quad r_{\nu} = 2^{N} P(E = \nu)$$

the numbers $m_\nu$ of configurations with energy $\nu$ are randomly distributed according to Poisson distributions

$$\text{Pro}(m_{\nu} = m) = \frac{(r_{\nu})^{m}}{m!} e^{-r_{\nu}}. \qquad (12)$$

Moreover these degeneracies $m_\nu$ are independent.

Then the generating function (see (10)) of the partition function $Z$ is given by

$$\langle e^{-tZ} \rangle = \prod_{\nu=-\infty}^{\infty} \left( \sum_{m=0}^{\infty} \frac{r_{\nu}^{m}}{m!} e^{-r_{\nu}} \exp[-mt e^{-\beta\nu}] \right) = e^{\phi(t)} \qquad (13)$$

with $\phi(t)$ given by

$$\phi(t) = \sum_{\nu} r_{\nu} (\exp[-t e^{-\beta \nu}] - 1). \qquad (14)$$

To obtain the average overlaps $\langle Y_k \rangle$ defined in (2) one first starts by computing the following average over disorder (i.e. over the $m_\nu$'s)

$$\langle Z(k\beta) e^{-tZ(\beta)} \rangle = \left\langle \sum_{\nu} m_{\nu} e^{-k\beta\nu} \exp \left[ -\sum_{\nu'} m_{\nu'} t e^{-\beta\nu'} \right] \right\rangle. \qquad (15)$$

Following the same procedure as for the derivation of (14) from (13) one gets

$$\langle Z(k\beta) e^{-tZ(\beta)} \rangle = (-)^k e^{\phi(t)} \frac{d^k \phi(t)}{dt^k}. \qquad (16)$$

Then rewriting (2) as

$$\langle Y_k(\beta) \rangle = \frac{1}{(k-1)!} \int_0^\infty dt \, t^{k-1} \langle Z(k\beta) e^{-tZ(\beta)} \rangle$$

leads to (9).
The validity of (9) can easily be checked in the case of the original REM, i.e. when the distribution of energies is continuous. The only difference is that (14) becomes

$$ \phi(t) = \log\langle e^{-tZ} \rangle = \int dE \rho(E) (\exp[-t e^{-\beta E}] - 1) \quad (17) $$

i.e. that the sum in (14) is replaced by an integral.

### 3. Continuous versus discrete energies

In the rest of this paper we will consider that each energy $E(C)$ is the sum of $N$ i.i.d. random variables as in (3) or (6). Therefore the generating function of each of these energies is of the form

$$ \langle e^{-\beta E(C)} \rangle = e^{N\Psi(\beta)}. \quad (18) $$

For example for the distributions of energies (3) and (6) one has respectively

$$ \Psi(\beta) = \frac{\beta^2}{4} \quad (19) $$

$$ \beta) = \log(p + (1 - p) e^{-\beta}). \quad (20) $$

From (18), writing that $\int dE P(E)e^{-\beta E} = e^{N\Psi(\beta)}$ one can see that for large $N$

$$ P(E) \approx \sqrt{\frac{G''(\frac{E}{N})}{2N\pi}} e^{-NG(\frac{E}{N})}, \quad (21) $$

where the function $G(\epsilon)$ can be expressed in a parametric form in terms of $\Psi(\beta)$

$$ \epsilon = -\Psi'(\beta); \quad G(\epsilon) = \beta\Psi'(\beta) - \Psi(\beta). \quad (22) $$

Since $G$ and $\Psi$ are Legendre transforms of each other, one also has

$$ \beta = -G'(\epsilon); \quad \Psi(\beta) = \epsilon G'(\epsilon) - G(\epsilon); \quad G''(\epsilon) = \frac{1}{\Psi''(\beta)}. \quad (23) $$

Here we are assuming that both $G(\epsilon)$ and $\Psi(\beta)$ are $C^2$, i.e. they have at least two continuous derivatives (which they have in the cases (3) and (6)).

Because the $2^N$ energy levels are i.i.d., one has that, close to the ground state energy, $2^N P(E)$ is neither exponentially large in $N$ nor exponentially small in $N$. Therefore if $\epsilon_c$ is the minimal solution of

$$ G(\epsilon_c) = \log 2, \quad (24) $$

the extensive part of the ground state energy is $N\epsilon_c$.

In the REM the low temperature phase is dominated by the configurations which have an energy at distance of order 1 from the ground state energy. Thus one can replace the distribution $P(E)$ by an exponential distribution [30] which approximates the true $P(E)$ near $E = N\epsilon_c$. Therefore

$$ \rho(E) = 2^N P(E) \approx A e^{\beta_c(E - N\epsilon_c)} \quad \text{where } A = \sqrt{\frac{G''(\epsilon_c)}{2\pi N}} = \sqrt{\frac{1}{2\pi N \Psi''(\beta_c)}} \quad (25) $$
and $\beta_c$ is defined by

$$ \beta_c = -G'(\epsilon_c). \qquad (26) $$

Note that, from (22)–(24) and (26),

$$ \beta_c \Psi'(\beta_c) - \Psi(\beta_c) = \log 2. \qquad (27) $$

This relation will appear in section 4 when we will look at the replica approach.

Let us now discuss the cases of continuous energies and discrete energies separately:

### 3.1. Continuous energies as in (3)

For continuous energies, the integral in (17) can be performed exactly for $\rho(E)$ given by (25) and one gets

$$ \phi(t) = - \frac{A e^{-N\beta_c \epsilon_c}}{\beta_c} t^{\frac{\beta_c}{\beta}} \Gamma\left(1 - \frac{\beta_c}{\beta}\right) \qquad (28) $$

and one gets for the average free energy (7) and for the negative moments of the partition function (8)

$$ \langle \log Z \rangle = -N\beta\epsilon_c + \frac{\beta}{\beta_c} \log \left( \frac{A}{\beta_c} \right) + \frac{\beta}{\beta_c} \log \left( \Gamma \left( 1 - \frac{\beta_c}{\beta} \right) \right) + \left( 1 - \frac{\beta}{\beta_c} \right) \Gamma'(1) + o(1) \qquad (29) $$

$$ \langle Z^n \rangle \approx \frac{\Gamma\left(1 - n\frac{\beta}{\beta_c}\right)}{\Gamma(1 - n)} \left[ \frac{A}{\beta_c} \Gamma\left(1 - \frac{\beta_c}{\beta}\right) \right]^{\frac{n\beta}{\beta_c}} e^{-Nn\beta \epsilon_c}. \qquad (30) $$

It is easy to check using (25) that for $G(\epsilon) = -\epsilon^2$, one has $\epsilon_c = -\sqrt{\log 2}$, $\beta_c = 2\sqrt{\log 2}$ and that (29) reduces to (4). (Also, not surprisingly one can recover (29) from (30) in the limit $n \to 0^-$.) Therefore (29) gives the generalization of (4) to other continuous distributions of energies $G(\epsilon)$.

For the overlaps, using (28) into (9) leads to

$$ \langle Y_k(\beta) \rangle = \frac{\Gamma\left(k - \frac{\beta_c}{\beta}\right)}{\Gamma\left(1 - \frac{\beta_c}{\beta}\right) (k - 1)!} \qquad (31) $$

which is identical to (5).

It is remarkable, though expected from the fact that $\rho(E)$ can be replaced by an exponential (25) that the average overlaps, in the large $N$ limit, do not depend on the details of the distribution $P(E)$.

### 3.2. Discrete energies as in (6)

When the energies take only integer values, the parameters $r_\nu$ in (12) are given by

$$ r_\nu = 2^N P(E = \nu) = A e^{-N\beta_c \epsilon_c} e^{\beta_c \nu}. $$
<figure><img src="image_4.png" /><figcaption>Figure 1. $w(t)$ versus $(\log t)/\beta$ for $\beta_c = \beta/2$ and $\beta = 2.5, 5, 7.5, \dots, 17.5$. The amplitude of the oscillations increase with $\beta$.</figcaption></figure>

Then $\phi(t)$ in (14) can be written as

$$ \phi(t) = -\frac{A e^{-N\beta_c\epsilon_c}}{\beta_c} t^{\frac{\beta_c}{\beta}} w(t), \quad (32) $$

where

$$ w(t) = \beta_c \sum_{\nu=-\infty}^{\infty} t^{-\frac{\beta_c}{\beta}} e^{\beta_c \nu} (1 - \exp[-t e^{-\beta_c \nu}]). \quad (33) $$

A priori the sum in (33) should run from 0, ..., N. However, in the low temperature phase, since the values of $\nu$ that contribute to the sum are those close to the ground state value $\nu = N\epsilon_c$, one can extend the sum to $-\infty$ and $+\infty$ with errors that are exponentially small in N.

Using the Poisson summation formula

$$ \sum_n f(n) = \sum_q \int dz f(z) e^{i2\pi qz} $$

one can also write $w(t)$ as

$$ w(t) = \sum_{q=-\infty}^{\infty} w_q t^{\frac{2i\pi q}{\beta}} \quad \text{where} \quad w_q = -\frac{\beta_c}{\beta} \Gamma\left(-\frac{\beta_c + 2i\pi q}{\beta}\right). \quad (34) $$

The main difference with the case of continuous energies (28) is that the function $w(t)$ in (32) is no longer constant. It is in fact easy to check in (33) that

$$ w(t) = w(e^{\beta} t) \quad (35) $$

i.e. is $w(t)$ a periodic function of $\log t/\beta$ (see figure 1).

Starting from (7)–(9) and making the change of variable

$$ t \rightarrow \left(\frac{A}{\beta_c}\right)^{-\frac{\beta}{\beta_c}} e^{N\beta_c\epsilon_c} \tau $$
one gets for the free energy

$$ \langle \log Z \rangle = -N\beta\epsilon_c + \frac{\beta}{\beta_c} \log \left( \frac{A}{\beta_c} \right) + \Gamma'(1) \\ + \int_0^\infty d\tau \log(\tau) \frac{d}{d\tau} \left[ \exp \left( -\tau^{\frac{\beta_c}{\beta}} w(B\tau) \right) \right] + o(1) \quad (36) $$

for the negative moments of Z

$$ \langle Z^n \rangle \approx \frac{1}{\Gamma(1-n)} e^{-n\beta\epsilon_c N} \left(\frac{A}{\beta_c}\right)^{\frac{\beta}{\beta_c}n} \int_0^\infty d\tau \tau^{-n} \frac{d}{d\tau} \left[ \exp \left[ -\tau^{\frac{\beta_c}{\beta}} w(B\tau) \right] \right] \quad (37) $$

and for the overlaps

$$ \langle Y_k(\beta) \rangle = \frac{(-1)^{k-1}}{(k-1)!} \int_0^\infty d\tau \tau^{k-1} \exp \left[ -\tau^{\frac{\beta_c}{\beta}} w(B\tau) \right] \frac{d^k \left( \tau^{\frac{\beta_c}{\beta}} w(B\tau) \right)}{d\tau^k} \quad (38) $$

where

$$ B = \left( \frac{\beta_c}{A} \right)^{\frac{\beta}{\beta_c}} e^{N\beta\epsilon_c} \quad (39) $$

with A given in (25). These expressions (36)–(39) are to be compared to (29)–(31). Because (see (35)) the function $w(t)$ is a periodic function of $\frac{\log t}{\beta}$, the overlaps (38) as well as the contributions of order 1 (i.e. the last terms in (36) and the last factor in (37)) are periodic functions of $(\log B)/\beta$.

The fact that the overlaps remain size dependent through the parameter B in (39) as the system size increases can be illustrated by the fact that in the plane $\langle Y_2 \rangle, \langle Y_3 \rangle$, expression (38) leads to closed loops instead of a single point for a fixed choice of $\beta_c$ and $\beta$ (see figure 2).

# 4. The replica approach

In this section we try to see how the exact expressions (36)–(38) in the discrete energy case and the corresponding expressions (29)–(31) in the continuous energy case could be recovered by Parisi’s replica approach.

## 4.1. The moments of the partition function

In the replica calculation one usually starts by the calculation of the positive integer moments of the partition function. To do so, one has (see (14) or (17))

$$ \phi(t) = \sum_{\mu \ge 1} \frac{(-t)^{\mu}}{\mu!} \langle Z(\beta\mu) \rangle \quad (40) $$

and from (10)

$$ \langle e^{-tZ} \rangle = \exp \left[ \sum_{\mu \ge 1} \frac{(-t)^{\mu}}{\mu!} \langle Z(\beta\mu) \rangle \right] $$
<figure><img src="image_5.png" /><figcaption>**Figure 2.** $\langle Y_3 \rangle$ versus $\langle Y_2 \rangle$ given by (38) when the system size increases at a fixed value of the ratio $\frac{\beta_c}{\beta} = \frac{1}{2}$. The various ellipses correspond to several choices of $\beta = 4, 5, \dots, 10$. Similar curves could be obtained by plotting quantities such as $\langle Y_4 \rangle$ vs $\langle Y_2 \rangle$ or $\langle Y_2^2 \rangle$ vs $\langle Y_2 \rangle$. As the system size increases, $\langle Y_2 \rangle$ and $\langle Y_3 \rangle$ vary when the energies are discrete (38). On the contrary for continuous energies, the plots would reduce to a single point: $\langle Y_2 \rangle = \frac{1}{2}, \langle Y_3 \rangle = \frac{3}{8}$ (see(5)) independent of $\beta$ and of the system size (in the limit of large system sizes).</figcaption></figure>

so that [25, 34]

$$ \begin{aligned} \langle Z^n(\beta) \rangle &= \sum_{r \ge 1} \frac{n!}{r!} \sum_{\mu_1 \ge 1} \frac{\langle Z(\beta \mu_1) \rangle}{\mu_1!} \cdots \sum_{\mu_r \ge 1} \frac{\langle Z(\beta \mu_r) \rangle}{\mu_r!} \delta_{\mu_1 + \cdots + \mu_r - n} \\ &= \sum_{r \ge 1} \frac{n!}{r!} 2^{Nr} \sum_{\mu_1 \ge 1} \cdots \sum_{\mu_r \ge 1} \frac{e^{N[\Psi(\beta\mu_1) + \cdots + \Psi(\beta\mu_r)]}}{\mu_1! \times \cdots \times \mu_r!} \delta_{\mu_1 + \cdots + \mu_r - n}. \end{aligned} \quad (41) $$

Each term represents a partition of the $n$ replicas into $r$ blocks of sizes $\mu_1, \dots, \mu_r$.

## 4.2. The overlaps

To obtain the overlap one has from (40)

$$ \frac{d^k \phi(t)}{dt^k} = (-)^k \sum_{n=0}^{\infty} \frac{(-t)^n}{n!} \langle Z((n+k)\beta) \rangle $$

and from (10)

$$ e^{\phi(t)} = \sum_{n=0}^{\infty} \frac{(-t)^n}{n!} \langle Z(\beta)^n \rangle. $$

Then as in (16)

$$ \langle Z(k\beta)e^{-tZ(\beta)} \rangle = (-)^k \left( \frac{d^k \phi(t)}{dt^k} \right) e^{\phi(t)} \quad (42) $$
one gets by expanding in powers of t

$$ \begin{aligned} \langle Z(k\beta) Z(\beta)^{n-k} \rangle &= (n-k)! \sum_{\mu=k}^{n} \frac{\langle Z(\mu\beta) \rangle \langle Z(\beta)^{n-\mu} \rangle}{(\mu-k)!(n-\mu)!} \\ &= (n-k)! \sum_{\mu=0}^{\infty} \frac{\langle Z(\mu\beta) \rangle \langle Z(\beta)^{n-\mu} \rangle}{(\mu-k)!(n-\mu)!}, \end{aligned} \quad (43) $$

where in the last line we simply have extended the sum to the range 0 to infinity, having in mind that factorials of negative integers are infinite.

This expression is also valid in the special case $k = 1$ so that one can write

$$ \frac{\langle Z(k\beta) Z^{n-k}(\beta) \rangle}{\langle Z^n(\beta) \rangle} = \sum_{\mu} W(\mu) \frac{(n-k)!(\mu-1)!}{(n-1)!(\mu-k)!}, \quad (44) $$

where the weights $W(\mu)$ are defined as

$$ W(\mu) = \frac{\langle Z(\mu\beta) \rangle \langle Z(\beta)^{n-\mu} \rangle}{(\mu-1)!(n-\mu)!} \left( \sum_{\mu'} \frac{\langle Z(\mu'\beta) \rangle \langle Z(\beta)^{n-\mu'} \rangle}{(\mu'-1)!(n-\mu')!} \right)^{-1} $$

These expressions will allow us to obtain the overlaps (2) using,

$$ \langle Y_k(\beta) \rangle = \lim_{n \to 0} \frac{\langle Z(k\beta) Z^{n-k}(\beta) \rangle}{\langle Z^n(\beta) \rangle} $$

but as usual in the replica calculation, the question is to analytically continue expressions like (41) or (44) to non-integer (positive or negative) values of $n$. To do so we will need to allow the sizes of the blocks $\mu$ in (41) and (44) to take non integer values.

### 4.3. The continuous energies case

Following Parisi's original approach [2] we assume that, for non-integer $n$, in the limit of a large system size ($N \to \infty$), the multiple sum in (41) is dominated by a single term, with all the $\mu_i$ equal to a single value $\mu$ and $r = \frac{n}{\mu}$ so that all blocks have the same size $\mu$. As he also did, we let this common value $\mu$ be a real number (instead of an integer) and we also assume that the contribution to the sum in (41) is minimum (in contrast with standard saddle point methods where one would need to find the maximum) at this value of $\mu$. In other words, for large $N$ we replace (41) by

$$ \langle Z^n \rangle = \mathcal{P} \times \min_{\mu} \left[ \frac{n!}{\left(\frac{n}{\mu}\right)! (\mu!)^{\frac{n}{\mu}}} \exp \left[ N \frac{n}{\mu} (\log 2 + \Psi(\beta\mu)) \right] \right] \quad (45) $$

where the prefactor $\mathcal{P}$ represents the 'putative' contributions to the multiple sums due to the neighbourhood of the extremum where all the $\mu_i$ are equal to $\mu$. The difficulty in determining the prefactor comes from the fact that once we let $r$ and the $\mu_i$ take real instead of integer values in (41), it is not clear how to estimate the contribution of the fluctuations in the neighbourhood of the extremum.

The extremum over $\mu$ is given by the solution $\mu_0$ of

$$ \Psi(\beta\mu) + \log 2 - \beta\mu\Psi'(\beta\mu) = 0. \quad (46) $$
Then from (27), one can see that

$$ \mu_0 = \frac{\beta_c}{\beta}. \qquad (47) $$

On the other hand, using (47), the expression of A (25) and the fact that $\epsilon_c = -\Psi'(\beta_c) = -(\log 2 + \Psi(\beta_c))/\beta_c$ (see (22) and (24)) the exact expression (30) can be rewritten as

$$ \langle Z^n \rangle \approx \left[ \frac{1}{2\pi N \Psi''(\beta\mu_0)\beta^2} \right]^{\frac{n}{2\mu_0}} \frac{\Gamma\left(1-\frac{n}{\mu_0}\right)}{\Gamma(1-n)} \left[ \frac{\Gamma(1-\mu_0)}{\mu_0} \right]^{\frac{n}{\mu_0}} \times \exp \left[ N \frac{n}{\mu_0} (\log 2 + \Psi(\beta\mu_0)) \right]. \qquad (48) $$

Clearly the exponential N dependence is the same in the exact (48) and the replica (45) expressions. Understanding the prefactor in (48) or the O(1) corrections in (29) using the replica approach is more problematic. There are however some similarities between the exact and the replica expressions. For example one could assume that factorials in (45) can be replaced by gamma functions in (48): $m! \rightarrow \frac{1}{\Gamma(1-m)}$. Then one would need to argue that the remaining factors come from Gaussian integrals near the saddle point so that

$$ \mathcal{P} = \left[ \frac{1}{2\pi N \Psi''(\beta\mu_0)\mu_0^2\beta^2} \right]^{\frac{n}{2\mu_0}}. \qquad (49) $$

In fact there is a way to recover this prefactor by using an integral representation [33, 35], but this is not strictly speaking a replica calculation because from the very start one tries to compute non-integer moments of the partition function. We will not discuss it here.

To obtain the average overlaps in the replica approach one can rewrite (44) as

$$ \langle Y_k(\beta) \rangle = \frac{\langle Z(k\beta) Z^{n-k}(\beta) \rangle}{\langle Z^n(\beta) \rangle} = \sum_{\mu} W(\mu) \frac{\Gamma(1-n)\Gamma(k-\mu)}{\Gamma(k-n)\Gamma(1-\mu)}, \qquad (50) $$

where we have used the fact that a generalization of the combinatorial factors in (44) to non-integer values of n and μ can be written in terms of gamma functions:

$$ \frac{(n-k)!(\mu-1)!}{(n-1)!(\mu-k)!} = \prod_{p=1}^{k-1} \left( \frac{p-\mu}{p-n} \right) = \frac{\Gamma(1-n)\Gamma(k-\mu)}{\Gamma(k-n)\Gamma(1-\mu)}. $$

Comparing the exact result (31) and (2) in the continuous energy case

$$ \langle Y_k(\beta) \rangle = \frac{\Gamma\left(k - \frac{\beta_c}{\beta}\right)}{\Gamma\left(1 - \frac{\beta_c}{\beta}\right) (k-1)!} $$

with the replica expression (50), in the $n \to 0$ limit, we see that they agree provided that we choose

$$ W(\mu) = \delta\left(\mu - \frac{\beta_c}{\beta}\right). $$

In other words the size μ of the blocks is fixed and takes the value $\mu_0$ already obtained in the calculation of the non-integer moments of Z (see (47)).
<figure><img src="image_6.png" /><figcaption>Figure 3. We plot $U_{N+1}$ versus $U_N$, where $U_N = E_N - E_{N-1} - \frac{\gamma}{N\beta_c}$ for the DREM on the left and for the directed polymer on the tree on the right. The parameters chosen are $p = 0.1$, $\epsilon_c = 0.42251$, $\beta_c = 2.5097$ and $500 < N < 1500$. We see clearly the oscillations. In absence of oscillations, the two graphs would reduce to a single value at $\epsilon_c$. One can notice that for the tree the amplitude of the oscillations is much smaller than for the DREM.</figcaption></figure>

### 4.4. The discrete energies case

For discrete energies the simplest procedure would be to repeat the replica analysis we did above in section 4.3 (i.e. consider that the sum in (41) is dominated by a single term where all the $\mu_i$'s are equal). This would lead to the expression (45) allowing to recover the right N-dependence in the exponential of (37) and therefore the right extensive part of the free energy. However considering that the sum over $\mu$ in the replica expression for the overlaps (50) reduces to a single term would lead to the same result as the continuous case (9) and (31) in contradiction with the exact expression (38).

On the other hand, as $\Psi(\beta)$ in (20) is a periodic function ($\Psi(\beta) = \Psi(\beta + 2i\pi)$), it is clear that if $\mu_0 = \frac{\beta_c}{\beta}$ is solution of (46) then all

$$ \mu_q = \frac{\beta_c}{\beta} + \frac{2i\pi q}{\beta} \quad \text{with } q \in \mathbb{Z} \qquad (51) $$

are also saddle points at the same height as $\mu_0$. Therefore in (41), one could keep the same exponential N-dependence by allowing the $\mu_i$'s to take all possible values $\mu_q$. This gives the correct extensive part of the free energy but, as in the continuous energy case, we do not have a clear way of understanding the prefactor in (48) using replicas.

For the overlap, however, we are going to see now that one can reproduce the exact expressions (38) if one allows the values of $\mu$ in (50) to fluctuate and to take complex values. From (34) we have

$$ \frac{d^k \left( \tau^{\frac{\beta_c}{\beta}} w(B\tau) \right)}{d\tau^k} = (-)^{k+1} \sum_{q=-\infty}^{\infty} \tau^{\frac{\beta_c + 2i\pi q}{\beta} - k} B^{\frac{2i\pi q}{\beta}} \Gamma \left( k - \frac{\beta_c + 2i\pi q}{\beta} \right). $$
Then inserting this expression into (38) and using (51) one gets the exact expression

$$ \langle Y_k(\beta) \rangle = \sum_{q=-\infty}^{\infty} W_q \frac{\Gamma(k - \mu_q)}{\Gamma(k) \Gamma(1 - \mu_q)}, \qquad (52) $$

where

$$ W_q = B^{\frac{2i\pi q}{\beta}} \frac{\beta_c}{\beta} \Gamma(1 - \mu_q) \int_0^{\infty} \tau^{\mu_q - 1} \exp \left[ -\tau^{\frac{\beta_c}{\beta}} w(B\tau) \right] d\tau. \qquad (53) $$

This is precisely the same form as the replica expression (50) (in the $n \to 0$ limit). Therefore (52) shows that one way to interpret the exact expression (38) within the replica approach is to allow the sizes of the blocks to fluctuate and to take complex values (with also complex weights $W_q$).

## 5. Conclusion

In the present paper we have obtained exact expressions for the overlaps and for the finite size corrections for the random energy model with discrete energies. Although the standard replica broken symmetry calculation does give the correct extensive part of the free energy, it fails to predict the overlaps. In fact, in the thermodynamic limit, the exact expression (38) of the overlaps oscillates as the system size $N$ increases. Trying to interpret these oscillations using replicas, we saw that one needs to allow the sizes of the blocks in the Parisi ansatz to fluctuate and to take complex values.

That $\langle Y_k \rangle \neq 1$ at zero temperature is not a surprise because it is expected that the ground state may be degenerate at zero temperature [36, 37]. What we have shown here is that even at finite temperature, the standard one step RSB does give the correct free energy but does not give the correct overlaps. That the oscillations appear both in overlaps and in the $O(1)$ corrections of the free energy makes a lot of sense because in the multi-valley picture of the RSB, a variation of $O(1)$ of the free energy of a valley would change the overlaps. It would be interesting to see whether a similar phenomenon could appear in more sophisticated models like the binary perceptron [38–42] or the K-sat problem [43, 44] for which the energy spectrum is discrete and the ground state is degenerate. One case at least for which the oscillations of the $O(1)$ in the free energy is visible (see figure 3) is the case of the directed polymer on a tree [45] where the energies $\epsilon_b$ of the bonds have a binary distribution $\epsilon_b = 1$ with probability $1-p$ and $\epsilon_b = 0$ with probability $p$. For the DREM as well for the tree the ground state energy $E_N$ is of the form

$$ E_N = N\epsilon_c + \frac{\gamma}{\beta_c} \log N + e_N + o(1) \qquad (54) $$

with $\gamma = \frac{1}{2}$ for the DREM and $\gamma = \frac{3}{2}$ for the tree [46–48]. According to (36), $e_N$ oscillates with $N$ for the DREM. What figure 3 shows is that it also oscillates in the case of the tree.

That the sizes of the blocks in the replica approach may fluctuate is a rather general phenomenon. It could even occur for continuous distributions of energies. In fact the expression (5) and (31) is a direct consequence of the fact that the energies close to the ground state energy have a distribution well approximated by a single exponential distribution (25). If for some reason, for example near a phase transition, this distribution would instead be well approximated by the sum of two exponentials

$$ \rho(E) \approx A e^{\beta_c(E - N\epsilon_c)} + A' e^{\beta'_c(E - N\epsilon_c)} $$
then the function $\phi(t)$ in (17) would become

$$\phi(t) = -\frac{A e^{-N\beta_c \epsilon_c}}{\beta_c} t^{\frac{\beta_c}{\beta}} \Gamma\left(1-\frac{\beta_c}{\beta}\right) - \frac{A' e^{-N\beta'_c \epsilon_c}}{\beta'_c} t^{\frac{\beta'_c}{\beta}} \Gamma\left(1-\frac{\beta'_c}{\beta}\right).$$

Then one could repeat the calculation of section 4.4 and obtain that

$$\langle Y_k(\beta) \rangle = W \frac{\Gamma(k-\mu)}{\Gamma(k) \Gamma(1-\mu)} + W' \frac{\Gamma(k-\mu')}{\Gamma(k) \Gamma(1-\mu')}, \quad (55)$$

with $\mu = \frac{\beta_c}{\beta}$, $\mu' = \frac{\beta'_c}{\beta}$ and the weights given by

$$W = \frac{A e^{-N\beta_c \epsilon_c}}{\beta} \Gamma\left(1 - \frac{\beta_c}{\beta}\right) \int_0^\infty t^{\frac{\beta_c}{\beta} - 1} e^{\phi(t)} dt;$$

$$W' = \frac{A' e^{-N\beta'_c \epsilon_c}}{\beta} \Gamma\left(1 - \frac{\beta'_c}{\beta}\right) \int_0^\infty t^{\frac{\beta'_c}{\beta} - 1} e^{\phi(t)} dt.$$

Clearly (55) shows that the sizes of the blocks in the replica calculation take the value $\mu$ with probability $W$ and $\mu'$ with probability $W'$ (it is easy to check that $W$ and $W'$ are positive numbers and that $W + W' = 1$ using the fact that $\phi(0) = 0$ and $\phi(\infty) = -\infty$).

As the discrete energy case can be viewed as a situation where the density of energies is a sum of exponentials with complex growth rates, it is not surprising that we found boxes with complex sizes in (52).

## Data availability statement

Any data that support the findings of this study are included within the article.

## 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 1980 A sequence of approximated solutions to the S-K model for spin glasses J. Phys. A: Math. Gen. 13 L115

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

[4] 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

[5] 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. France 45 843–54

[6] 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)

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

[8] Sherrington D and Kirkpatrick S 1975 Solvable model of a spin-glass Phys. Rev. Lett. 35 1792–6
[9] Guerra F 2003 Broken replica symmetry bounds in the mean field spin glass model *Commun. Math. Phys.* **233** 1–12

[10] Guerra F and Toninelli F L 2002 The thermodynamic limit in mean field spin glass models *Commun. Math. Phys.* **230** 71–9

[11] Guerra F 1996 About the overlap distribution in mean field spin glass models *Int. J. Mod. Phys. B* **10** 1675–84

[12] Talagrand M 2003 *Spin Glasses: A Challenge for Mathematicians: Cavity and Mean Field Models* vol 46 (Berlin: Springer)

[13] Talagrand M 2010 *Mean Field Models for Spin Glasses: Volume I: Basic Examples* vol 54 (Berlin: Springer)

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

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

[16] 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

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

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

[19] Moukarzel C and Parga N 1991 Numerical complex zeros of the random energy model *Physica A* **177** 24–30

[20] Ogure K and Kabashima Y 2005 An exact analytic continuation to complex replica number in the discrete random energy model of finite system size *Prog. Theor. Phys. Suppl.* **157** 103–6

[21] Jana N K 2007 Contributions to random energy models *PhD Thesis* Indian Statistical Institute-Kolkata arXiv:0711.1249 [math.PR]

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

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

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

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

[26] Bovier A and Kurkova I A 2003 Rigorous results on some simple spin glass models *Markov Proc. Related Fields* **9** 209 http://math-mprf.org/journal/articles/id964/

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

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

[29] Derrida B and Flyvbjerg H 1987 Statistical properties of randomly broken objects and of multivalley structures in disordered systems *J. Phys. A: Math. Gen.* **20** 5273–88

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

[31] Gardner E and Derrida B 1989 The probability distribution of the partition function of the random energy model *J. Phys. A: Math. Gen.* **22** 1975–81

[32] Ogure K and Kabashima Y 2009 On analyticity with respect to the replica number in random energy models: I. An exact expression for the moment of the partition function *J. Stat. Mech. P03010*

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

[34] Dotsenko V 2011 Replica solution of the random energy model *Europhys. Lett.* **95** 50006

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

[36] Sasaki M and Martin O C 2002 Temperature chaos in a replica-symmetry-broken spin glass model: a hierarchical model with temperature chaos *Europhys. Lett.* **60** 316–22

[37] Krzakala F and Martin O C 2001 Discrete energy landscapes and replica symmetry breaking at zero temperature *Europhys. Lett.* **53** 749–55

[38] Gardner E and Derrida B 1988 Optimal storage properties of neural network models *J. Phys. A: Math. Gen.* **21** 271–84

[39] Gardner E and Derrida B 1989 Three unfinished works on the optimal storage capacity of networks *J. Phys. A: Math. Gen.* **22** 1983–94
[40] Krauth W and Mézard M 1989 Storage capacity of memory networks with binary couplings *J. Phys. France* **50** 3057–66

[41] Huang H and Kabashima Y 2014 Origin of the computational hardness for learning with binary synapses *Phys. Rev. E* **90** 052813

[42] Ding J and Sun N 2019 Capacity lower bound for the Ising perceptron *Proc. 51st Annual ACM SIGACT Symp. on Theory of Computing (STOC '19)* ed M Charikar and E Cohen pp 816–27

[43] Monasson R and Zecchina R 1997 Statistical mechanics of the random K-satisfiability model *Phys. Rev. E* **56** 1357–70

[44] Mézard M, Parisi G and Zecchina R 2002 Analytic and algorithmic solution of random satisfiability problems *Science* **297** 812–5

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

[46] Bramson M D 1983 *Convergence of Solutions of the Kolmogorov Equation to Travelling Waves* vol 285 (Providence, RI: American Mathematical Society)

[47] Bramson M D 1978 Maximal displacement of branching Brownian motion *Commun. Pure Appl. Math.* **31** 531–81

[48] Majumdar S N and Krapivsky P L 2000 Extremal paths on a random Cayley tree *Phys. Rev. E* **62** 7735–42