Phase Transitions in the Coloring of Random Graphs
Lenka Zdeborová1 and Florent Krza̧ka la2
LPTMS, UMR 8626 CNRS et Univ. Paris-Sud,
91405 Orsay CEDEX, France
PCT, UMR Gulliver 7083 CNRS-ESPCI,
10 rue Vauquelin, 75231 Paris, France
We consider the problem of coloring the vertices of a large sparse random graph with a given
number of colors so that no adjacent vertices have the same color. Using the cavity method, we
present a detailed and systematic analytical study of the space of proper colorings (solutions).
We show that for a fixed number of colors and as the average vertex degree (number of constraints)
increases, the set of solutions undergoes several phase transitions similar to those observed in the
mean field theory of glasses. First, at the clustering transition, the entropically dominant part of the
phase space decomposes into an exponential number of pure states so that beyond this transition
a uniform sampling of solutions becomes hard. Afterward, the space of solutions condenses over a
finite number of the largest states and consequently the total entropy of solutions becomes smaller
than the annealed one. Another transition takes place when in all the entropically dominant states
a finite fraction of nodes freezes so that each of these nodes is allowed a single color in all the
solutions inside the state. Eventually, above the coloring threshold, no more solutions are available.
We compute all the critical connectivities for Erdős-Rényi and regular random graphs and determine
their asymptotic values for large number of colors.
Finally, we discuss the algorithmic consequences of our findings. We argue that the onset of
computational hardness is not associated with the clustering transition and we suggest instead that
the freezing transition might be the relevant phenomenon. We also discuss the performance of a
simple local Walk-COL algorithm and of the belief propagation algorithm in the light of our results.
PACS numbers: 89.20.Ff, 75.10.Nr, 05.70.Fh, 02.70.-c
I. INTRODUCTION
Graph coloring is a famous yet basic problem in combinatorics. Given a graph and q colors, the problem consists
in coloring the vertices in such a way that no connected vertices have the same color [1]. The celebrated four-
colors theorem assures that this is always possible for planar graphs using only four colors [2]. For general graphs,
however, the problem can be extremely hard to solve and is known to be NP-complete [3], so that it is widely
believed that no algorithm can decide in a polynomial time (with respect to the size of the graph) if a given arbitrary
instance is colorable or not. Indeed, the problem is often taken as a benchmark for the evaluation of the performance
of algorithms in computer science. It has also important practical application as timetabling, scheduling, register
allocation in compilers or frequency assignment in mobile radios.
In this paper, we study colorings of sparse random graphs [4, 5]. Random graphs are one of the most fundamental
source of challenging problems in graph theory since the seminal work of Erdős and Rényi [6] in 1959. Concerning
the coloring problem, a crucial observation was made by focusing on typical instances drawn from the ensemble
of random graphs with a given average vertex connectivity c, as c increases a threshold phenomenon is observed.
Bellow a critical value cs a proper coloring of the graph with q colors exists with a probability going to one in
the large size limit, while beyond cs it does not exist in the same sense. This sharp transition also appears in
other Constraint Satisfaction Problems (CSPs) such as the satisfiability of Boolean formulae [1]. The existence
of the sharp COLorable/UNCOLorable (COL/UNCOL) transition was partially [85] proven in [7], and computing
rigorously its precise location is a major open problem in graph theory. Many upper and lower bounds were established
[8, 9, 10, 11, 12, 13, 14, 15] for Erdős-Rényi and regular random graphs.
It was also observed empirically [16, 17] that deciding colorability becomes on average much harder near to the
coloring threshold cs than far away from it. This onset of computational hardness cannot be explained only by the
simple fact that near to the colorable threshold the number of proper colorings is small [18]. Some progress in the
theoretical understanding has been done by the analysis of search algorithms [19, 20]. For the coloring problem, it
was proven [21] that a simple algorithm q-colors almost surely in linear time random graphs of average connectivity
c ≤ q log q − 3q/2 for all q ≥ 3 (see [21] for references on previous works). For 3-coloring the best algorithmic lower
bound is c = 4.03 [10]. An important and interesting open question [22] is the existence of an ǫ > 0 and of a polynomial
algorithm which q-colors almost surely a random graph of connectivity c = (1 + ε)q log q for arbitrary large q.
The sharp coloring threshold and the onset of hardness in its vicinity has also triggered a lot of interest within the
http://arXiv.org/abs/0704.1269v2
COL/UNCOLRigidityClustering Condensation
cd cc cr cs
Average connectivity
FIG. 1: A sketch of the set of solutions when the average connectivity (degree) is increased. At low connectivities (on the
left), all solutions are in a single cluster. For larger c, clusters of solutions appear but the single giant cluster still exists and
dominates the measure. At the dynamic/clustering transition cd, the phase space slits in an exponential number of clusters. At
the condensation/Kauzmann transition cc, the measure condenses over the largest of them. Finally, no solutions exist above the
COL/UNCOL transition cs. The rigidity/freezing transition cr (which might come before or after the condensation transition)
takes place when the dominating clusters start to contain frozen variables (dominating clusters is a minimal set of clusters such
that it covers almost all proper colorings). The clusters containing frozen variables are colored in black and those which do not
are colored in grey.
statistical physics community following the discovery of a close relation between constraint satisfaction problems and
spin glasses [23, 24]. In physical terms, solving a CSP consists in finding groundstates of zero energy. The limit of
infinitely large graphs corresponds to the thermodynamic limit. In the case of the q−coloring problem, for instance,
one studies the zero temperature behavior of the anti-ferromagnetic Potts model [25]. Using this correspondence,
a powerful heuristic tool called the cavity method has been developed [24, 26, 27, 28]; it allows an exact analytical
study of the CSPs on sparse random graphs. Unfortunately, some pieces are still missing to make the cavity method
a rigorous tool although many of its results were rigorously proven. The cavity method is equivalent to the famous
replica method [29] (and unfortunately for clarity, it has also inherited some of its notations, as we shall see).
The cavity method and the statistical physics approach have been used to study the q-coloring of random graphs
in [30, 31, 32]. The coloring threshold cs was calculated [30], the self-consistency of the solution checked [31] and the
large q asymptotics of the coloring threshold computed [31]. These results are believed to be exact although proving
their validity rigorously remains a major subject in the field. Nevertheless —as the results obtained for the random
satisfiability problem [28, 33, 34]— they agree perfectly with rigorous mathematical bounds [8, 9, 10, 11, 12, 13, 15],
and with numerical simulations. The coloring threshold is thus fairly well understood, at least at the level of cavity
method.
A maybe more interesting outcome of the statistical physics analysis of the CSPs was the identification of a new
transition which concerns the structure of the set of solutions, and that appears before the coloring threshold [24,
28, 35]: while at low connectivities all solutions are in a single pure state (cluster) [86], the set of solutions splits
in an exponential number of different states (clusters) at a connectivity strictly smaller than cs. Roughly speaking,
clusters are groups of near-by solutions that are in some sense disconnected from each other. Recently, the existence
of the clustered phase was proven rigorously in some cases for the satisfiability problem [36, 37]. A major step was
made by applying the cavity equations on a single instance: this led to the development of a very efficient message-
passing algorithm called Survey Propagation (SP) that was originally used for the satisfaction problem in [24] and later
adapted for the coloring problem in [30]. Survey propagation allows one to find solutions of large random instances
even in the clustered phase and very near to the coloring threshold.
Despite all this success, the cavity description of the clustered phase was not complete in many aspects, and a
first improvement has been made with the introduction of a refined zero temperature cavity formalism that allows
a more detailed description of the geometrical properties of the clusters [38, 39]. We pursue in this direction and
give for the first time a detailed characterization of the structure of the set of solutions. We observed in particular
that the clustering threshold was not correctly computed, that other important transitions were overlooked and
the global picture was mixed up. The corrected picture that we describe in this paper is the following: when the
connectivity is increased, the set of solutions undergoes several phase transitions similar to those observed in mean
field structural glasses (we sketch these successive transitions in fig. 1). First, the phase space of solutions decomposes
into an exponential number of states which are entropically negligible with respect to one large cluster. Then, at the
clustering threshold cd, even this large state decomposes into an exponential number of smaller states. Subsequently,
above the condensation threshold cc, most of the solutions are found in a finite number of the largest states. Eventually,
the connectivity cs is reached beyond which no more solutions exist. Another important transition, that we refer to
as rigidity, takes place at cr when a finite fraction of frozen variables appears inside the dominant pure states (those
containing almost all the solutions). All those transitions are sharp and we computed the values of the corresponding
critical connectivities.
A nontrivial ergodicity breaking takes place at the clustering transition, in consequence uniform sampling of solutions
becomes hard. On the other hand, clustering itself is not responsible for the hardness of finding a solution. Moreover,
until the condensation transition many results obtained by neglecting the clustering effect are correct. In particular for
all c < cc: i) the number of solutions is correctly given by the annealed entropy (and, for general CSP, by the replica
symmetric entropy), and ii) simple message passing algorithm such as Belief Propagation (BP) [40, 41] converges to
a set of exact marginals (i.e. the probability that a given node takes a given color). In consequence we can use BP
plus a decimation-like strategy to find proper colorings on a given graph. Finally we give some arguments to explain
why the rigidity transition is a better candidate for the onset of computational hardness in finding solutions.
Our results are obtained within the one-step replica symmetry breaking approach, and we believe (and argue
partially in section IV), that our results would not be modified by considering further steps of RSB (as opposed to
previous conclusions [31]).
A shorter and partial version of our results, together with a study of similar issues in the satisfiability problem, was
already published in [42]. We refer to [43] for a detailed discussion of the satisfiability problem. The paper is organized
as follows: In section II, we present the model. In section III, we introduce the cavity formalism at the so-called
replica symmetric level, and discuss in detail why and where this approach fails. In section IV we take into account the
existence of clusters of solutions and employ the so called one-step replica-symmetry breaking formalism to describe
the properties of clusters. The results for several ensembles of random graphs are then presented in section V. We
finally discuss the algorithmic implications of our findings in section VI and conclude by a general discussion. Some
appendixes with detailed computations complete the paper.
II. THE MODEL
A. Definition of the model
For the statistical physics analysis of the q−coloring problem [30, 32, 44] we consider a Potts [25] spin model with
anti-ferromagnetic interactions where each variable si (spin, node, vertex) is in one of the q different states (colors)
s = 1, . . . , q. Consider a graph G = (V , E) defined by its vertices V = {1, . . . , N} and edges (i, j) ∈ E which connect
pairs of vertices i, j ∈ V ; we write the Hamiltonian as
H({s}) =
(i,j)∈E
δ(si, sj) . (1)
With this choice there is no energy contribution for neighbors with different colors, but a positive contribution
otherwise. The ground state energy is thus zero if and only if the graph is q-colorable. The Hamiltonian leads to a
Gibbs measure [45] over configurations (where β is the inverse temperature) :
µ({s}) = 1
e−βH({s}) , (2)
In the zero temperature limit, where β → ∞, this measure becomes uniform over all the proper colorings of the graph.
B. Ensembles of Random Graphs
We will consider ensembles of graphs that are given by a degree distribution Q(k). The required property of Q(k)
is that its parameters should not depend on the size of the graph. All the analytical results will concern only very
large sparse graphs (N → ∞). Provided the second and higher moments of Q(k) do not diverge, such graphs are
locally tree-like in this limit [4, 5]. More precisely, call a d-neighborhood of a node i the set of nodes which are at
distance at most d from i. For d arbitrary but finite the d-neighborhood is almost surely a tree graph when N → ∞.
This property is connected with the fact that the length ℓ of the shortest loops (up to a finite number of them) scales
with the graph size as ℓ ∼ log(N). We will consider the two following canonical degree distribution functions:
(i) Uniform degree distribution, Q(k) = δ(k− c), corresponding to the c-regular random graphs, where every vertex
has exactly c neighbors.
(ii) Poissonian degree distribution, Q(k) = e−cck/k!, corresponding to the Erdős and Rényi (ER) random graphs [6].
A simple way to generate graphs with N vertices from this ensemble is to consider that each link is present with
probability c/(N − 1). The binomial degree distribution converges to the Poissonian in the large size limit.
It will turn out that the cavity technics simplify considerably for the regular graphs. However, ER graphs have the
advantage that their average connectivity is a real number that can be continuously tuned, which is obviously very
convenient when one wants to study phase transitions. It is thus useful to introduce a third ensemble, which still has
the computational advantage of regular graphs, but that at the same time gives more freedom to vary the connectivity:
(iii) bi-regular random graphs, where nodes with connectivity c1 are all connected to nodes with connectivity c2, and
vice-versa. There are thus two sets of nodes with degree distributions Q(k) = δ(k − c1) and Q(k) = δ(k − c2).
Notice that these graphs are bipartite by definition and therefore have always a trivial 2-coloring which we will have
to dismiss in the following. This can be easily done within the cavity formalism (it is equivalent to neglecting the
ordered “crystal” phase in glass models [46]).
In the first two cases, the parameter c plays the role of the average connectivity, c =
kQ(k). In the cavity
approach, a very important quantity is the excess degree distribution Q1(k), i.e. the distribution of the number of
neighbors, different from j, of a vertex i adjacent to a random edge (ij):
Q1(k) =
(k + 1)Q(k + 1)
. (3)
This distribution remains Poissonian for Erdős-Rényi graphs, whereas the excess degree is equal to c− 1 in the case
of regular graphs. In the bi-regular case, there are two sets of nodes with excess degrees c1 − 1 and c2 − 1.
III. THE CAVITY FORMALISM AT THE REPLICA SYMMETRIC LEVEL
We start by reviewing the replica symmetric (RS) version of the cavity method [26, 27] and its implications for the
coloring problem. In the last part of the section we show when, and why, the RS approach fails.
A. The replica symmetric cavity equations
The coloring problem on a tree is solved exactly by an iterative method called the belief propagation algorithm [40]
(note some boundary conditions have to be imposed, otherwise a tree is always 2-colorable) that is equivalent to the
replica symmetric cavity method [41]. At this level, the method is in fact a classical tool in statistical physics to deal
with that tree structure that dates back to the original ideas of Bethe, Peierls and Onsager [47]. It allows one to
compute the marginal probabilities that a given node takes a given color and, in the language of statistical physics,
observables like energy, entropy, average magnetization, etc. The applicability of the method goes however beyond
tree graphs and we will discuss when it is correct for random locally tree-like graphs in section III C.
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
���k k2
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
�������
������
������
������
������
������
������
������
������
������
������
������
������
������
������
������
������
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
�����
������������
FIG. 2: Iterative construction of a tree by adding a new Potts spin.
Let us now describe the RS cavity equations. Denote ψi→jsi the probability that the spin i has color si when the
edge (ij) is not present and consider the iterative construction of a tree in fig. 2. The probability follows a recursion:
ψi→jsi =
k∈i−j
e−βδsiskψk→isk =
k∈i−j
1 − (1 − e−β)ψk→isi
. (4)
where Z
0 is a normalization constant (the cavity partition sum) and β the inverse temperature. The notation
k ∈ i − j means the set of neighbors of i except j. The normalization Zi→j0 is related to the free energy shift after
the addition of the node i and the edges around it, except (ij), as
0 = e
−β∆F i→j . (5)
In the same way, the free energy shift after the addition of the node i and all the edges around it is
e−β∆F
= Zi0 =
1 − (1 − e−β)ψk→is
. (6)
The Zi0 is at the same time the normalization of the total probability that a node i takes color si (the marginal of i):
ψisi =
1 − (1 − e−β)ψk→isi
. (7)
The free energy shift after the addition of an edge (ij) is
e−β∆F
0 = 1 − (1 − e−β)
ψj→is ψ
s . (8)
The free energy density in the thermodynamic limit is then given by (see for instance [26])
f(β) =
∆F i −
∆F ij
 . (9)
Note that this relation for the free energy is variational, i.e. that if one differentiates with respect to β, then only the
explicit dependence needs to be considered (see [26] for details). The energy (the number of contradictions) and the
corresponding entropy (the logarithm of the number of solutions with a given energy) densities can be then computed
from the Legendre transform
− βf(β) = −βe+ s(e) , (10)
were f = F/N , e = E/B and s = S/N are intensive variables.
The learned reader will notice that in some previous works using the cavity method [30, 31], these equations were
often written for a different object than the probabilities ψ. Instead, the so-called cavity magnetic fields h and biases
u where considered. The two approaches are related via
−βhi→jsi ≡ Zi→j0 ψi→jsi ≡
k∈i−j
−βuk→isi . (11)
Strictly speaking the ψ are “cavity probabilities” while h and u are “cavity fields”, however ψ are sometimes also
referred to as cavity fields (or messages) in the literature, and the reader will thus forgive us if we do so.
Note that each of the two notations is suitable for performing the zero temperature limit in a different way: in the
first one we fix the energy zero and we obtain the zero temperature BP recursion which gives the marginal probabilities
for each variable, while in the second case, we obtain a simpler recursion called Warning Propagation (49) that deals
with the energetic contributions but neglects the entropic ones [28, 30]. This is the origin of the discrepancy between
the RS results of refs. [32] and [30]. As we shall see, the differences between these two limits will be an important
point in the paper.
B. Average over the ensemble of graphs and the RS solution
To compute (quenched) averages of observables over the considered ensemble of random graphs given by the degree
distribution Q(k) we need to solve the self-consistent cavity functional equation
P(ψ) =
Q1(k)
dψiP(ψi)δ(ψ −F({ψi})), (12)
where Q1(k) is the distribution of the number of neighbors given that there is already one neighbor, and the function
F({ψi}) is given by eq.(4). Beware that ψ here is a q-components vector while we omit the vector notation to lighten
the reading. This equation is quite complicated since the order parameter is nontrivial, but we can solve it numerically
using the population dynamics method described in [26, 27].
Throughout this paper we will search only for the color symmetric functions P(ψ), i.e. invariant under permutation
of colors. Clearly with this assumption we might miss some solutions of (12). Consider for example q > 2 colors and
the ensemble of random bi-regular graphs. Since every bipartite graph is 2-colorable there are q(q−1)/2 corresponding
color asymmetric solutions for P(ψ). For the ensemble of random graphs considered here, we later argue that this
assumption is however justified.
Another important observation is that for regular graphs the equation (12) crucially simplifies: the solution factor-
izes [48] in the sense that the order parameter ψ is the same for each edge. This is due to the fact that, locally, every
edge in such a regular graph has the same environment. All edges are therefore equivalent and thus the distribution
P(ψ) has to be a delta function. For the bi-regular graphs, the solution of (12) also factorizes, but the two sets of
nodes of connectivity c1 and c2 (each of them being connected to the other) have to be considered separately.
It is immediate to observe that P(ψ) = δ(ψ−1/q) (i.e. each of the q components of each cavity field ψ equals 1/q),
is always a solution of (12). By analogy with magnetic systems we shall call this solution paramagnetic. Numerically,
we do not find any other solution in the colorable phase. For regular random graphs the paramagnetic solution
is actually the only factorized one. The number of proper colorings predicted by the RS approach is thus easy to
compute. Since all messages are of the type P(ψ) = δ(ψ − 1/q), the free energy density simply reads
− βfRS = log q +
1 − 1 − e
. (13)
The entropy density at zero temperature thus follows
sRS = log q +
1 − 1
. (14)
It coincides precisely with the annealed (first moment) entropy. We will see in the following that, surprisingly, the
validity of this formula goes well beyond the RS phase (actually until the so-called condensation transition).
C. Validity conditions of the replica symmetric solution
We used the main assumption of the replica symmetric approach when we wrote eq. (4): we supposed that the
cavity probabilities ψk→isi for the neighbors k of the node i are “sufficiently” independent in absence of the node i,
because only then the joint probability factorizes. This assumption would be true if the lattice were a tree with
non-correlated boundary conditions, but loops, or correlations in the boundaries, may create correlations between the
neighbors of node i (in absence of i) and the RS cavity assumption might thus cease to be valid in a general graph.
The aim here is to precise and quantify this statement both from a rigorous and heuristic point of view.
1. The Gibbs measure uniqueness condition
Proving rigorously the correctness of the RS cavity assumption for random graphs is a crucial step that has not yet
been successfully overcome in most cases. The only success so far was obtained by proving a far too strong condition:
the Gibbs measure uniqueness [49, 50, 51, 52]. Roughly speaking: the Gibbs measure (2) is unique if the behavior
of a spin i is totally independent from the boundary conditions (i.e. very distant spins) for any possible boundary
conditions. More precisely, let us define {sl} colors of all the spins at distance at least l from the spin i. The Gibbs
measure µ is unique if and only if the following condition holds for every i (and in the limit N → ∞)
{sl},{s
|µ(si|{sl}) − µ(si|{s′l})|
0 , (15)
where the average is over the ensemble of graphs. In [50, 51], it was proven that the Gibbs measure in the coloring
problem on random regular graphs is unique only for graphs of degree c < q. To the best of our knowledge, this has
not been computed for Erdős-Rényi graphs (later in this section, we argue on the basis of a physical argument that
it should be c < q − 1 in this case).
2. The Gibbs measure extremality condition
In many cases, the RS approach is observed to be correct beyond the uniqueness threshold. It was thus suggested
in [53] (see also [42]) that the Gibbs measure extremality provides a proper criterion for the correctness of the replica
symmetric assumption. Roughly speaking, the difference between uniqueness and extremality of a Gibbs measure is
that although there may exist boundary conditions for which the spin i is behaving differently than for others, such
boundary conditions have a null measure if the extremality condition is fulfilled. Formally (and keeping the notations
from the previous section), the extremality corresponds to
µ({sl})
|µ(si|{sl}) − µ(si)|
0 . (16)
In mathematics the “extremal Gibbs measure” is often used as a synonym for a “pure state”. Recently, the authors
of [53] provided rigorous bounds for the Gibbs measure extremality of the coloring problem on trees.
There exist two heuristic equivalent approaches to check the extremality condition. In the first one, one studies
the divergence of the so called “point-to-set” correlation length [54, 55]. The second one is directly related to the
cavity formalism: one should check for the existence of a nontrivial solution of the one-step replica symmetry breaking
equations (1RSB) at m = 1 (see section V). Both these analogies were remarked in [53] and exploited in [42]. We will
show in section V that the extremality condition ceases to be valid at the clustering threshold cd, beyond which the
1RSB formalism will be needed.
3. The local stability: a simple self-consistency check
A necessary, simple to compute but not sufficient, validity condition for the RS assumption is the non-divergence
of the spin glass susceptibility (see for instance [56]). If it diverges, a spin glass transition occurs, and the replica
symmetry has to be broken [29]. The local stability analysis thus gives an upper bound to the Gibbs extremality
condition, which remarkably coincides with the rigorous upper bound of [57]. This susceptibility is defined as
χSG =
〈sisj〉2c . (17)
The connectivities above which it diverges at zero temperature can be computed exactly within the cavity formalism
(we refer to appendix A for the derivation). It follows for regular and Erdős-Rényi graphs:
RS = q
2 − 2q + 2 , cERRS = q2 − 2q + 1 , (18)
while the stability of the bi-regular graphs of connectivities c1, c2 is equivalent to regular graphs with c = 1 +
(c1 − 1)(c2 − 1).
Note that for regular and ER graphs the RS instability threshold is in the colorable phase only for q = 3. Indeed the
5-regular graphs are 3-colorable [31] (and rigorous results in [13, 14]) and exactly critical since c
RS(3) = 5, and for ER
graphs the COL/UNCOL transition appears at cs ≈ 4.69 [30] while cERRS (3) = 4. This means that the replica symmetry
breaking transition appears continuously at the point cRS so that above it the RS approach is not valid anymore. For
all q ≥ 4, however, the local stability point is found beyond the best upper bound on the coloring threshold for both
regular and ER graphs. In this case, the extremality condition will not be violated by the continuous mechanism, but
we will see that, instead, a discontinuous phase transition, as happens in mean field structural glasses, will take place.
Interestingly enough, a similar computation can be made for the ferromagnetic susceptibility χF =
〈sisj〉c
(see again appendix A). It diverges at c = q for regular graph and c = q− 1 for Erdős-Rényi graphs. This divergence
(called modulation instability in [56]) would announce the transition towards an anti-ferromagnetic ordering on a tree,
which is however incompatible with the frustrating loops in a random graph (although such order might exist on the
bi-regular graphs). This is precisely the solution we dismiss when considering only the color symmetric solution of
(12). Note however that the presence of this instability shows that the problem ceases to a have a unique Gibbs
state (although it is still extremal) as for some specific (and well-chosen) boundary conditions, an anti-ferromagnetic
solution may appear. Indeed it coincides perfectly with the rigorous uniqueness condition c = q for regular graph,
and suggests strongly that the uniqueness threshold (or at least an upper bound) for Erdős-Rényi graphs is c = q− 1.
IV. ONE-STEP REPLICA SYMMETRY BREAKING FRAMEWORK
So far we described the RS cavity method for coloring random graphs and explained that the extremality of the Gibbs
measure gives a validity criterion. We now describe the one-step replica symmetry breaking cavity solution [26, 27].
In this approach, the non-extremality of the Gibbs measure is cured, by decomposing it into several parts (pure states,
clusters) in such a way that within each of the states the Gibbs measure becomes again extremal.
This decomposition has many elements/states, not just a finite numbers like the q-states of the usual ferromagnetic
Potts models. It is actually found that the number of pure states is growing exponentially with the size of the system.
Let us define the state-entropy function Σ(f) —called the complexity— which is just the logarithm of the number of
states with internal free energy density f , i.e. N (f) = exp[NΣ(f)]. In the glass transition formalism, this complexity
is usually referred to as the configurational entropy. Dealing with exponentially many pure states is obviously a
nightmare for all known rigorous approaches to the thermodynamic limit. The heuristic cavity method overcomes
this problem elegantly, as was shown originally in the seminal work of [26, 27].
Another very useful intuition about the 1RSB cavity method comes from the identification of states α with the
fixed points {ψ} of the belief propagation equations (4). The goal is thus to compute the statistical properties of these
fixed points. Each of the states is weighted by the corresponding free energy (9) to the power m, where m is just
a parameter analogous to the inverse temperature β (in the literature m is often called the Parisi replica symmetry
breaking parameter [29, 58]). The probability measure over states {ψ} is then
µ̃({ψ}) = Z0({ψ})
e−βmNf({ψ}) , (19)
where Z1 is just the normalization constant. To write the analog of the belief propagation equations we need to define
the probability (distribution) P i→j(ψi→j) of the fields ψi→j . This can be obtained from those of incoming fields as
P i→j(ψi→jsi ) =
k∈i−j
dψk→isi P
k→i(ψk→isi )δ(ψ
−F({ψk→isi }))
δ(ψ −F)
. (20)
The function F is given by eq. (4) and the delta function ensures that the set of fields ψi→j is a fixed point of the
belief propagation (4). The re-weighting term
takes into account the change of the free energy of a state
after the addition of a cavity spin i and its adjacent edges except (ij), as defined in eq. (5). This term appears for the
same reason as a Boltzmann factor e−βδsi,sk in eq. (4): it ensures that the state α is weighted by (Zα)
m in the same
way a configuration {s} is weighted by e−βH({s}) in (2). Finally Zi→j1 is a normalization constant. In the second line
of (20) we introduced an abbreviation that will be used from now on to make the equations more easily readable.
The notation POP comes from “population dynamics” which refers to the numerical method we use to solve eq. (20).
The probability distribution P (ψ) can be represented numerically by a set of fields taken from P (ψ), and then the
probability measure P (ψ)dψ becomes uniform sampling from this set, for more details see appendix D.
We define the “replicated free energy” and compute it in analogy with eq. (9) as
Φ(β,m) ≡ − 1
log(Z1) =
∆Φi −
 , (21)
where
e−βm∆Φ
e−βm∆F
, e−βm∆Φ
e−βm∆F
. (22)
Putting together (19) and (21) we have
Z1 = e
−βmNΦ(β,m) =
e−βmNf({ψ}) =
df e−N [βmf(β)−Σ(f)], (23)
where the sum over {ψ} is over all states (or BP fixed points). In the interpretation of [59] m is the number of replicas
of the system, thus the name “replicated free energy” for Φ(β,m). Note that we are using the word “replica” only
to refer to the established terminology as no replicas are needed within the cavity formalism. From the saddle point
method, it follows that the Legendre transform of complexity function Σ(f) gives the replicated free energy Φ(m)
− βmΦ(β,m) = −βmf(β) + Σ(f). (24)
Notice that this equation is correct only in the highest order in the system size N , i.e. in densities and at the
thermodynamic limit. From the properties of the Legendre transform we have
Σ = βm2∂mΦ(β,m) , f = ∂m[mΦ(β,m)] , βm = ∂fΣ(f). (25)
Thus, from eq.(21), the free energy reads
f(β) =
∆F ie−βm∆F
e−βm∆F
∆F ije−βm∆F
e−βm∆F
. (26)
When the parameter m is equal to one (the number of replicas is actually one in the approach of [59]), then
−βΦ(β, 1) = −βf(β) + Σ(f) reduces to the usual free energy function considered in the RS approximation
Φ(β, 1) = e− Σ + s
= e− Tstot , (27)
where s in the internal entropy density of the corresponding clusters and stot the total entropy density of the system.
A. Analyzing the 1RSB equations
Combining (21) and (26) we can compute Σ and f for each value m, that gives us implicitly Σ(f). To compute
the thermodynamic observables in the model we have to minimize the total free energy ftot = f(β) − Σ/β over such
values of f where the complexity Σ(f) is non-negative (so that the states exist in the thermodynamic limit). The
minimum of the total free energy corresponds to a value of parameter m = m∗ and states with the free energy f∗
dominate the thermodynamics. Three different cases are then observed:
a) If there is only the trivial (replica symmetric) solution at m = 1, then the Gibbs measure (2) is extremal and
the replica symmetric approach is correct. If at the same time a nontrivial solution exists for some m 6= 1, then
the clusters corresponding to this solution have no influence on the thermodynamics.
b) If there is a nontrivial solution at m = 1 with a positive complexity, then m∗ = 1 minimizes the total free
energy. The system is in a clustered phase with an exponential number of dominating states.
c) If however the complexity is negative at m = 1, then the corresponding states are absent with probability
one in the thermodynamic limit. Instead the total entropy is dominated by clusters corresponding to m∗ such
that Σ(m∗) = 0: the system is in a condensed phase. Note that the condition Σ(m∗) = 0 corresponds to the
maximum of the replicated free energy (21).
The transitions between these cases are well known in structural glass phenomenology where they appear when the
temperature is lowered [60, 61]. The transition from the paramagnetic phase to the clustered one is usually referred
to as the dynamical transition [62] or the clustering transition. It is not a true thermodynamic transition as the total
free energy of the system at m∗ = 1 is still equal to the replica symmetric one (9) (see appendix C) and thus is an
analytical function of connectivity. However, the phase space is broken into exponentially many components and, as
a consequence, the dynamics fall out-of-equilibrium beyond this transition.
The second transition from the clustered to the condensed phase is, however, a genuine thermodynamic transition
(the free energy has a discontinuity in the second derivative at cc) and is called the replica symmetry breaking
transition, or the static glass transition. At this point the measure condenses into few clusters, and we shall call it
the condensation transition. In structural glasses, it corresponds to the well known Kauzmann transition [63]. The
sizes of the clusters in the condensed phase follow the so called Poisson-Dirichlet process which is discussed shortly
in appendix B.
The procedure to compute the replicated free-energy (21) and the related observables was described above for a
single large random graph. To compute the averages over the ensemble of random graphs, we need to solve an equation
analogous to eq. (12)
P [P (ψ)] =
Q1(k)
dP i(ψi)P [P i(ψi)] δ(P (ψ) −F2({P i(ψi)})), (28)
where the functional F2 is given by eq. (20). Solving this equation for a general ensemble of random graphs and a
general parameter m is a numerically quite tedious problem. In the population dynamics algorithm [26, 27] we need
to deal with a population of populations of q-components fields. It is much more convenient to look at the ensemble
of random regular graphs where a factorized solution P [P (ψ)] = δ(P (ψ) − P0(ψ)) must exists. Then we are left with
only one functional equation (20).
Before discussing the zero temperature limit, we would like to point out that there exists another very important
case in which eq. (28) simplifies. For m = 1, as first remarked and proved in [53], when dealing with the problem of
reconstruction on trees, the equations can be written (and numerically solved) in a much simpler way. We again refer
to the appendix C for details. Especially for the Poissonian random graphs this simplification is very useful.
B. Zero temperature limit
We now consider the zero temperature limit β → ∞ of the 1RSB equations (24)-(28) to study the coloring problem.
In most of the previous works [27, 28, 30] the energetic zero temperature limit was employed. The β → ∞ limit of eq.
(24) was taken in such a way that mβ = y remains constant. The replicated free energy (24) then becomes
− yΦe(y) = −ye+ Σ(e) . (29)
It is within this approach that the survey propagation (SP) algorithm was derived. The connectivity at which the
complexity function Σ(e = 0) becomes negative is the coloring threshold. Above this connectivity Σ(e) was used to
compute the minimal number of violated constraints (the ground state energy). The reweighting in eq. (20) becomes
e−y∆E
, and when y → ∞ all the configurations with positive energy are forbidden.
In this paper we adopt the entropic zero temperature limit, suggested originally in [38, 39]. The difference in the
two approaches was already underlined in sec. III A. We want to study the structure of proper colorings, i.e. the
configurations of zero energy, and we thus fix the energy to zero. Then we obtain the entropy by considering −βf = s
and introduce a free entropy —or in replica term a “replicated entropy”— as Φs(m) = −βmΦ(β,m)|β→∞. Eq. (24)
then becomes
Φs(m) = ms+ Σ(s). (30)
The belief propagation update (4) becomes
ψi→jsi =
k∈i−j
1 − ψk→isi
, (31)
while the 1RSB equation (20) keeps the same expression (and thus the same computational complexity).
The partition sum Z0 in (2) becomes in this limit the number of proper colorings or solutions. The clusters are now
sets of such solutions, and are weighted by their size to the power m. The free entropy Φs(m) is then computed as
Φs(m) =
∆Φis −
∆Φijs
 , (32)
where ∆Zi and ∆Zij are given by eqs. (5) and (6). The analysis from the previous section is valid also for the entropic
zero temperature limit. The information extracted from the number of clusters of a given size, Σ(s), is one of the
main results of this paper and will be discussed and interpreted in section V.
C. The role of frozen variables
In this section we discuss the presence and the role of the frozen variables and explain the connection between the
energetic and the entropic zero temperature limits. This allows us to revisit (and extend) the survey propagation
equations. Remember that the components of the cavity field ψi→jsi are the probabilities that the node i takes the
color si when the constraint on the edge ij is not present. In the zero temperature limit we can classify them in two
categories:
(i) A hard field corresponds to the case when all components of ψi→j are zero except one, s. Then only that color
is allowed for the spin i, in absence of edge (ij).
(ii) A soft field corresponds to the case when more than one component of ψi→jsi is nonzero. The variable i is thus
not frozen in absence of edge (ij), and the colors of all the nonzero components are allowed.
This distinction is also meaningful for the full probabilities ψisi (4), if ψ
is a hard field then the variable i is
frozen. In the colorable region there cannot exist a finite fraction of frozen variables (even if we consider properly
the permutational symmetry of colors) since by adding a link the connectivity changes by 1/N but the probability
of becoming uncolorable would be finite. On the contrary, in the 1RSB picture, we observe that a finite fraction of
variables can be frozen within a single cluster. In other words, in all the solutions that belong to this given cluster a
finite fraction of variables can take one color only. By adding a link into the graph, the connectivity grows by 1/N ,
and there is a finite probability that a cluster with frozen variables disappears. The distinction between hard and soft
fields is useful not only for the intuition about clusters, but also for the analysis of the cavity equations and it also
leads to the survey propagation algorithm.
The distribution of fields over states P i→j(ψi→j) (20) can be decomposed into the hard- and soft-field parts
P i→j(ψi→j) =
ηi→js δ(ψ
i→j − rs) + ηi→j0 P̃ i→j(ψi→j) , (33)
where P̃ i→j is the distribution of the soft fields and the normalization is
s=0 η
s = 1.
Interestingly, the presence of frozen variables in the entropically dominating clusters is connected to the divergence
of the size of average minimal rearrangement [55, 64]. Precisely, choose a random proper coloring {s} and a random
node i in the graph. The average minimal rearrangement is the Hamming distance to the nearest solution in which
node i has a color different from si averaged over the nodes i, the proper colorings, and graphs in the ensemble.
Another interesting role of the frozen variables arises within the whitening procedure, introduced in [65] and studied,
between others, for the satisfiability problem in [66, 67]. This procedure is equivalent to the warning propagation
(WP) update (49) which we outlined in sec. III A. Whitening is able to identify if a solution belongs to a cluster with
frozen variables or not. Particularly, the result of the whitening is a set of hard cavity fields.
Since the survey propagation algorithm is computing statistics over the states that contain hard fields, then the
solution found after decimating the survey propagation result should a priori also contain hard fields. However, recent
works show that if one applies the whitening procedure starting from solutions found by SP on large graphs, whitening
converges every time to the trivial fixed point (see detailed studies for K-SAT in [66, 67]). A possible solution to this
apparent paradox is discussed in sec. VI.
1. Hard fields in the simplest case, m = 0
Let us now consider the survey propagation equations originally derived in [30] from the energetic zero temperature
limit (29) when y → ∞. For simplicity we will write them only for the 3-coloring. We consider the 1RSB cavity
equation (20) for m → 0, then the reweighting factor (Zi→j0 )m is equal to zero when the arriving hard fields are
contradictory, and equal to one otherwise. The update of probability ηs that a field is frozen in direction s is then
written from eq.(20):
ηi→js =
k∈i−j(1 − ηk→is ) −
k∈i−j(η
0 + η
p ) +
k∈i−j η
k∈i−j(1 − ηk→ip ) −
k∈i−j(η
0 + η
p ) +
k∈i−j η
. (34)
In the numerator there is a telescopic sum counting the probability that color s and only color s is not forbidden by
the incoming fields. In the denominator the telescopic sum is counting the probability that there is at least one color
which is not forbidden. If we do not want to actually find a proper coloring on a single graph but just to compute
the replicated free energy/entropy, we can further simplify eq. (34) by imposing the color symmetry. Indeed, the
probability that in a given state a field is hard in direction of a color s has to be independent of s (except s = 0 which
corresponds to a soft field). Then (34) becomes, now for general number of colors q:
ηi→j = w({ηk→i}) =
l=0 (−1)l
k∈i−j
1 − (l + 1)ηk→i
l=0 (−1)l
k∈i−j [1 − (l + 1)ηk→i]
. (35)
Note that since ∂Σ(s)/∂s = −m, the value m = 0 corresponds to the point where the function Σ(s) has a zero
slope. If a nontrivial solution of (35) exists, then Σ(s)|m=0 is the maximum of the curve Σ(s) and is counting the
total log-number of clusters of size s, which is due to the exponential dependence, also the total log-number of all
clusters, regardless their size. There are two points that we want to emphasize:
• Suppose that a nontrivial solution of (35) exists, i.e. many clusters exist and their number can be computed
with the energetic zero temperature limit calculations. Then the clusters might be very small and contain very
few solution in comparison to bigger less numerous clusters; or in comparison to a giant single cluster which
might still exist. This situation cannot be decided by the energetic formalism that weights clusters equally
independently of their size.
• Suppose, on contrary, that a nontrivial solution of (35) does not exist. It might still well be that many clusters
exist, but the Σ(s) curve has no part with zero slope.
We will see that these two cases are actually observed. The energetic method, that can locate the coloring threshold
and from which the survey propagation can be derived, is therefore not a good tool to study the clustering transition.
2. Generalized survey propagation recursion
Let us compute how the fraction of hard fields η evolves after one iteration of equation (20) at general m. There
are two steps in each iterations of (20). In the first step, η iterates via eq. (35). In the second step we re-weight the
fields. Writing P hardm (Z) the —unknown— distribution of the reweightings Z
m for the hard fields, one gets
ηi→j =
dZ P hardm (Z)Z
mw({ηk→i}) = w({η
k→i})
dZ P hardm (Z)Z
w({ηk→i})
hard. (36)
A similar equation can formally be written for the soft fields
1 − qηi→j = 1 − qw({η
k→i})
soft. (37)
Writing explicitly the normalization N , we finally obtain the generalized survey propagation equations:
ηi→j =
w({ηk→i})
qw({ηk→i}) + [1 − qw({ηk→i})] r({ηk→i}) , with r({η
k→i}) =
Zmsoft
Zmhard
. (38)
In order to do this recursion, the only information needed is the ratio r between between soft- and hard-field reweight-
ings, which is in general difficult to compute since it depends on the full distribution of soft fields.
There are two cases where eq. (38) simplifies so that the hard-field recursion become independent from the soft-field
distribution. The first case is, of course, m = 0 then r = 1 independently of the edge (ij), and the equation reduces
to the original SP. The second case arise for m = 1, where one can use the so-called reconstruction formalism and
obtain again a closed set of equations. The computation is done in appendix C, and the SP equations at m = 1 read
ηi→js =
(−1)l
s1,...,sl 6=s
k∈i−j
1 − q
q − 1
ηk→isα
 . (39)
It would be interesting in the future to use eq. (38) in an algorithm to find proper graph colorings, as it has been
done with the original SP equation [30]. As an approximation one might also use a value r independent of the edge
(ij), but different from one.
For the purpose of the present work, it is important to notice that it is also possible to use eq. (38) in the population
dynamics to simplify the numerical evaluation of the 1RSB solution by separating the hard-field and the soft-field
contributions. Indeed, it gives the exact density of hard fields provided the ratio r is calculated, which is doable
numerically (see appendix D). This allows us to monitor precisely the hard-field density and only the soft-field part
is given by the population dynamics. This turns out to greatly improve the precision of the numerical solution of the
cavity equations and to considerably fasten the code.
3. The presence of frozen variables
A natural question is: “When are the hard fields present?” or more precisely: “When does eq. (38) have a nontrivial
solution η > 0?” First notice that in order to constrain a node into one color, one needs at least q− 1 incoming fields
that forbids all the other colors. It means that function w({ηk→i}) defined in eq. (35) is identically zero for k < q− 1
and might be non-zero only for k ≥ q − 1, where k is the number of incoming fields.
In the limit r → 0 (which corresponds to m → −∞) eq. (38) gives η = 1/q if w({ηk→i}) is positive, and η = 0
if w({ηk→i}) is zero. Updating eq. (38) on a given graph, from initial conditions η = 1/q everywhere, is equivalent
to recursive removing of all the nodes of connectivity smaller than q. This shows that the first nontrivial solution
with hard fields exists if and only if the q-core [68] of the graph is extensive. For regular graphs it is simply at
connectivity c = q while for Erdős-Rényi graphs these critical connectivities can be computed exactly and read, for
small q, c3 = 3.35, c4 = 5.14, c5 = 6.81 [68]. Indeed we see that the first nontrivial solution to the 1RSB equation
appears much before those of the original SP equation at m = 0.
On a regular graph, the equations further simplify as η factorizes (is edge independent) and follows a simple
self-consistent equation
η = w(η)
qw(η) + [1 − qw(η)] r . (40)
This equation can be solved for every possible ratio r so that for all c ≥ q, we can compute and plot the curve η(r).
We show the results in fig. 3 for different numbers of colors q. On this plot, we observe that η = 1/q, as predicted,
for r = 0. It then gets smaller for larger value of the ratio and, at a critical value rcrit, the solution disappears
discontinuously and only the (trivial) solution η = 0 exists. The values rcrit correspond to a critical value of mr. For
all m > mr no solution containing frozen variables can exist.
 0.65
 0.75
 0.85
 0.95
 0  0.5  1  1.5  2  2.5  3  3.5  4
3-coloring of regular graphsc=3
SP, general
SP m=0
SP m=1
 0.65
 0.75
 0.85
 0.95
 0  0.2  0.4  0.6  0.8  1  1.2  1.4  1.6  1.8
4-coloring of regular graphs
c=6 c=7
SP, general
SP m=0
SP m=1
FIG. 3: The lines are solutions of eq. (40) and give the total fraction qη of hard fields for a given value of ratio r = Zmsoft/Z
for q = 3 (left) and q = 4 (right) in regular random graphs. There is a critical value of the ratio (full point) beyond which only
the trivial solution η = 0 exists. Note that the solutions at m = 0 and m = 1 only exist for a connectivity large enough.
D. Validity conditions of the 1RSB solution
Now that we have discussed in detail the 1RSB formalism the next question is: “Is this approach correct?” To
answer this question, one has to test if the Gibbs measure is extremal within the thermodynamically dominating
pure states. This is equivalent to checking if the two-step Replica Symmetry Breaking (2RSB) solution is non trivial.
Computing explicitly the 2RSB solution is however very complicated numerically, especially for Erdős-Rényi graphs.
Instead, the local stability of the 1RSB solution towards 2RSB should be checked, in analogy with the RS stability in
sec. III C 3. It is indeed a usual feature in spin glass physics to observe that the 1RSB glass phase become unstable
at low temperatures towards a more complex RSB phase and this phenomenon is called the Gardner transition [69].
To perform the stability analysis [31, 34, 56, 70], one first writes the 2RSB recursion —where the order parameter is
a distribution of distributions of fields on every edge P1(P2(ψ))— and then two types of 1RSB instabilities have to be
considered depending on the way the 2RSB arises from the 1RSB solution. The first type of instability —called states
aggregation— corresponds to δ(P (ψ)) → P1(P2(ψ)) while the second type —called states splitting— corresponds to
P (δ(ψ)) → P1(P2(ψ)). A complete stability analysis is left for future works, but it is worth to discuss the relevance
of the results derived over the last few years [31, 34, 70].
The 1RSB stability was studied for the coloring problem in [31] but only for the energetic zero temperature limit
(29). In this case the parameter y = βm is conjugated to the energy. The results derived in [31] —as well as those
previously derived for other problems [34, 70]— thus concern only the clusters of sizes corresponding to m = 0 at zero
(for y = ∞) or at positive (for finite y) energy. The main result of [31, 34, 70] was that the 1RSB approach was stable
in vicinity of the coloring threshold cs. As we shall see the clusters corresponding to m = 0 are those dominating
the total entropy at the coloring threshold and as a consequence its location is thus exact within the cavity approach.
The states of the lowest energy (the ground states) in the uncolorable phase also correspond to m = 0, and thus
the conclusions of [31, 34, 70] concerning the uncolorable phase are also correct. In particular, a Gardner transition
towards further steps of RSB appears in the uncolorable phase beyond a connectivity denoted cG in [31].
On the other hand in the colorable phase the stability of the entropically dominating clusters that correspond to
m > 0 should be investigated. Some more relevant information can be, however, already drawn from known results. It
was indeed found that the 1RSB approach at m = 0 is type I stable for all y, and type II stable for all y > yI in vicinity
of the coloring threshold. These results concerns the states of positive energy, but keeping in mind the interpretation
of y as a slope in T,m diagram, we see that the clusters of zero energy corresponding to small but nonzero positive m
and zero temperature are also stable with respect to both types of stabilities. Near the colorable threshold, the value
of m∗ which describes the dominating clusters is close to zero and as a consequence all the dominating clusters are
1RSB stable in vicinity of the coloring threshold. Far from the coloring threshold, the stability analysis of [31, 34] is
irrelevant. In particular, the predictions of a full-RSB colorable phase made in [31, 34, 70] is not correct. Quite the
contrary, our preliminary results indicate that all the dominating clusters are 1RSB stable for q > 3.
The 3−coloring is however a special case as the clustering transition is continuous. Although the type II instability
seems irrelevant in this case as well, we cannot at the moment dismiss a type I instability close to the clustering
transition. Indeed the entropically relevant clusters correspond to values of m∗ close to one in this case, and it is
simple to show that the clusters at m = 1 are type I unstable in the case there is a continuous transition: this is
because the type I stability is equivalent to the convergence of the 1RSB update on a single graph. Since for m = 1
the averages of the 1RSB fields satisfy the RS belief propagation equations, and since we know from the RS stability
analysis in section III C 3 that those equations do not converge in the RS unstable region (i.e. for c > cRS = 4 in
3-coloring of Erdős-Rényi graphs), it then follows that the 3-coloring is unstable against state aggregation at m = 1 for
all connectivities c > 4. Therefore, it is possible that the 1RSB result for 3-coloring are only approximative for what
concerns the number and the structure of solutions close to the clustering transition (note that the critical values for
the phase transition are however correct and do not depend on that). This, and related issues [71], will be hopefully
clarified in future works.
To conclude, we believe that all the transition points we discuss in this paper (and those computed in the K-SAT
problem in [42]) as well as the overall picture, are exact and would not be modified by considering further steps of
replica symmetry breaking.
V. THE COLORING OF RANDOM GRAPHS: CAVITY RESULTS
We now solve the 1RSB equations, discuss and interpret the results. We solve the equation (28) by the population
dynamics technique, the technical difficulties and the precision of the method are discussed in appendix D. Let
us stress at this point that the correctness of eq. (28) is guaranteed only in the limit of large graphs (N → ∞),
unfortunately the cavity method does not give us any direct hint about the finite graph-size corrections. We start by
the results for the regular random graphs, then consider ensemble of bi-regular graphs and after that we turn towards
Erdős-Rényi graphs. Finally, we discuss the limit of large number of colors.
A. Regular random graphs
Let us fix the number of colors q, vary the connectivity, and identify successively all the transitions that we shall
encounter. For the sake of the discussion, we choose as a typical example the 6-coloring and we discuss later in details
the cases, for different number of colors, where some transitions are missing or are arriving in a different order. We
solved the 1RSB equation (20) for regular graphs, where the distribution P i→j(ψ) is the same for every edge (ij) (see
appendix D) and plot the curves for Σ(s) we obtained doing so in fig. 4. We now describe the phase space of solutions
when the connectivity is increased:
1) At very low connectivities c < q, only the paramagnetic RS solution is found at all m. i.e. P (ψ) = δ(ψ − 1/q).
The phase space is made of a single RS cluster.
2) For larger connectivities c ≥ q, we saw in section IV C 3 that the 1RSB equations start to have nontrivial solutions
with hard fields in an interval [−∞,mr]. Interestingly, another nontrivial solution, without hard fields, can be
found numerically in an interval [ms,∞], and we shall call this one the soft-field solution. As the connectivity
-0.04
-0.02
 0.02
 0.04
 0.06
 0.08
 0  0.02  0.04  0.06  0.08  0.1  0.12  0.14  0.16  0.18
q=4, c=9
q=5, c=14
q=5, c=13
-0.05
 0.05
 0.15
 0  0.02  0.04  0.06  0.08  0.1  0.12  0.14
 0  0.1  0.2  0.3
FIG. 4: Complexity as a function of the internal state entropy for the q-coloring problem on random regular graphs of
connectivity c. The full line corresponds to the clusters where a finite fraction of hard fields (frozen variables) is present and
the dotted line to the clusters without hard fields. The circle signs the entropically dominating clusters. Left: (q = 4, c = 9) is
in the clustered phase; (q = 5, c = 13) is in a simple replica symmetric phase and (q = 5, c = 14) is in the condensed clustered
phase. Right: results for 6-coloring for connectivities 17 (RS), 18 (clustered), 19 (condensed) and 20 (uncolorable). For 4-,5-
and 6-coloring all the smaller connectivities are in the RS phase while all the larger one are uncolorable.
increases, we find that mr increases while ms decreases, so that the gap [mr,ms] where no nontrivial solution
exists it getting continuously smaller.
However, there is no nontrivial solution at m = 1 for connectivities smaller then cd (see fig. 4 for the example
of the 6-coloring at c = 17). This means that the Gibbs measure (2) is still extremal. In other words the large
RS state still exists and is entropically dominant (its entropy (14) is noted by a circle in fig. 4). Despite the
fact that an exponential number of clusters of solutions exist and that the SP equations converge to a nontrivial
result, a random proper coloring will almost surely belong to the large RS cluster.
3) If the connectivity is increased at and above the clustering threshold cd, a nontrivial solution with positive
complexity Σ is found at m = 1. In fig. 4, we see that this happens at cd = 18 for the regular 6-coloring.
At this point, the RS Gibbs measure (2) ceases to be extremal and the single large RS cluster splits into
exponentially numerous components. To cover almost all proper colorings we need to consider exponentially
many clusters N ∼ eNΣ(m∗=1). The probability that two random proper colorings belong to the same cluster
is going exponentially to zero with the system size. The connectivity cd is thus the true clustering (dynamic)
transition. This is not, however, a thermodynamic phase transition because the 1RSB total entropy reduces to
the RS entropy (14) at m = 1 which is analytical in c. Thus the RS approach gives a correct number of solution
and correct marginals as long as the complexity function at m = 1 is non-negative.
4) For even larger connectivities c ≥ cc, the complexity at Σ(m = 1) becomes negative, e.g. cc = 19 for 6-coloring.
It means that the clusters corresponding to m = 1 are absent with probability one. The total entropy is then
smaller than the RS/annealed one and is dominated by clusters corresponding to m∗ < 1 such that Σ(m∗) = 0.
The ordered weights of the entropically dominating clusters follow the Poisson-Dirichlet process (explained in
appendix B). As a consequence, the probability that two random proper colorings belong to the same cluster is
finite in the thermodynamic limit. Another way to describe the situation is that the entropy condenses into a
finite number of clusters. This condensation is a true thermodynamic transition, since the total entropy is non-
analytical at cc (there is a discontinuity in its second derivative with respect to connectivity). The condensation
is analogous to the static (Kauzmann) glass transition observed in mean field models of glasses [60, 61].
5) For connectivities c ≥ cs (cs = 20 for 6-coloring) even the maximum of the complexity Σ(m = 0) becomes
negative. In this case proper colorings are absent with probability going to one exponentially fast with the size
of the graph, and we are in the uncolorable phase.
It is useful to think of the growing connectivity as additions of the constraints into a fixed set of nodes. From this
point of view the set of solutions which exists at connectivity c gets smaller when new edges are introduced and the
connectivity increased. This translates into the cartoon in the introduction (fig. 1) where all the successive transitions
are represented. Finally, another important transition has to be considered:
6) There is a connectivity cr beyond which the measure is dominated by clusters that contain a finite fraction of
frozen variables. For the regular 6−coloring, cr = 19. We refer to this as the rigidity transition, by analogy
with [72, 73].
The presence or the absence of hard fields inside a given cluster is crucial: if a cluster contains only soft fields, then
after the addition of a small but finite fraction of new constraints, its size will get smaller (or it will split). If, however,
a cluster contains a finite fraction of frozen variable, then after adding a small but finite fraction of links the cluster
will almost surely disappear.
Since the connectivities of regular graphs are integer numbers, we define the dynamical threshold cd as the smallest
connectivity where a nontrivial 1RSB solution exists at m = 1, the condensation transition cc as the smallest connec-
tivity where complexity at m = 1 is negative, cr the smallest connectivity where hard fields are present at m
∗ and
the coloring threshold cs as the first uncolorable case. The scenario described here is observed for all cases of the
regular ensemble, although, since connectivities are integer, the transitions are not very well separated at small q. We
summarize the results in table I.
Note that for q > 3, the local RS stability discussed in section III C 3 is irrelevant in the colorable regime. The
only subtle case being for 3−coloring of 5-regular graphs where the RS solution is only marginally stable, i.e. the
spin glass correlation function goes to zero only algebraically instead of exponentially (from this point of view c = 5
would correspond to the critical point well known in the second order phase transitions). More interesting cases will
arise in the other ensembles of random graphs.
q cSP [31] cd [53] cr cc cs [31]
3 5 5+ - 6 6
4 9 9 - 10 10
5 13 14 14 14 15
6 17 18 19 19 20
7 21 23 - 25 25
8 26 29 30 31 31
9 31 34 36 37 37
10 36 39 42 43 44
20 91 101 105 116 117
q c m∗ mr ms
5 3 RS+ 0.12 1.2(1)
4 8 RS -0.03 3.4(1)
4 9 1 0.41 0.41
5 12 RS -0.02 3.7(1)
5 13 RS 0.20 2.0(1)
5 14 0.50 0.90 0.90
6 16 RS -0.02 4.3(1)
6 17 RS 0.05 3.2(1)
6 18 1 0.40 0.40
6 19 0.92 0.96 0.96
7 21 RS 0.01 4.7(1)
7 22 RS 0.17 3.2(1)
7 23 1 0.60 0.60
7 24 1 0.95 0.95
TABLE I: Left: The transition thresholds for regular random graphs: cSP is the smallest connectivity with a nontrivial solution
at m = 0; the clustering threshold cd is the smallest connectivity with a nontrivial solution at m = 1; the rigidity threshold
cr is the smallest connectivity at which hard fields are present in the dominant states, the condensation cc is the smallest
connectivity for which the complexity at m = 1 is negative and cs the smallest UNCOL connectivity. Note that 3−coloring
of 5−regular graphs is exactly critical for that cd = 5
+. The rigidity transition may not exist due to the discreteness of the
connectivities. Right: Values of m∗ (corresponding to the dominating clusters), and in the range of [−∞,mr] the hard-field
solution exists, in the range [ms,∞] the soft-field solution exists.
B. Results for the bi-regular ensemble
The bi-regular ensemble allows us to fine-tune the connectivity while preserving the factorization of the 1RSB
solution, which is crucial for the numerical precision. It is actually more correct to say that the solution is “bi-
factorized”, as all the messages going from the nodes with connectivity c1 to c2 are the same and the other way
around. The bi-regular ensemble allows us to describe with large precision two interesting cases, which reappear in
the Erdős-Rényi ensemble and which are not present in the regular ensemble (again, due to the discrete nature of
the connectivity). Let us remind here that bipartite graphs are always 2-colorable, but we consider only the color
symmetric cavity solutions and that is why we get a nontrivial result from this ensemble.
-0.02
-0.01
 0.01
 0.02
 0.03
 0.04
 0.05
 0.06
 0.1  0.12  0.14  0.16  0.18  0.2  0.22  0.24
4-coloring of bi-regular graphs
(5-21)
-0.01
 0.01
 0.02
 0.03
 0.04
 0.05
 0.2  0.22  0.24  0.26  0.28  0.3  0.32  0.34  0.36
(4,36)
(4,39)
(4,42)
(4,45)
4-coloring of bi-regular graphs
 0.003
 0.006
 0.33  0.34
(4,39)
FIG. 5: The complexity as a function of entropy for 4-coloring or bi-regular graphs. Left: 5-21-bi-regular graph, an example
where the entropy is dominated by clusters with soft fields while the gap in the curve Σ(s) still exists. Right: 4-c-bi-regular
graphs for c=36, 39, 42, 45. In all these cases the replica symmetric solution is locally unstable. In the dependence Σ(s) we
see an unphysical branch of the complexity which is zoomed in the inset for c = 39.
In fig. 5, the left picture is the result for the complexity as a function of entropy Σ(s) for 4-coloring of 5-21-bi-regular
graphs. The replica symmetric solution on this case is locally stable. We see clearly the gap between the hard-field
and the soft-field solution, and yet we are already beyond the clustering transition cd; actually the system is in the
condensed phase. This example is similar to what happens for the 4-coloring of Erdős-Rényi graphs.
The second interesting case, the right hand side of fig. 5, is given by the results for Σ(s) for the 4-coloring of
4-c-bi-regular graphs, which are RS unstable for c > 28. Both the clustering and the condensation transitions coincide
with the RS instability cd = cc = 28. The survey propagation equations have a nontrivial solution starting from
cSP = 37. The rigidity transition is at cr = 49. Finally the coloring threshold is cs = 57. Qualitatively, the results for
this 4-c-bi-regular ensemble are the same as those for the 3-coloring of Erdős-Rényi random graphs.
We see that for c ≤ 42 the gap between the hard-field (full line) and soft-field (dotted line) solution exists. For
m > ms there is a non-physical nontrivial soft-field solution, the convex part of the line in the figure, zoomed in the
inset. It means that for m < mr we actually can find two solutions depending if we start or not with a population
containing enough hard fields. The unphysical branch survives even when the gap [mr,ms] closes, see the example of
c = 45 in the figure.
We would like to stress at this point the enormous similarity of the soft-field part of the curve Σ(s) to the one in
fig. 4 in ref. [35]. Actually the variational results of [35] should be very precise and relevant near to the continuous
clustering transition (this is also case for the 3-coloring of Erdős-Rényi graphs or for 3-SAT).
C. Results for Erdős-Rényi random graphs
For Erdős-Rényi random graphs obtaining the solution of eq. (28) is computationally more involved as the solution is
no longer factorized. In the population dynamics a population of populations has to be updated, which is numerically
possible only for small populations, and so one has to be careful that the finite population-size corrections are small
enough, see details in appendix D. However, all the computations can be done with the same computational complexity
as for the regular graphs for m = 0, the energetic zero temperature limit (section IV C 1), and for m = 1 (appendix C).
That is enough to obtain the SP, clustering, condensation and COL/UNCOL transitions (from which the first and
last one was computed in [30]). We can also compute exactly when hard fields appear for m = 1, eq. (C9), this
transition is further studied in [64]. Finally, using the generalized survey propagation equation introduced in section
IV C 2, the rigidity transition can be computed quite precisely.
1. The general case for q > 3, discontinuous clustering transition
The phase transitions in q-coloring of random Erdős-Rényi graphs for q > 3 are qualitatively identical to those
discussed in the case of random regular graphs. We plot the results for the total entropy (number of solutions) and
complexity (number of clusters which dominate the entropy) in the 4− and 5− coloring in fig. 6.
At the clustering transition cd the complexity becomes discontinuously positive, the large RS cluster suddenly
splits in an exponential number of smaller ones. The total entropy Σ∗ + s∗ is given by the RS formula (14) up to the
condensation transition cc. At the condensation transition the complexity of the dominating clusters becomes zero,
the total entropy stot = s∗ < sRS is given by the point where Σ(s
∗) = 0. The function stot(c) is non-analytical at
the point cc, it has a discontinuity in the second derivative. At the coloring threshold cs all the clusters of solutions
disappear, note, however, that the total entropy of the last existing clusters is strictly positive (about a half of the
total entropy at the condensation transition). That means that the COL/UNCOL transition is not only sharp but
also discontinuous in terms of entropy of solutions. Note that the positive entropy has two contributions: the trivial
and smaller one coming from presence of leaves and other small subgraphs, and the nontrivial and more important
one connected with the fact, that the ground state entropy is positive, even in the uncolorable phase or for the random
regular graphs.
Finally we located the rigidity transition, when frozen variables appears in the dominating clusters. For 3 ≤ q ≤ 8
this transition appears in the condensed phase. As the number of colors grows it approaches the clustering transition.
All the four critical values cd, cr, cc and cs are summarized in table II, values of cSP and cr(m = 1) are given for
comparison.
 0.02
 0.04
 0.06
 0.08
 0.12
 0.14
 0.16
 0.18
 8.3  8.4  8.5  8.6  8.7  8.8  8.9  9
Σ*(c)
cd cc cs
Poissonian graphs, q=4
 0.02
 0.04
 0.06
 0.08
 0.12
 0.14
 0.16
 0.18
 12.8  12.9  13  13.1  13.2  13.3  13.4  13.5  13.6  13.7
Σ*(c)
cd cc cs
Poissonian graphs, q=5
FIG. 6: The 1RSB total entropy and complexity of the dominating clusters for 4- and 5-coloring of Erdős-Rényi random
graphs. The complexity jumps discontinuously at the clustering transition cd while the total entropy stays analytical. The
complexity disappears at the condensation transition cc causing a non-analyticity in the total entropy. Finally the total entropy
discontinuously disappears at the coloring threshold. Dashed is the RS entropy left for comparison.
2. The special case of 3−coloring, continuous clustering transition
The only case which is left to be discussed is the 3-coloring of Erdős-Rényi graphs. It is different from q > 3
because the replica symmetric solution is locally unstable in the colorable phase (see section III C 3). The extremality
condition underlying the RS assumption ceases to be valid because of the mechanism discussed in section III C 3, with
a divergence of the spin glass correlation length: the main difference with the previous cases is therefore that the
clustering transition is continuous and coincide with the condensation transition.
However, the phenomenology does not differ too much from the other cases: cRS = cd = cc = 4; the phase where
the entropy is dominated by exponential number of states is thus missing and the complexity corresponding to m = 1
is always negative (see fig. 7 left together with the dependence of the total entropy on the connectivity). Note that
the curves Σ(s) for the 3-coloring have been already studied in [38, 39] where the authors considered however only
the range of connectivities c = [4.42, 4.69] = [cSP, cs].
q cd cr cc cs cSP cr(m=1)
3 4 4.66(1) 4 4.687(2) 4.42(1) 4.911
4 8.353(3) 8.83(2) 8.46(1) 8.901(2) 8.09(1) 9.267
5 12.837(3) 13.55(2) 13.23(1) 13.669(2) 12.11(2) 14.036
6 17.645(5) 18.68(2) 18.44(1) 18.880(2) 16.42(2) 19.112
7 22.705(5) 24.16(2) 24.01(1) 24.455(5) 20.97(2) 24.435
8 27.95(5) 29.93(3) 29.90(1) 30.335(5) 25.71(2) 29.960
9 33.45(5) 35.658 36.08(5) 36.490(5) 30.62(2) 35.658
10 39.0(1) 41.508 42.50(5) 42.93(1) 35.69(3) 41.508
TABLE II: Critical connectivities cd (dynamical, clustering), cr (rigidity, rearrangments), cc (condensation, Kauzmann) and cs
(COL/UNCOL) for the phase transitions in the coloring problem on Erdős-Rényi graphs. The connectivities cSP (where the
first non trivial solution of SP appears) and cr(m=1) (where hard fields appear at m = 1) are also given. The error bars consist
of the numerical precision on evaluation of the critical connectivities by the population dynamics technique, details are given
in appendix D.
All the results derived for the 4-coloring of 4-c-bi-regular bipartite graphs are quantitatively valid also here. We are
thus not surprised by the fact that in interval c = [4, 4.42] the survey propagation algorithm gives us a trivial result:
simply the maximum of the curve Σ(s) does not exist yet there is no nontrivial solution at m = 0. Yet, the entropy
is dominated by finite number of largest clusters which do not contain hard fields. The two solutions (hard-field and
soft-field) join at a connectivity around 4.55. Finally at cr = 4.66 the hard fields arrive to the dominating states (and
in consequence to all others).
 0.05
 0.15
 0.25
 3.9  4  4.1  4.2  4.3  4.4  4.5  4.6  4.7  4.8
Σm=1(c)
cd=cc cs
Poissonian graphs, q=3
 0  0.2  0.4  0.6  0.8  1
(c-cc)/(cs-cc)
 0  0.2  0.4  0.6  0.8  1
(c-cc)/(cs-cc)
large q
FIG. 7: Left: The total entropy for 3-coloring of Erdős-Rényi random graphs. The dashed line is the replica symmetric (also
the annealed) entropy, left for comparison. The complexity at m = 1 is shown, it is negative for c > 4, however, for connectivity
near to four it is very near to zero. Right: The values of parameter m∗ (Σ(m∗) = 0) as a function of connectivity for q = 3, 4, 5
and in the large q limit. The connectivity c is rescaled as (c− cc)/(cs − cc). It is striking that for q > 3 the curves are so well
fitted by the large q limit one. We are even not able to see the difference due to the error bars which are roughly of the point
size.
3. The overlap structure
We now give some results about the overlap structure in the random coloring to elaborate the intuition about
clusters. First, consider marginal probabilities ψi,αsi within a cluster α. Note that due to the color symmetry there
exist another q! − 1 clusters different only in the permutation of colors. We define the intra-cluster overlap of two
solutions (averaged over states) as
〈(ψi,αsi )
2〉α . (41)
In the paramagnetic phase δ = 1/q, otherwise we have to compute it from the fixed point of equation (20). The
overlap between two solutions which lie in two clusters, which differ just by permutation π of colors is
δj = δ
j − 1
q − 1 +
q − j
q(q − 1) , (42)
where j is the number of fixed positions in the permutation π (in particular δq = δ, and δ1 = 1/q). In fig. 8 we show
the overlap structure for 3- and 4-coloring. The probabilities that two random solutions have one of the overlaps can
be computed from the Poisson-Dirichlet process described in appendix B, in fact this is not a self-averaging quantity
[29].
 0.25
 0.75
 4.8 4.687 4.4 4.2 4 3.8 3.6
connectivity
cd=cc cs
 0.25
 0.75
 9 8.9 8.8 8.6 8.46 8.35 8.2 8
connectivity
cd cc cs
FIG. 8: Left: Overlaps structure in 3-coloring of random graphs as a function of connectivity. The intra-cluster overlap (upper
curve) grows continuously from 1/3 at the clustering transition c = 4. In the figure from up there are δ = δ3, δ1 and δ0. Right:
Overlaps structure in 4-coloring of random graphs as a function of connectivity. The intra-cluster overlap (upper curve) jumps
discontinuously from 1/4 at the clustering transition c = 8.35. The probability that two random solutions belong to the same
cluster, however, is zero between the clustering and condensation transition [8.35, 8.46]. In the figure from up there are δ = δ4,
δ2, δ1 and δ0.
D. Large q Asymptotics
We give here the exact analytical large q expansion of the previous results. In the asymptotic computations the
regular and Erdős-Rényi ensembles are equivalent (the corrections are of smaller order in q that the orders we give).
We refer to the appendix E for the explicit derivation of the formulae.
At large q a first set of transitions arises for connectivities scaling as q log q:
cSP = cr(m = 0) = q [log q + log log q + 1 − log 2 + o(1)] , (43)
cr(m = 1) = q[log q + log log q + 1 + o(1)]. (44)
cSP was already computed in [31] and cr is the rigidity transition. The clustering transition has to appear before the
rigidity one cd < cr. For all the finite q cases we looked at, cd was between cSP and cr.
A second set of transitions arises for connectivities scaling as 2q log q:
cc = 2q log q − log q − 2 log 2 + o(1) , (45)
cs = 2q log q − log q − 1 + o(1) . (46)
The condensation thus appears very close the COL/UNCOL transition and both are very far from the clustering and
rigidity transitions (those are on a half way in the phase diagram).
We show also in appendix E that for connectivity c = 2q log q − log q + α, one has
2qs(m) ≃ 2m log 2 , (47)
2qΣ(m) ≃ 2m − 2 −m2m log 2 − α . (48)
Since the RS free energy is correct until cc, which differers just by constant from cs, that means that for all connec-
tivities bellow cc the number of solutions is correctly given by the replica symmetric entropy (14). Indeed, the value
s(m = 1) can be obtain by a large expansion of eq. (14).
In fig. 9 we plot the complexity of dominating clusters Σ∗ = Σ(m∗), the total entropy stot = Σ∗ + s∗, and the
physical value of m∗ as a function of connectivity c = 2q log q − log q+ α. Note that the properly scaled values of the
total number of solutions at cc and cs, and the values cc, cs themselves, are already very close to those at q = 3, 4, 5
(see figs. 6 and fig. 7 left). The closeness is particularly striking for the values m∗ for q = 4 and q = 5 (see fig. 7
right).
These formulae show that in the large q limit, near to the coloring threshold, it is the number of clusters which
change with connectivity (i.e. α), and not their internal entropy (size). In the leading order, adding constraints near
to the COL/UNCOL transition thus destroys clusters of solutions, but do not make them smaller: this is due to the
fact that these clusters are dominated by frozen variables so that adding a link kills them most of the time. We also
computed the entropy value at the condensation transition, and found s(m = 1) = log 2/q. The entropy of the last
cluster (exactly at the COL/UNCOL transition) is s = log 2/2q.
 0  0.5  1  1.5  2
α=-1.75
αg=-2log2
α=-1.2
αq=-1
-2 -1.8 -1.6 -1.4 -1.2 -1 -0.8
2q Σ*(α)
αg αq
FIG. 9: Analytical result for the large q asymptotics close to the COL/UNCOL transition for c = 2q log q − log q + α. Left:
(rescaled) complexity versus (rescaled) internal entropy for different connectivities. The condensation transition appears for
α = −2 log 2. The maximum of the complexity becomes zero at cs for α = −1. Right: Total entropy (s
tot), complexity (Σ∗)
and the parameter m∗ versus α. Notice how the values for the total number of solution are already very close to those for finite
low q in figs. 6, 7.
VI. ALGORITHMIC CONSEQUENCES
In this section we give some algorithmic consequences of our findings. First, we discuss the whitening procedure. We
then introduce a random walk algorithm adapted from the Walk-SAT strategy and study its performance. We show
in particular that the clustering/dynamical transition does not correspond to the onset of hardness in the problem and
argue that it is instead the rigidity transition. Finally, we discuss the performance of the belief propagation algorithm
in counting and finding solutions, and show that is works much better than previously anticipated.
A. The whitening procedure
The whitening procedure as introduced in [65] can distinguish between solutions which belong to a cluster containing
hard fields and those which do not. Generally whitening is equivalent to the warning propagation (version of belief
propagation which distinguish only if a field is hard or not). Warning propagation for coloring was derived in [30].
Let us call ui→j = (1, 0, 0, . . . , 0) the hard field in the direction of the first color, i.e. in absence of node j the node i
takes only the first color in all the colorings belonging to the cluster in consideration, and similarly for other colors.
Denote ui→j = (0, 0, 0, . . . , 0) if ψi→j is not frozen in the cluster, we say that the oriented edge i→ j is then “white”.
The update for u’s follows from (4)
ui→js = min
k∈i−j
uk→ir + δr,s
− min
k∈i−j
uk→ir
 . (49)
To see if a solution {si} belongs to a cluster with frozen variables or not we initialize warning propagation with
ui→js = δs,si , and update iteratively according to (49) until a fixed point is reached (the update every time converge,
because starting from a solution we are only adding white edges). In the fixed point or all edges are white, then the
solution {si} does not belong to a frozen cluster, or some of the edges stay colored (non-white), then the solution {si}
belongs to a frozen cluster. Note that in the K-SAT problem (but not in general), whitening is equivalent to a more
intuitive procedure, where the directed edged are not considered [66, 67].
We wish to offer here an explanation of a paradox observed in [66, 67]. The SP algorithm gives information
on the frozen variables in the most numerous clusters (m = 0). Yet, the solutions which are found by the standard
implementation (decimation and SP plus Walk-SAT) do not belong to clusters with frozen variables, since they always
give a trivial whitening result (all directed edges are white)[66, 67]. We suggest that the decimation strategy drives
the system towards a solution belonging to a large cluster, which does not contain frozen variables. In this case, it
is logical that the result of the whitening is trivial, as it is observed. We believe this is reason why no nontrivial
whitenings are observed so far in the study of the K-SAT problem on large graphs.
Note that beyond the rigidity transition this argument does not work anymore, since there all the clusters (for all
m such that Σ(m) > 0) contain frozen variables. More precisely, for q ≥ 9 we could in principle end up in soft-clusters
even beyond the rigidity transition (since that one concerns only the dominant states), if this is possible is let for
further investigation. Interestingly, in the coloring problem we have not been able to find solutions beyond the rigidity
transition even with survey propagation algorithm (compare cr with the performance of SP in [30]). Further, more
systematic, investigations have to be done about these issues, employing other strategies for the use of the survey
propagation equations (for example the reinforcement [74]).
B. A Walk-COL algorithm to color random graph
In this paper, we have computed the correct clustering transition cd for the random coloring problem. Beyond this
transition, Monte Carlo algorithms are proven not to reach equilibrium as their time of equilibration diverges [55, 62].
It was often claimed, or assumed, that this point corresponds to the onset of hardness of the problem. However,
the fact that the physical dynamics does not equilibrate just means that the complete set of solutions will not be
correctly sampled —indeed Monte-Carlo experiments clearly display slow relaxation [75]— but not that no solutions
can be eventually found. This simple fact explains the results of [32] where a simple annealing procedure was shown
to 3-color a ER graph beyond cd = 4.
In this section, we use a local search strategy which does not satisfy the detailed balance condition. Therefore, we
do not expect to be able to find typical solutions, however it might be possible to find some solutions to the problem.
The Walk-COL algorithm [87] is a simple adaptation of the celebrated Walk-SAT [76]. More precisely, we adapted
the method designed for satisfiability in [77]. Given a graph, and starting from an initial random configuration, we
recursively apply the following procedure:
1) Choose at random a spin which is not satisfied (i.e. at least one of its neighbors has the same color).
2) Change randomly its color. Accept this change with probability one if the number of unsatisfied spin has been
lowered, otherwise accept it with probability p.
3) If there still are unsatisfied nodes, go to step 1) unless the maximum running time is reached
The probability p has to be tuned in each different case for a better efficiency of the algorithm. Typically, values
between 0.01− 0.05 give good results. We shall now briefly discuss the performance of the algorithm, to illustrate the
two following points: (a) When the phase space is RS, we observe that Walk-COL finds a solution in linear time. (b)
Even in the “complex” phase for c > cd, the algorithm can find in some cases solutions in linear time.
Concerning the first point, we tested the algorithm in the RS phase of regular random graphs for q = 3, 4, 5, 6, 7.
In all these cases, we were able to color in linear time all the graphs of connectivities that correspond to a replica
symmetric solution. In particular, the cases (q = 3, c = 5), (q = 5, c = 13), (q = 6, c = 17), (q = 7, c = 21),
(q = 7, c = 22) are found to be colorable with the Walk-COL algorithm even if a nontrivial solution to the SP
equations exists.
Concerning the second point, we considered the 3− and 4−coloring of Erdős-Rényi random graphs. The results
are shown in fig. 10 where the percentage of unsatisfied spins versus the number of attempted flips (averaged over
5 different realizations) divided by N is plotted. We observe that the curves corresponding to different values of N
superpose quite well (and that actually the results for N = 2 · 105 are systematically lower than those for N = 5 · 104)
so that an estimation of the time needed to color a graph can be obtained. The connectivities of these graphs are
beyond the dynamical transition (cd = 4 for 3-coloring and cd = 8.35 for 4-coloring). It would be interesting to
systematically test Walk-COL, as it has been done for Walk-SAT in [77], to derive the precise connectivity at which
it ceases to be linear.
 1e-05
 1e-04
 0.001
 0.01
 100  1000  10000  100000  1e+06
c=4.1
c=4.3
c=4.4 c=4.5
3-coloring
N=50 000
N=200 000
 1e-05
 1e-04
 0.001
 0.01
 100  1000  10000  100000  1e+06
c=8.0
c=8.3
c=8.4
c=8.54-coloring
N=50 000
N=200 000
FIG. 10: Performance of the Walk-COL algorithm in coloring random graphs for 3−coloring (left) and 4−coloring (right). We
plot the rescaled time (averaged over 5 instances) needed to color a graph of connectivity c. The strategy allows one to go
beyond the clustering transition (cd = 4 for 3-coloring and cd = 8.35 for 4-coloring) in linear time with respect to the size of
the graph.
Already these results show that the dynamical transition is not a problem for the algorithms. This can also be
observed in a number of numerical experiments for the satisfiability [77, 78] and the coloring [73, 79] problems.
We believe, however, that the rigidity transition plays a fundamental role for the average computational complexity.
A first argument for this is that, for large graphs, it seems that all the known algorithms are only able to find solutions
with a trivial whitening, i.e. solutions that belong to clusters without hard fields. Beyond the rigidity transition
however, the clusters without hard fields become very rare (in the sense that the dominating clusters and all the
smaller, more numerous ones, contain hard fields). For q ≥ 9 maybe the connectivity where hard fields appear in
clusters corresponding to Σ(m) = 0, m > 1 should be considered. This suggests that the known algorithms will not
be able to find a solution beyond this point.
A second argument is the following: local search algorithms are either attracted into a solution or stucked in a
metastable state. These metastable states, in order to be able to trap the dynamics, have to contain a finite fraction of
hard fields. Given an algorithm, determining which of these two situations happens is not only a question of existence
of states, but also a question of basins of attraction and a theoretical analysis of such basins is a very difficult task
so that the precise analysis of the behavior of local algorithms remains a hard problem. However, the metastable
states are known, from the cavity formalism, to be much more numerous than the zero-energy states. Moreover the
basin of attraction of a zero-energy state that contains hard fields does not probably differ much from those of the
metastable state (while, on the other hand, the basin of attraction of a zero-energy state which does not contain hard
fields might be slightly different and arguably relatively larger). It thus seems to us reasonable that local algorithms
will get trapped by the metastable states beyond the rigidity transition.
A similar conclusion was reached recently in [73] where the recursive implementation of the Walk-COL algorithm
 4  4.1  4.2  4.3  4.4  4.5  4.6  4.7  4.8
cd=cc cr cs
N=2000
N=4000
N=8000
N=16000
N=32000
 8  8.2  8.4  8.6  8.8  9
cd cc cr csN=2000
N=4000
N=8000
N=16000
N=32000
FIG. 11: Performance of the BP algorithm plus decimation in coloring random graphs for 3−coloring (left) and 4−coloring
(right). The strategy described in the text allows to color random graphs beyond the clustering and even the condensation
transitions.
was studied and found to be somehow simpler to analyze. Again the strategy was found to be efficient (with linear
time with respect to the size of the graph) beyond cd but bellow cr. The precise algorithmic implications of the rigidity
transition thus require further investigation, maybe in the lines of [73, 78].
C. A belief propagation algorithm to color random graph
Another consequence of our results, that we already discussed shortly in [42], is that the standard belief propagation
(BP) algorithm gives correct marginals until the condensation transition. It is actually a simple algebraic fact that the
1RSB approach at m = 1 gives the same results for the marginals (average probabilities over all clusters) as the simple
RS approach (see appendix C). Moreover the log-number of solutions in clusters corresponding to m = 1 is also equal
to the RS one. This suggests to use the BP marginals (as was already suggested in [40]) and a decimation procedure
to find proper colorings. Compared with the SP algorithm which has computational complexity proportional to q!
(factorial) the BP is only q. We have seen moreover that for large numbers of colors, the condensation point is very
close to the COL/UNCOL transition, so that BP could be used in a large interval of connectivities.
As a simple application, we tested how the straightforward implementation of the BP algorithm plus decimation
allows one to find solutions of 3− and 4−coloring of random Erdős-Rényi graphs. Note that in the 3−coloring the
clustered but not condensed phase is missing, so the argumentation above does not concern this case. The algorithm
works by iterating the following procedure:
(i) Run BP on the graph for a given number l of iterations.
(ii) Consider the most biased variable, and color it with its most probable color.
Two problems have to be mentioned. The first one is rather trivial: since at the beginning all colors are symmetric,
the first color had to be put at random. The second one is more serious and concerns the convergence of BP. Indeed,
we saw that there is local instability in the BP (replica symmetric) equations at connectivity c = 4 for the 3−coloring
of random graphs, so that the BP equations do not converge. This seems to be a problem restricted to the 3−coloring,
but even in the case of 4− or more coloring, the BP equations do not converge on the decimated graph when a finite
fraction (typically few percent) of variables is fixed. The reason or that is to be understood.
Nevertheless, since we merely want to design an effective tool to solve the coloring problem, we choose to avoid this
problem by fixing the number of iteration l at each step and thus ignore the non-convergence. We tried the method
on both the 3− and 4−coloring and obtained unexpectedly good results. We used the following protocol in the code:
We first try to find a solution with l = 10. If we do not succeed, we restart with l = 20 and once more with l = 40.
We tried that on 10 different samples for different connectivities. The probability to find a proper coloring with these
conditions is shown in fig.11.
We thus observe that the BP strategy is able to find solutions, even beyond the condensation transition. This shows
clearly that the decimation procedure is a nontrivial one, and that the problem is not really hard in that region of
connectivities. Note that the SP algorithm plus decimation has been shown to work in the 3−coloring very well until
about 4.60 [30]: our results are thus very close to those obtained using SP. This rises again the question of the rigidity
transition cr = 4.66, which might also be problematic for the decimated survey propagation solver.
VII. CONCLUSIONS
Let us summarize the results. They are perhaps best illustrated looking back to the cartoon in fig. 1, where the
importance of the size of clusters is evidenced. We find that the set of solutions of the q-coloring problem undergoes
the following transitions as the connectivity is increased:
(i) At low connectivity, c < cd, many clusters might exist but they are very small and the measure over the set of
solutions is dominated by the single giant cluster described by the replica symmetric approach.
(ii) Only at the dynamical transition cd the giant cluster decomposes abruptly into an exponentially large number
of clusters (pure states). For connectivities cd < c < cc, the measure is dominated by an exponential number of
clusters. Yet, the total number of solutions is given by the replica symmetric entropy (14), and the marginals
are given by the fixed point of the replica symmetric equations (belief propagation) (4). Starting from this
transition the uniform sampling of solution becomes hard.
(iii) At connectivity cc the set of solutions undergoes a condensation transition, similar to the one appearing in mean
field spin glasses. In the condensed phase the measure is dominated by finite number of the largest clusters. The
total entropy is strictly smaller than the replica symmetric one and has a discontinuity in the second derivative
at cc.
(iv) When connectivity cs is reached, no more clusters exist: this is the COL/UNCOL transition. Note that the
entropy of last existing clusters is strictly positive, and not given only by the contribution of the isolated nodes,
leaves and other small subgraphs, the COL/UNCOL transition is thus discontinuous in entropy.
This picture is very similar to the well-known scenario of the glass transition in temperature, with the dynamical
and glass (Kauzmann) transition [61]. In some cases, the main one being the 3-coloring of Erdős-Rényi graphs,
the clustering and the condensation transition merge and a continuous transition take place at cd = cc, which then
coincide with the local instability of the replica symmetric solution. Interestingly the variational approach of [35] is
very precise near to the continuous clustering transition. Since the 3-SAT problem behaves in the same manner, this
solves the apparent contradictions between the results of [35] and [24].
In addition to the transitions describing the organization of clusters, another important phenomenon concerning
the internal structure of clusters takes place. A finite fraction of frozen variables can appear in the clusters (a frozen
variable takes the same color in all the solutions that belong to the cluster). We found that the fraction of such
variables in each cluster undergoes a first order transition and jumps from zero to a finite fraction at a connectivity
that depends on the size of the cluster. In particular:
(v) There exists a critical connectivity cr (rigidity/freezing) at which the thermodynamically relevant clusters —
those that dominate the Gibbs measure— start to contain a finite fraction of frozen variables.
The results above were obtained within the 1RSB scheme, but should not change when considering further steps
of RSB (an exception might be the 3-coloring near to the clustering transition).
We discussed some algorithmic consequences of these transitions. First, the belief propagation algorithm is efficient
in counting solution and estimating marginals until the condensation transition. More interestingly, it can also be
used, just like survey propagation, together with a decimation procedure in order to find solutions as we numerically
demonstrated. Secondly, the dynamical transition is not the one at which simple algorithms fail as we illustrated using
the Walk-COL strategy. For the 3-coloring of ER graph, there is even a rigorous proof of algorithmic performance
beyond cd = 4 and until c = 4.03 [10]. We argued that, instead, the rigidity phenomenon is responsible for the onset
of computational hardness. This is a major point that we hope to see more investigated in the future.
Our study opens a way to many new and promising investigations and developments. For instance, we wrote the
equivalent of the survey propagation equations for general value of m, which has a particularly simple form for m = 1
(39). It would be interesting to use these equations to find solutions. The behavior at finite temperature and the
performance of the annealing procedure are also of interest. It would furthermore be interesting to re-discuss other
finite connectivity spin glass models like for instance the lattice glass models [56] in the light of our findings. The
stability towards more steps of replica symmetry breaking, or the super-symmetric approach [71], should be further
investigated. Finally, it would be interesting to combine the entropic and energetic approach to investigate the frozen
variables in the meta-stable states. We hope that our results will stimulate the activity in these lines of thoughts.
Acknowledgments
We thank Jorge Kurchan, Marc Mézard, Andrea Montanari, Federico Ricci-Tersenghi, Guilhem Semerjian and
Riccardo Zecchina for cheerful and very fruitful discussions concerning these issues. The numerical computations
were done on the cluster EVERGROW (EU consortium FP6 IST) at LPTMS, Orsay, and on the cluster DALTON
at ESPCI, Paris. This work has been partially supported by EVERGROW (EU consortium FP6 IST).
APPENDIX A: STABILITY OF THE PARAMAGNETIC SOLUTION
In this appendix, we show how to compute the stability of the paramagnetic solution towards the continuous
appearance of a 1RSB solution. This happens, as usual for continuous transition, when the spin glass correlation
length, or equivalently, the spin glass susceptibility, diverges. Obviously, the presence of the diverging correlation
length invalidate the premise of the RS cavity method. Recall that the spin glass susceptibility is defined as
χSG =
〈sisj〉2c . (A1)
and can be rewritten for the present purpose as
χSG =
γdE(〈s0sd〉2c) , (A2)
where we consider the average over graphs, in the thermodynamic limit, where spins s0 and sd are at distance d. The
factor γd stands for the average number of neighbors at distance d, when d≪ logN . Assuming that the limit for large
d of the summands in (A2) exists (with the limit N → ∞ performed first), we relate it to the stability parameter :
λ = lim
E(〈s0sd〉2c)
. (A3)
Then the series in (A2) is essentially geometric, and converges if and only if λ < 1.
Using the fluctuation-dissipation theorem we relate the correlation 〈s0sd〉c to the variation of magnetization in s0,
caused by an infinitesimal field in sd. Finally, using the fact that we perform the large-N limit first, the variation
above is dominated by the direct influence through the length-d path between the two nodes, and this induces a
“chain” relation: if the path involves the nodes (d, d− 1, . . . , 0) we have
E(〈s0sd〉2c) = C ·
∂ψa→0
∂ψb→d
= C · E
∂ψl→l−1
∂ψl+1→l
. (A4)
The stability parameter of the paramagnetic solution of the cavity equations towards small perturbations can be
computed from the following Jacobian
T τσ =
∂ψ1→0τ
∂ψ2→1σ
, (A5)
which gives the infinitesimal probability that a change in the input probability ψ2→1σ will change the output probability
ψ1→0τ . The index RS says that the expression has to be evaluated at the RS paramagnetic solution.
This matrix has only two different entries, all the diagonal elements are equal, and all the non-diagonal elements
are also equal. As an immediate consequence all Jacobians commute and are thus simultaneously diagonalizable so
that it will be sufficient to study the effect after one cavity iteration (one step in the chain). The matrix T has only
two distinct eigenvalues,
∂ψ2→11
∂ψ2→12
∂ψ2→11
+ (q − 1) ∂ψ
∂ψ2→12
The second eigenvalue corresponds to the homogeneous eigenvector (1, 1, ..., 1) and describes a fluctuation changing
all ψ2→1τ , τ = 1, ..., q, by the same amount, and maintains the color symmetry. It is thus not likely to be the relevant
one and we will see that indeed λ2 = 0. The first eigenvalue, however, is (q − 1)-fold degenerate and its eigenvectors
are spanned by (1,−1, 0, ..., 0), (0, 1,−1, 0, ..., 0), . . . , (0, ..., 0, 1,−1). The corresponding fluctuations explicitly break
the color symmetry, and are in fact the critical ones. Using the cavity recursion (4), the two derivatives simply read
∂ψ2→12
= (1 − e−β) (ψ
1−(1−e−β)ψ21
∂ψ2→11
= (1 − e−β)
(ψ1→01 )
1−(1−e−β)ψ2→11
1−(1−e−β)ψ2→11
so that the values of the two eigenvalues evaluated at the RS solution, where all ψ are equal to 1/q, are
1 − q
1−e−β
, λ2 = 0. (A8)
The stability parameter (A3) is thus λ = γλ21 and the critical temperatures bellow which the instability sets in are
T regc (q, c) = −1/log
1 − q√
c− 1 + 1
, TERc (q, c) = −1/log
1 − q√
. (A9)
For regular and Erdős-Rényi graphs respectively. Thus at zero temperature the critical connectivities reads
RS stab = q
2 − 2q + 2 , cERRS stab = q2 − 2q + 1. (A10)
These results coincide perfectly with the numerical simulations of the cavity recursion of [32]. The analytical ex-
pressions equivalent to (A9) were in fact first obtained in [57] in the context of the reconstruction problem on trees
as an upper bound for the Gibbs measure extremality, and its connection with the statistical physics approach was
explained in [53]. The case of bi-regular random graphs can be easily understood by noticing that two recursions
should be considered, one with γ = c1−1, and one with γ = c2−1. As a consequence, the stability point is equivalent
in this case to the one of a regular random graph with an effective connectivity equal to c = 1 +
(c1 − 1) (c2 − 1).
Another instability appears when γ|λ1| > 1. This has been refered to as the modulation instability in [56]. Actually,
this is the continuous instability towards the appearance of the anti-ferromagnetic order. Since at zero temperature
λ1 = (1 − q)−1, then for connectivities larger than cmod = q for random regular graphs (and cmod = q − 1 for Erdős-
Rényi) the paramagnetic solution becomes unstable towards the anti-ferromagnetic order. However, this is correct if
we study a tree with some given (and well chosen) boundary condition, but as noted in [56], the anti-ferromagnetic
solution in impossible on random graphs because of the existence of frustrating loops of arbitrary length. The cavity
equations (4) can actually never converge towards an anti-ferromagnetic solution of a random graph. Instead, when
iterating, the fields oscillate between different solutions (thus the name modulation).
In other words, although on a random tree with special boundary conditions there exists for c > cmod a nontrivial
solution to the cavity recursion (for the Gibbs state is no longer unique (15)), this solution does not exist on a
random graph (and the Gibbs state is still extremal (16)). Note that this instability could anyway be a source of
numerical problems that can be overcome considering that the distribution of cavity fields P(ψ) over the ensemble of
random graphs has to be symmetric in the color permutation. Another possibility is to randomly mix the new and
old populations in the population dynamics so that the anti-ferromagnetic oscilations are destroyed.
APPENDIX B: THE RELATIVE SIZES OF CLUSTERS IN THE CONDENSED PHASE
In this section, we introduce the Poisson-Dirichlet point process and we shortly review some of its important
properties. We also sketch its deep connection with the size of clusters in the condensed phase. Poisson-Dirichlet
(PD) point process is a set of points {xi}, i = 1, . . . ,∞ such that x1 > x2 > x3 > . . . and
i=1 xi = 1. To construct
these points we consider a Poisson process {yi}, i = 1, . . . ,∞ of intensity measure y−1−m
, 0 < m∗ < 1 (note that
this measure is not a probability measure). We order the sequence {yi} in such a way that y1 > y2 > y3 > . . . and
define the PD point process as
i=1 yi
. (B1)
If we identify the parameter m∗ with the value for which the complexity is zero Σ(m∗) = 0 then yi is proportional
to the number of solutions in cluster i (or to e−βF for non-zero temperature), and xi is the size on that cluster
relative to the total number of solutions. This connection was (on a non-rigorous level) understood in [80], for more
mathematical review see [81]. Note that that due to the permutation symmetry in graph coloring there are every
time q! copies of one clusters (different in the color permutation).
To get feeling about the PD statistics let us answer in fig. 12 to the following question: Given the value m∗ how
many clusters do we need to cover fraction r of solutions, in other words what is the smallest k such that
i=1 xi > r?
The mathematical properties of the PD process are very clearly reviewed in [82]. To avoid confusion, note at this
point that the PD process we are interested in is the PD(m∗, 0) in the notation of [82]. In the mathematical literature,
it is often referred to the PD(0, θ) without indexing by the two parameters. Let us remind two useful results. Any
moment of any xi can be computed from the generating function
E[exp (−λ/xi)] = e−λφm∗(λ)i−1ψm∗(λ)−i , (B2)
where λ ≥ 0 and the functions φm∗ and ψm∗ are defined as
φm∗(λ) = m
e−λxx−1−m
dx , (B3)
ψm∗(λ) = 1 +m
(1 − e−λx)x−1−m
dx . (B4)
Another relation is that the ratio of two consequent points Ri = xi+1/xi, i = 1, 2, . . . is distributed as im
In particular its expectation is E[Ri] = im
∗/(1 + im∗) and the random variables Ri are mutually independent. We
used these relation to obtain data in figure 12.
 0 0.2 0.4 0.6 0.8 1
 0 0.2 0.4 0.6 0.8 1
FIG. 12: The sketch of size of the largest clusters for given value of parameter m∗. The lower curve is related to the average
size of the largest clusters as 1/E[1/x1] = 1−m
∗. The following curves are related to the size of i largest clusters, their distances
are E[Ri]E[Ri−1] . . .E[R1](1 −m
APPENDIX C: THE 1RSB FORMALISM AT m = 1 AND THE RECONSTRUCTION EQUATIONS
In this appendix we discuss the considerable simplification of eq. (28) that is obtained by working directly at m = 1.
This was first remarked and proved in [53] when dealing with the tree reconstruction problem (for a discussion of a
case where the RS solution is not paramagnetic see [43]). We first introduce the probability distribution of fields (20)
averaged over the graph
P (ψ) ≡
dP P [P (ψ)]P (ψ) =
Q1(k)
dψi P
(ψi) δ
ψ −F({ψi})
Z0, (C1)
where Z1 is computed from (20) as
dψi P i(ψi)Z0 =
1 − ψi
, (C2)
where ψ =
dψ P (ψ)ψ. Generally, ψ is a solution of the RS equation (12), which is easily seen from (20), (28). Since
the RS solution is the paramagnetic one ψ = 1/q the form of Z1 is particularly simple.
In the next step, we want to get rid of the term Z0 in eq. (C1). We thus introduce q distribution functions P s
P s(ψ) = qψsP (ψ) . (C3)
It is then easy to show that if ψ = 1/q then P s(ψ) satisfies
P s(ψ) =
Q1(k)
s1...sk
π(si|s)
ψ −F({ψi})
dψiP si(ψ
i) , (C4)
where
π(si|s) =
1 − (1 − e−β)δ(si, s)
q − (1 − e−β) . (C5)
We solve eq. (C4) by population dynamics. In order to do this, one needs to deal with q populations of q-component
fields, and to update them according to (C4). Is is only a functional equation and not a double-functional as the
general 1RSB equation (28). Moreover the absence of the reweighting term Zm0 simplifies the population dynamics
algorithm significantly. Finally, it is important to note that the computational complexity here is the same as the one
for regular and ER random graphs.
A crucial theorem is also proven in [53]: the population dynamics of eq. (C4) has a nontrivial solution if and only
if it converge to a nontrivial solution starting from initial conditions:
r (ψs) = δ(r, s). (C6)
This shows that when a paramagnetic solution is found, then no other solutions exist.
Similar manipulations allow us to obtain the replicated free energy (21) which is equal in this case to the replica
symmetric free energy (9), and the free energy (26) inside the corresponding states as
− βf(β) =
s1...sk
π(si|s)
logZi0
dψiP si(ψ
i) , (C7)
s1,s2
π(s1|s2)
logZ120
dψ1dψ2 P s1(ψ
1)P s2(ψ
where the normalization factors Zi0, Z
0 are defined by (6) and (8). The complexity follows from (24). Since the
replicated free energy Φ(β, 1) is equal, according to (27), to the total free energy, we showed the statement used
several time in the paper, i.e. the total free energy (entropy) at m = 1 is equal to the replica symmetric free energy
Another important point is that one can write the recursion separating the hard and soft fields. In general, at zero
temperature, we can write the distribution P q(ψs) in eq. (C4) as
P r(ψ) =
µr,sδ(ψs − 1) + (1 −
µr,s)P̃r(ψ). (C8)
Plugging this to eq. (C4) and taking into account the initial condition (C6) and color symmetry, we see that µq,s =
qηδ(q, s), where η satisfies
Q1(k)
(−1)m
q − 1
1 − mq
q − 1η
. (C9)
On ER graphs, the sum can be performed analytically and one finds
1 − e−
. (C10)
This equation can be solved iteratively starting from ηinit = 1/q. It is a very simple equation, as the one obtained
for m = 0, which gives us a very efficient way to compute the fraction of hard fields at m = 1 for both regular and
ER graphs. Indeed η is larger that zero only above a certain average connectivity cr(m = 1).
APPENDIX D: NUMERICAL METHODS
In this section, we detail the numerical methods we used to solve the 1RSB equations (20,28), and the procedures
used to generate the data. We use a population dynamics method, as introduced in [26, 27], and model the distribution
P i→j(ψi→j) by a population of N vectors ψi→j . To compute P i→j(ψi→j) knowing the P k→i(ψk→i) for all incoming
k we perform the 1RSB recursion in eq. (20) in two steps: (i) first we compute the new vectors ψi→j using the simple
RS recursion in eq. (31) (this is the iterative step) and (ii) we take into account the weight (Z
m for each of the
vectors (this is the reweighting step). For the reweighting we tried different strategies, two of them perform very well.
a) For every field ψ in the population, we keep its weight Z0. We then compute the cumulative distribution of
weights Z0 and sample uniformly the incoming fields. Using dichotomy we generate a random fields with its
proper weight in O(log(N)) steps. A complete iteration thus takes O(N logN) steps.
b) We compute N new vectors and then we make a new population when we clone some of them while erasing
others so that in this new population each field is present according to its weight (in principle, one can even
change the size of the population, although we have not implemented this strategy). This second approach can
be implemented in linear time (generating an ordered list of random numbers is a linear problem, see [83]), but
is a bit less precise as we introduce redundancy in the population.
We finally choose to use the second strategy, as we observed that it performs almost as good at the first one (for a
given size of population) while it was much faster, so that, for given computer time, it allows a better representation
of the population. We also force the population to be color-symmetric by adding a random shift of colors in the
incoming messages. This is needed in order to avoid the anti-ferromagnetic solution. The learned reader will notice
that this is equivalent to solving a disordered Potts glass instead of a anti-ferromagnet model. Indeed the fact that
an Ising anti-ferromagnet on a random graph is equivalent to an Ising spin glass was already noticed [84].
Another important issue is the presence of hard fields. In fig. 13 are histograms of the first component of the vectors
in the population for 3- and 4-coloring of 5- and 9-regular random graphs respectively. It is interesting to see how
they peak around fractional values due to the presence of hard fields (see the three upper one). Maybe even more
interesting are the lower one where no hard fields are present. However, since there are soft fields with values 1 − ǫ,
where ǫ can be almost arbitrary small, one cannot see from these picture the absence of frozen variables. For the case
c = 9, q = 4,m = 0.8 for instance, the presence of the quasi-hard fields makes the distribution clearly concentrate on
values around one, zero and half (note however that the amplitude —on a logarithmic scale— is far less important).
 1e-05
 1e-04
 0.001
 0.01
 0  0.2  0.4  0.6  0.8  1
c=5, q=3, m=0
 1e-05
 1e-04
 0.001
 0.01
 0  0.2  0.4  0.6  0.8  1
c=5, q=3, m=2
 1e-05
 1e-04
 0.001
 0.01
 0  0.2  0.4  0.6  0.8  1
c=9, q=4, m=0
 1e-05
 1e-04
 0.001
 0.01
 0  0.2  0.4  0.6  0.8  1
c=9, q=4, m=0.4
 1e-05
 1e-04
 0.001
 0.01
 0  0.2  0.4  0.6  0.8  1
c=9, q=4, m=0.8
 1e-05
 1e-04
 0.001
 0.01
 0  0.2  0.4  0.6  0.8  1
c=9, q=4, m=1.6
FIG. 13: Histograms of the first component of the cavity field ψ1, i.e. the probability that a node takes color one. Notice
the logarithmic y-scale. Frozen fields (ψ1 = 1) are present in the solution for the three upper cases; there are delta-peaks on
0,1,1/2 and other simple fractions depending on q. Notice that even when frozen fields are not present, there are many almost
frozen fields (the distributions only concentrate around 0,1,1/2 and other fractions).
The quasi-hard fields are therefore very hard to distinguish numerically from the true hard ones. This is evidenced
on fig. 14 where we plot the fraction of hard fields computed using the expression (38) together with a numerical
estimate made by population dynamics without the separate hard/soft implementation. We show that the fraction
of fields of value 1 − ǫ is not zero in regions where we know that there are no hard fields even for ǫ = 10−20. This
demonstrates the presence of quasi-hard fields, with ǫ going to zero as the critical m is approached. This transition
is further studied in [64].
 0  0.1  0.2  0.3  0.4  0.5  0.6  0.7  0.8  0.9
c=9, q=4
mrε=10-4
ε=10-8
ε=10-12
ε=10-16
ε=10-20
FIG. 14: Fraction of the hard and the quasi-hard cavity fields qη (a field is quasi-hard if ψ > 1−ǫ) in the 4−coloring of 9−regular
graph. The bold line is obtained with the analytical computation of the fraction of hard fields and the dot corresponds to the
threshold mr.
An important simplification of the 1RSB update (20) arises when we consider the soft and the hard fields separately.
The fraction of hard fields can be computed using the generalized SP equation (38), provided the ratio Zmsoft/Z
hard is
computed. This considerably reduced the size of the population as only soft field has to be kept in memory. Another
way to further speed up the code is to generate directly soft fields with a uniform measure instead of waiting for
them to come out from the 1RSB recursion. Indeed they might be quite rare in the region of small m and one can
spend a considerable amount of time before being able to sample them correctly. Generating the soft fields with a
uniform weight turns out to be rather easy using the following method: (i) Choose two random colors q1 and q2. (ii)
Perform the usual recursion (31) in order to have a new vector but forbid incoming hard fields to q1 and q2. (iii) To
obtain a uniform soft field generator, the resulting field should be weighted by 1/
where s is the number of non-null
components in the vector. This is specially useful in the case of Erdős-Rényi graphs.
The formula for the free entropy Φs(m) (32) also simplifies in this case. Consider a given site i; the site free entropy
term can be split into three parts when (i) the field is hard, (ii) the total field is soft and (iii) the field is contradictory.
Φis = log
phard(Zhard)m + psoft(Zsoft)m
, (D1)
where phard/psoft are the probabilities that the total field is frozen/soft, and are given by the SP recursion. Indeed
the probability that the total field is not contradictory (phard + psoft) is the denominator in eq. (35) while phard is the
numerator of eq. (35). The link part can also be simplified using the fact that contradictions arise when two incoming
frozen messages of the same color are chosen, so that
Φijs = log
pno contr(Zno contr)m
, (D2)
where pno contr is simply (1 − qηi→jηj→i).
For m = 0 the formula further simplify as Zmhard = Z
soft = Z
no contr = 1 so that
Φis = log
(−1)l
1 − (l + 1)ηk→i
, (D3)
Φijs = log
1 − qηi→jηj→i
. (D4)
This is precisely what was obtained within the energetic cavity approach in [30]. The numerical population dynamics
implementation with mixed hard/soft strategy is therefore as precise as it could be since we obtain the exact evaluation
in the m = 0 case. This simple computation also demonstrates how one can recover the energetic zero temperature
limit from the generic formalism.
Finally, we obtain the function Φs(m). We fit this function using an ansatz Φs(m) = a + b2
m + c3m . . . and then
perform the Legendre transform to obtain the entropy and complexity. It is also possible to compute directly the
complexity from the population data using the expression of the derivative of the potential directly in the code. Both
methods lead to very good results. We show an example of the raw data and their fit in fig. D, where the data have
been obtained with relatively small population (N = 5000) but where the mixed strategy separating the hard and
soft fields have been used. For the purely soft-field branch, we used N = 50000. It took few hours up to few days to
generate these curves on present Intel PCs.
-0.02
 0.02
 0.04
 0.06
 0.08
 0.12
 0.14
 0.16
-1 -0.5  0  0.5  1  1.5  2
Φ(x) for q=6, c=19
-0.08
-0.06
-0.04
-0.02
 0.02
 0  0.05  0.1
Σ(s) for q=6, c=19
 0.05
 0.15
 0.25
-2 -1  0  1  2
Φ(x) for q=4, c=9
-0.12
-0.08
-0.06
-0.04
-0.02
 0.02
 0.04
 0  0.05  0.1  0.15
Σ(s) for q=4, c=9
FIG. 15: The numerical results for the free entropy (30) and its fit with a function a+b2x +c3x... for the 6-coloring of 19-regular
graphs and the 4-coloring of 9-regular graphs. Circles give the analytical results at m = 0 and m = 1. On the right parts,
we present the complexity versus internal entropy with the numerical points and the Legendre transform of the fit of the free
entropy. The analytical result for Σmax is also shown.
In the case of bi-regular random graphs, one needs two different populations: one for the fields going from nodes
with connectivity c1 and one for the fields going from nodes with connectivity c2. Then each iteration for population
1 (resp. 2) should be performed using as incoming messages the vectors of population 2 (resp. 1.). Again, one can
separately perform the recursion for the hard-field fractions in both population.
The case of Erdős-Rényi random graphs is more involved, as one needs a large number Npop of populations, each
of them of size N . In this case, using the separate hard/soft fields implementation and the formulae (D1,D2) for
complexity is crucial, as it allows a good precision even for smaller population sizes. We used typically 2Npop/c ≈
(1 − 3) · 103 and N ≈ (1 − 3) · 102. The error bars in table II are computed from several independent runs of the
population dynamics. In each case we were able to make the equilibration times and the population sizes large enough
such that by doubling the time or the population size we did not observed any significant systematic changes in the
average results.
APPENDIX E: HIGH-q ASYMPTOTICS
The quenched averages in the large q limit are the same for the regular and Erdős-Rényi graphs and we thus consider
directly the regular ensemble of connectivity c = k + 1. The appearance of a nontrivial 1RSB solution for m = 0,
which correspond to cSP, was already computed in [31] and reads
cr(m = 0) = cSP = q [log q + log log q + 1 − log 2 + o(1)] , (E1)
ηd(m = 0) =
1 − 1
log q
log q
, (E2)
while the coloring threshold is [31]
cs = 2q log q − log q − 1 + o(1). (E3)
We now show how the connectivity where a solution with hard fields at m = 1 first appears, and how the complete
free entropy Φs(m) (30) can be computed close to the COL/UNCOL transition
1. The appearance of hard fields at m = 1
We first show that the correct scaling for the appearance of hard fields at m = 1 is
k = q[log q + log log q + α]. (E4)
and compute the value of α. In the order O(q) we can write also k = (q − 1)[log(q − 1) + log log(q − 1) + α].
The starting point is the equation (C9), with Q1(x) = δ(x − k). In the large q limit the fraction of hard fields is
µ(q, k) = qη(q, k) = 1− θ(q, k), where θ(q, k) = o(1) is the fraction of soft fields. We check self-consistently at the end
of the computation that only the two first terms of (C9) are important. Then we have
µ(q, k) = 1 − (q − 1)
1 − 1
q − 1µ(q, k)
= 1 − (q − 1)e−
kµ(q,k)
q−1 . (E5)
A self-consistent equation for θ(q, k) follows
log(q − 1)θ(q, k) = (q − 1)θ(q,k)e−α , (E6)
which is solved by θ(q, k) = γ(α)/ log(q − 1) where
γ(α)e−γ(α) = e−α . (E7)
The maximum of the left hand side is 1/e for γ(α) = 1. It means that a solution of (E7) exists for α > 1. Finally the
hard fields appear in the 1RSB solution for m = 1 at connectivity
cr(m = 1) = q[log q + log log q + 1 + o(1)]. (E8)
The clustering transition cd should be between cSP and cr(m = 1) as this is what we observed for finite q. We see
that cSP and cr(m = 1) differs only in the third order and both are very far from the coloring threshold and also from
the condensation transition as we show in section E 2. It would be interesting to compute a large q expansion of the
connectivity at which the hard fields appears in all the clusters (for all m such that Σ(m) > 0). Together with our
conjecture about rigidity being responsible for the computational hardness that might give a hint about the answer on
the long-standing question [22]: “Is there a polynomial algorithm and ǫ such that the algorithm would color random
graphs of average connectivity (1 + ǫ)q log q for all large q?”
2. The condensation transition
To compute the large-q asymptotic of the condensation transition, we first need to derive the large-q expansion of
the free entropy (30) in the connectivity regime c = 2q log q. Let us show self-consistently that the following scaling
is relevant for the condensation transition in the large q limit
cs = 2q log q − γ log q + α , (E9)
, (E10)
and compute the constants γ, α, B. Using the above scaling, the function w(η) (35) is dominated by the first two
terms in numerator and denominator, and reads in the first two leading orders 1 − qw(η) ≈ qe
so that
w(η) =
log q
(E11)
independently of γ, α, and B. To take into account the reweighting we expand eq. (38) in the two leading orders
log q
. (E12)
Note that almost all the incoming fields are hard, i.e. have one component of value 1. Since there are on average only
2B log q incoming soft fields, the leading order of the hard-field reweighting (the normalization in eq. (4)) is different
from 1 with a probability only O(log q/q). Similarly, almost all the soft fields have two nonzero and equal components.
The normalization in eq. (4) is thus almost surely 2, thus the average reweighting factor of the soft fields is
Zms = 2
m + O
log q
. (E13)
Finally,
log q
. (E14)
Therefore the constant B in (E10) is B = 2m/2, independently of γ and α.
The computation of the complexity requires the next order in the hard-field reweighting. Indeed the normalization
in (4) might not be 1 but 1/2; and this happens when there is a soft field arriving of the color corresponding to the
hard field in consideration. The probability of this event is
2c(1−qη)
= O( q
log q
). The hard-field reweighting is thus
1 − 2c (1 − qη)
2c (1 − qη)
log q
. (E15)
We now expand the replicated free energy (30) in the large q limit and regime (E9). Remind that from (6, 8)
Φs(m) = log (Z
m − c
. (E16)
The averages are over the population in the sense of (21).
The site free energy is the logarithm of the average of the total field normalization. This average can be split into
three parts when (i) the total field is a hard field, (ii) the total field is a soft field and (iii) the total field is contradictory
(and its normalization zero). The probability that the total field is not contradictory is the denominator in eq. (35)
g(η) =
(−1)l
l + 1
[1 − (l + 1)η]c, (E17)
where again only the first two terms are relevant in the expansion. The site free energy is then
log (Zi0)
m = log g(η) + log
q w(η)Zm
+ (1 − q w(η))Zms
≃ log
q(1 − η)c − q(q − 1)
(1 − 2η)c
+ log
1 − 1
− 2c (1 − qη)
1 − 1
.(E18)
where
q [1 − η]c + q (q − 1)
[1 − 2η]c
= log q + c log [1 − η] + log
1 − q − 1
1 − 2η
1 − η
(E19)
≈ log q + c log
1 − 1
+ log
1 − 1
+ o(1/q)
. (E20)
To compute the link contribution in (E16) we need to consider two fields ψi→js and ψ
s and to compute
0 = 1 −
ψi→js ψ
s . (E21)
There are three different cases:
1. Two hard fields are chosen, then Z
0 = 0 with probability qη
2 (this is of order 1/q) and Z
0 = 1 with probability
q(q − 1)η2 (this is of order 1).
2. Two soft fields are chosen then Z
0 = 1 with probability (1 − qη)2 (this is of order 1/q2), all other situations
being O(1/q3). Let us remind that the dominant soft fields are two-component of type (1/2, 1/2, 0, 0, . . . ).
3. One hard and one soft field is chosen, then Z
0 = 1 with probability 2η(1 − qη)(q − 2)/q, and Z
0 = 1/2 with
probability 4η(1 − qη) (this is of order 1/q2).
On average, one thus obtain for the link contribution
= log
1 − qη2 − 4η (1 − qη)
1m + 4η (1 − qη) 1
. (E22)
Putting together the two pieces (E18) and (E22), expanding η according to (E14) and considering only the highest
order in c, we can finally write the free energy as
Φs(m) = log q −
2m − 2
. (E23)
The internal entropy s(m) and the complexity Σ = Φs(m) −ms(m) are then
s(m) =
∂Φs(m)
2m log 2
, (E24)
Σ(m) = log q − c
2m − 2 −m 2m log 2
, (E25)
and the complexity is thus zero for cΣ=0 = 2q log q − log q − 2 + 2m [1 −m log 2] + o(1). In particular, one has for the
coloring and the condensation thresholds
cΣ=0(m = 0) = 2q log q − log q − 1 + o(1) , (E26)
cΣ=0(m = 1) = 2q log q − log q − 2 log 2 + o(1) . (E27)
For connectivity c = 2q log q − log q + α, one gets
2qs(m) ≃ 2m log 2 , (E28)
2qΣ(m) ≃ 2m − 2 −m2m log 2 − α . (E29)
[1] S. Cook, Proceedings of the 3rd Annual ACM Symposium. on Theory of Computing (ACM, New York, 1971), pp. 151-158.
M. R. Garey and D. S. Johnson Computers and Intractability: A Guide to the Theory of NP-Completeness, (Freeman,
New York, 1979). C.H Papadimitrious, Computational complexity (Addison-Wesley, 1994).
[2] K. Appel and W. Haken, Illinois J. Math. 21, 421 (1977). K. Appel and W. Haken, Illinois J. Math. 21, 491 (1977).
[3] R. E. Miller and J. W. Thatcher (eds.), Complexity of Computer Computation, Plenum Press. New York, 85 (1972).
[4] B. Bollobás, Random Graphs (CUP, Cambridge, 2001).
[5] S.Janson, T. Luczak and A.Rucinski, Random Graphs, John Wiley and Sons, 2000.
[6] P. Erdős and A. Rényi, Publ. Math. (Debrecen) 6, 290 (1959).
[7] D. Achlioptas and E. Friedgut, Random Structures and Algorithms, Volume 14, Issue 1, Pages 63 - 70 (1999).
[8] D. Achlioptas, A. Naor, Y. Peres, Nature, 435 (2005), 759-764.
[9] S. Dubois, O. Mandler, unpublished.
[10] D. Achlioptas, C. Moore, Journal of Comp. & Sys. Sci., 67 (2), (2003), p.441-471.
[11] D. Achlioptas, C. Moore, Proc. RANDOM 2004.
[12] L. Shi, N.C. Wormald, Research Report CORR 2004-24, Faculty of Mathematics, University of Waterloo (September 2004).
[13] G. Kemkes, N.C. Wormald, Research Report CORR 2005-22, Faculty of Mathematics, University of Waterloo (2005).
[14] J. Dı́az, G. Grammatikopoulos, A.C. Kaporis, L.M. Kirousis, X. Perez, and D.G. Sotiropoulos; Proceedings of the 13th
Annual European Symposium on Algorithms, 2005.
[15] T. Luczak, Combinatorica 11, 45 (1991).
[16] P. Cheeseman, B. Kanefsky, W.M. Taylor, Proceedings of the Twelfth International Joint Conference on Artificial Intelli-
gence, IJCAI-91, Sydney, Australia, 1991.
[17] B. Selman, S. Kirkpatrick, Artificial Intelligence 81 (1-2): 273-295 MAR 1996. B. Selman, D. G. Mitchell, H. J. Levesque,
Artificial Intelligence 81 (1-2): 17-29 MAR 1996.
[18] D.L. Mammen, T. Hogg, J. Artif. Intell. Res. (JAIR) 7: 47-66 (1997).
[19] C. P. Williams, T. Hogg, Artificial Intelligence 70 (1-2): 73-117 OCT 1994, T. Hogg, B. A. Huberman and C. P. Williams,
Artificial Intelligence 81 (1-2): 1-15 MAR 1996.
[20] S. Cocco and R. Monasson, Phys. Rev. Lett. 86 (8): 1654-1657 FEB 19 2001.
[21] D. Achlioptas and M. Molloy, in Proc. of FOCS-97 (IEEE, New York, 1997).
[22] A. Frieze and C.J.H. McDiarmid, Rand. Struc. Alg. 10, 5 (1996).
[23] R. Monasson, R. Zecchina, S. Kirkpatrick, B. Selman, and L. Troyansky, Nature 400, 133 (1999).
[24] M. Mézard, Parisi G and Zecchina R, Science 297, 812 (2002).
[25] F.Y. Wu, Rev. Mod. Phys. 54, 235 (1982).
[26] Mézard M and Parisi G, 2001 Eur. Phys. J. B 20 217.
[27] Mézard M and Parisi G, 2003 J. Stat. Phys 111 no. 1-2, 1-34.
[28] M. Mézard and R. Zecchina, Phys. Rev. E, 66, 056126, (2002).
[29] Mézard M, Parisi G and Virasoro M A, Spin-glass Theory and Beyond, vol. 9 of Lecture notes in Physics. World Scientific,
Singapore, 1987
[30] R. Mulet, A. Pagnani, M. Weigt and R. Zecchina, Phys. Rev. Lett. 89, 268701 (2002). A. Braunstein, R. Mulet, A. Pagnani,
M. Weigt, and R. Zecchina, Phys. Rev. E 68, 036702 (2003).
[31] F. Krzakala, A. Pagnani, and M. Weigt, Phys. Rev. E 70, 046705 (2004).
[32] J. van Mourik and D. Saad, Phys. Rev. E 66, 056120 (2002).
[33] A. Braunstein, M. Mézard, and R. Zecchina, Random Struct. Algorithms 27, 201 (2005).
[34] Montanari A, Parisi G and Ricci-Tersenghi F, J. Phys. A: Math. Gen. 37 2073-2091.
[35] G. Biroli, R. Monasson, M. Weigt, Eur. Phys. J. B 14, 551 (2000).
[36] M. Mézard, T. Mora and R. Zecchina, Phys. Rev. Lett. 94, 197205 (2005).
[37] D. Achlioptas, F. Ricci-Tersenghi, Proc. 38th annual ACM symposium on Theory of computing (2006), p. 130.
[38] M. Mézard, M. Palassini and O. Rivoire, Phys. Rev. Lett. 95, 200202 (2005).
[39] O. Rivoire, These de doctorat, CNRS http://tel.ccsd.cnrs.fr/tel-00009956/en/ (2005).
[40] Pearl J, Probabilistic Reasoning in Intelligent Systems: Networks of Plausible Inference (2nd edition). Morgan Kaufmann
Publishers, San Francisco, 1988.
[41] Yedidia, J.S., Freeman, W.T. and Weiss, Y., Exploring Artificial Intelligence in the New Millennium, ISBN 1558608117,
Chap. 8, pp. 239-236, January 2003 (Science and Technology Books).
[42] F. Krzakala, A. Montanari, F. Ricci-Tersenghi, G. Semerjian, L. Zdeborová, Proc. Natl. Acad. Sci. 104, 10318 (2007).
[43] A. Montanari, F. Ricci-Tersenghi, G. Semerjian, in preparation (2007).
[44] D. J. Gross, I. Kanter and H. Sompolinsky Phys. Rev. Lett. 55, 304-307 (1985).
[45] H. O. Georgii, Gibbs Measures and Phase Transitions, Berlin (1988).
[46] S. Franz, M. Mézard, F. Ricci-Tersenghi, M. Weigt and R. Zecchina, Europhys. Lett. 55 (2001) 465.
[47] H. A. Bethe, Proc. R. Soc. Lond. 150, 552 (1935). R. Peierls, Proc. Roy. Soc. 154, 207 (1936). L Onsager - Journal of the
American Chemical Society, (1936). P. R. Weiss, Phys. Rev. 74, 1493 - 1504 (1948).
[48] K.Y. Wong and D. Sherrington, J. Phys. A 21 (1988) L459. Y.Y. Goldschmidt and P.-Y.Lai, J. Phys.A 23 (1990) L775.
[49] D. Aldous, The zeta(2) Limit in the Random Assignment Problem. Random Structures and Algorithms 18 (2001) 381-418.
[50] Bandyopadhyay A and Gamarnik D, 2006 ACM-SIAM Symp. on Discrete Algorithms, 890.
[51] J. Jonasson, Statistics and Probability Letters 57 (2002), 243-248.
[52] A. Montanari, Devavrat Shah, Proceedings of Symposium of Discrete Algorithms (SODA), New Orleans, September 2007.
[53] M. Mézard, A. Montanari, J. Stat. Phys. 124, 1317-1350 (2007).
[54] J.P. Bouchaud, G. Biroli, J. Chem. Phys. 121, 7347 (2004).
[55] A. Montanari, G. Semerjian, J. Stat. Phys. 124, 103-189 (2006); Phys. Rev. Lett. 94, 247201 (2005)
[56] O. Rivoire, G. Biroli, O. C. Martin, M. Mézard, Eur. Phys. J. B 37, 55-78 (2004).
[57] H. Kesten, B. P. Stigum, Ann Math Statist 37, 1463 (1966). H. Kesten, B. P. Stigum, J Math Anal Appl 17, 309 (1966)
[58] Parisi G, J. Phys. A13 L115-L121, (1980).
[59] R. Monasson. Phys. Rev. Lett. 75, 2847 (1995).
[60] B. Derrida, Phys. Rev. Lett. 45, 79-82 (1980). B. Derrida, Phys. Rev. B 24, 2613-2626 (1981).
[61] T. Kirkpatrick and D. Thirumalai, Phys. Rev. Lett. 58, 2091 (1987). T. Kirkpatrick, D. Thirumalai, and P. Wolynes, Phys.
Rev. A 40, 1045 (1989). M. Mézard and G. Parisi, Phys. Rev. Lett. 82 (1999) 747. M. Mézard, in ”More is different”, Ong
and Bhatt editors, Princeton (2002).
[62] A. Montanari, G. Semerjian, J. Stat. Phys. 125, 23 (2006).
[63] W Kauzmann, Chemical Reviews 43 (2), 1948.
[64] G. Semerjian, arXiv:0705.2147v1 (2007).
[65] G. Parisi, cs.CC/0212047.
[66] E. N. Maneva, E. Mossel and M. J. Wainwright, in Proceedings of SODA, 2005.
[67] A. Braunstein, R. Zecchina, J. Stat. Mech. P06007 (2004).
[68] B. Pittel, J. Spencer, N. Wormald, J. of Combinatorial Theory, Series B, 67, 1996, pp. 111–151.
[69] E. Gardner, Nucl. Phys. B297, [FS14] 74 (1985).
[70] A. Montanari and F. Ricci-Tersenghi, Eur. Phys. J. B 33, 339 (2003)
[71] G. Parisi and T. Rizzo, Phys. Rev. B 72, 184431 (2005). M. Müller, L. Leuzzi, and A. Crisanti, Phys. Rev. B 74, 134431
(2006).
[72] C.S. O’Hern, L. E. Silbert, A. J. Liu and S. R. Nagel, Phys. Rev. E 68, 011306 (2003).
[73] F. Krzakala and J. Kurchan, cond-mat/0702546.
[74] J. Chavas, C. Furtlehner, M. Mézard and R. Zecchina, J. Stat. Mech. (2005) P11016
[75] P. Svenson, G. Nordahl, Phys. Rev. E, 59, n. 4, 3983 (1999).
[76] B. Selman, H. A. Kautz and B. Cohen, proc. of AAAI-94, Seattle (1994);
[77] J. Ardelius, E. Aurell, Phys. Rev. E 74, 037702 (2006).
[78] S. Seitz, M. Alava, P. Orponen, J Stat. Mech. (2005) P06006.
[79] S. Boettcher and A. G. Percus, Phys. Rev. E 64,026114(2001).
http://tel.ccsd.cnrs.fr/tel-00009956/en/
http://arXiv.org/abs/0705.2147
http://arXiv.org/abs/cs/0212047
http://arXiv.org/abs/cond-mat/0702546
[80] Mézard M., Parisi G., Virasoro M.A., J. Physique Lett. 46 (1985) L-217.
[81] M. Talagrand , Spin Glasses: A challenge for Mathematicians, Springer, 2000.
[82] Pitman J., Yor M., Ann. Probab. 25, 855-900 (1997).
[83] L. Devroye, Non-Uniform Random Variate Generation, Springer-Verlag, New York, 1986
[84] T. Castellani, F. Krzakala and F. Ricci-Tersenghi; Eur. Phys. J. B 47, 99 (2005)
[85] More precisely [7] proves the existence of a sharp threshold for a possibly size dependent sequence of thresholds, whose
convergence is not proven.
[86] In this paper we use the words ”cluster” and ”state” as synonyms.
[87] A slightly different adaptation, closer to [76], was performed few years ago by Andrea Pagnani and Martin Weigt (private
communication).
	Introduction
	The Model
	Definition of the model
	Ensembles of Random Graphs
	The cavity formalism at the replica symmetric level
	The replica symmetric cavity equations
	Average over the ensemble of graphs and the RS solution
	Validity conditions of the replica symmetric solution
	The Gibbs measure uniqueness condition
	The Gibbs measure extremality condition
	The local stability: a simple self-consistency check
	One-step replica symmetry breaking framework
	Analyzing the 1RSB equations
	Zero temperature limit
	The role of frozen variables
	Hard fields in the simplest case, m=0
	Generalized survey propagation recursion
	The presence of frozen variables
	Validity conditions of the 1RSB solution
	The coloring of random graphs: cavity results
	Regular random graphs
	Results for the bi-regular ensemble
	Results for Erdos-Rényi random graphs
	The general case for q>3, discontinuous clustering transition
	The special case of 3-coloring, continuous clustering transition
	The overlap structure
	Large q Asymptotics
	Algorithmic consequences
	The whitening procedure
	A Walk-COL algorithm to color random graph
	A belief propagation algorithm to color random graph
	Conclusions
	Acknowledgments
	Stability of the paramagnetic solution
	The relative sizes of clusters in the condensed phase 
	The 1RSB formalism at m=1 and the reconstruction equations
	Numerical Methods
	High-q asymptotics
	The appearance of hard fields at m=1
	The condensation transition
	References

  We consider the problem of coloring the vertices of a large sparse random
graph with a given number of colors so that no adjacent vertices have the same
color. Using the cavity method, we present a detailed and systematic analytical
study of the space of proper colorings (solutions).
  We show that for a fixed number of colors and as the average vertex degree
(number of constraints) increases, the set of solutions undergoes several phase
transitions similar to those observed in the mean field theory of glasses.
First, at the clustering transition, the entropically dominant part of the
phase space decomposes into an exponential number of pure states so that beyond
this transition a uniform sampling of solutions becomes hard. Afterward, the
space of solutions condenses over a finite number of the largest states and
consequently the total entropy of solutions becomes smaller than the annealed
one. Another transition takes place when in all the entropically dominant
states a finite fraction of nodes freezes so that each of these nodes is
allowed a single color in all the solutions inside the state. Eventually, above
the coloring threshold, no more solutions are available. We compute all the
critical connectivities for Erdos-Renyi and regular random graphs and determine
their asymptotic values for large number of colors.
  Finally, we discuss the algorithmic consequences of our findings. We argue
that the onset of computational hardness is not associated with the clustering
transition and we suggest instead that the freezing transition might be the
relevant phenomenon. We also discuss the performance of a simple local Walk-COL
algorithm and of the belief propagation algorithm in the light of our results.

Introduction
	The Model
	Definition of the model
	Ensembles of Random Graphs
	The cavity formalism at the replica symmetric level
	The replica symmetric cavity equations
	Average over the ensemble of graphs and the RS solution
	Validity conditions of the replica symmetric solution
	The Gibbs measure uniqueness condition
	The Gibbs measure extremality condition
	The local stability: a simple self-consistency check
	One-step replica symmetry breaking framework
	Analyzing the 1RSB equations
	Zero temperature limit
	The role of frozen variables
	Hard fields in the simplest case, m=0
	Generalized survey propagation recursion
	The presence of frozen variables
	Validity conditions of the 1RSB solution
	The coloring of random graphs: cavity results
	Regular random graphs
	Results for the bi-regular ensemble
	Results for Erdos-Rényi random graphs
	The general case for q>3, discontinuous clustering transition
	The special case of 3-coloring, continuous clustering transition
	The overlap structure
	Large q Asymptotics
	Algorithmic consequences
	The whitening procedure
	A Walk-COL algorithm to color random graph
	A belief propagation algorithm to color random graph
	Conclusions
	Acknowledgments
	Stability of the paramagnetic solution
	The relative sizes of clusters in the condensed phase 
	The 1RSB formalism at m=1 and the reconstruction equations
	Numerical Methods
	High-q asymptotics
	The appearance of hard fields at m=1
	The condensation transition
	References
