Statistical Physics of Hard Optimization Problems

Acknowledgment

First of all I would like to express my thanks to my advisor Marc Mézard from whom I learned a lot. He showed me how to combine enthusiasm, patience, computations, and intuition in the correct proportion to enjoy the delicious taste of the process of discovery. I thank as well to my advisor Václav Janiš who guided my scientific steps mainly in the earlier stages of my work.

This work would never be possible without the contact and discussions with my colleagues and collaborators all over the world. Without them I would feel lost in the vast world of unknown. I am grateful to all the organizers of workshops, conferences and summer schools where I had participated. I also thank for invitations on visits and seminars which were always very inspiring.

I am very thankful to the whole Laboratoire de Physique Théorique et Modèles Statistiques for a very warm reception and to all its members for helping me whenever I needed. I also owe a lot to my professors from the Charles University in Prague, and to my colleagues from the Institute of Physics of the Academy of Sciences in Prague. I also thank the referees of my thesis and other members of the thesis committee for their interest in my work and for accepting this task.

I value profoundly the scholarship granted by the French Government which covered the largest part of my stay in France. Further, I appreciated greatly the subvention from the French ministry for higher education and research ”cotutelles internationales de thése”. I also acknowledge gratefully the support from the FP6 European network EVERGROW.

My deepest thanks go to my parents for their constant support, encouragement, and love. Finally, thank you Flo for all the items above and many more. You taught me well, and yes I am the Jedi now.

Title: Statistical Physics of Hard Optimization Problems

Abstract: Optimization is fundamental in many areas of science, from computer science and information theory to engineering and statistical physics, as well as to biology or social sciences. It typically involves a large number of variables and a cost function depending on these variables. Optimization problems in the NP-complete class are particularly difficult, it is believed that the number of operations required to minimize the cost function is in the most difficult cases exponential in the system size. However, even in an NP-complete problem the practically arising instances might, in fact, be easy to solve. The principal question we address in this thesis is: How to recognize if an NP-complete constraint satisfaction problem is typically hard and what are the main reasons for this? We adopt approaches from the statistical physics of disordered systems, in particular the cavity method developed originally to describe glassy systems. We describe new properties of the space of solutions in two of the most studied constraint satisfaction problems - random satisfiability and random graph coloring. We suggest a relation between the existence of the so-called frozen variables and the algorithmic hardness of a problem. Based on these insights, we introduce a new class of problems which we named ”locked” constraint satisfaction, where the statistical description is easily solvable, but from the algorithmic point of view they are even more challenging than the canonical satisfiability.

Keywords: Constraint satisfaction problems, combinatorial optimization, random coloring problem, average computational complexity, cavity method, spin glasses, replica symmetry breaking, Bethe approximation, clustering of solutions, phase transitions, message passing, belief propagation, satisfiability threshold, reconstruction on trees.

Titre: Physique statistique des problèmes d’optimisation

Résumé: L’optimisation est un concept fondamental dans beaucoup de domaines scientifiques comme l’informatique, la théorie de l’information, les sciences de l’ingénieur et la physique statistique, ainsi que pour la biologie et les sciences sociales. Un problème d’optimisation met typiquement en jeu un nombre important de variables et une fonction de coût qui dépend de ces variables. La classe des problèmes NP-complets est particulièrement difficile, et il est communément admis que, dans le pire des cas, un nombre d’opérations exponentiel dans la taille du problème est nécessaire pour minimiser la fonction de coût. Cependant, même ces problèmes peuveut être faciles à résoudre en pratique. La principale question considérée dans cette thèse est comment reconnaître si un problème de satisfaction de contraintes NP-complet est ”typiquement” difficile et quelles sont les raisons pour cela ? Nous suivons une approche inspirée par la physique statistique des systèmes desordonnés, en particulier la méthode de la cavité développée originalement pour les systèmes vitreux. Nous décrivons les propriétés de l’espace des solutions dans deux des problèmes de satisfaction les plus étudiés : la satisfiabilité et le coloriage aléatoire. Nous suggérons une relation entre l’existence de variables dites ”gelées” et la difficulté algorithmique d’un problème donné. Nous introduisons aussi une nouvelle classe de problèmes, que nous appelons ”problèmes verrouillés”, qui présentent l’avantage d’être à la fois facilement résoluble analytiquement, du point de vue du comportement moyen, mais également extrêmement difficiles du point de vue de la recherche de solutions dans un cas donné.

Les mots clefs: Problèmes d’optimisation de contraintes, optimisation combinatoire, problèmes de coloriage, complexité de calcul moyenne, méthode de la cavité, verres de spins, brisure de la symétrie des répliques, approximation de Bethe, regroupement des solutions en amas, transitions de phases, passage de messages, propagation des convictions, seuil de satisfiabilité, reconstruction sur des arbres.

Název: Statistická fyzika těžkých optimalizačních úloh

Abstrakt: Optimalizace je fundamentální koncept v mnoha vědních oborech, počínaje počítačovou vědou a teorií informace, přes inženýrství a statistickou fyziku, až po biologii či ekonomii. Optimalizační úloha se typicky skládá z minimalizace funkce závisející na velkém množství proměnných. Problémy z takzvané NP-úplné třídy jsou obzvláště složité, věří se, že počet operací potřebný k nalezení řešení v tom nejtěžším případě roste exponenciálně s počtem proměných. Nicméně i pro NP-úplné úlohy platí, že praktické případy mohou být jednoduché. Hlavní otázka, kterou se zabývá tato práce, je: Jak rozpoznat, zda je NP-úplný problém splnitelnosti podmínek v typickém případě těžký a čím je tato složitost způsobena? K této otázce přistupujeme s využitím znalostí ze statistické fyziky neuspořádaných, a zejména skelných, systémů. Popíšeme nové vlastnosti prostoru řešení ve dvou z nejvíce studovaných optimalizačních problémů – splnitelnosti náhodných Booleovských formulí a barvení náhodných grafů. Navrhneme existenci vztahu mezi typickou algoritmickou složitostí a existencí takzvaně zamrzlých proměnných. Na základě těchto poznatků zkonstruujeme novou třídu problémů, které jsme nazvali ”uzamknuté”, zde je statistický popis množiny všech řešení poměrně jednoduchý, ale z algoritmického pohledu jsou tyto typické případy těchto problémů ještě težši než v kanonickém problému splitelnosti Booleovských formulí.

Klíčová slova: Problémy splnitelnosti podmínek, kombinatorická optimalizace, barvení náhodných grafů, průměrná algoritmická složitost, metoda kavity, spinová skla, narušení symetrie replik, Betheho aproximace, shlukování řešení, fázové přechody, posílání zpráv, propagace domněnek, práh splnitelnosti, rekonstrukce na stromech.

Foreword

P.-G. de Gennes in his foreword to the book ”Stealing the gold – A celebration of the pioneering physics of Sam Edwards” wrote:

But he {meaning S. Edwards} also has another passion, which I {meaning P.-G. de Gennes} call ”The search for unicorns.” To chase unicorns is a delicate enterprise. Medieval Britons practised it with great enthusiasm (and this still holds up to now: read Harry Potter). Sir Samuel Edwards is not far from the gallant knights of the twelfth century. Discovering a strange animal, approaching it without fear, then not necessarily harnessing the creature, but rapidly drawing a plausible sketch of its main features.

One beautiful unicorn prancing in the magic garden of Physics has been names ”Spin glass.” It is rare: not many pure breeds of Spin glasses have been found in Nature. But we have all watched the unpredictable jumps of this beast. And we have loved its story – initiated by Edwards and Anderson.

Unicorn is a mythical animal, described in the book of Job, together with another strange and fascinating creature which is less peaceful: the leviathan. Leviathans are described as immense terrible monsters, invincible beasts. Most people prefer not to even think about them.

This thesis tells a story about what happens when the fierce and mysterious beauty of a unicorn meets with the invincibility of a leviathan.

Chapter 1 Hard optimization problems

In this opening chapter we introduce the constraint satisfaction problems and discuss briefly the computer science approach to the computational complexity. We review the studies of the random satisfiability problem in the context of average computational complexity investigations. We describe the connection between spin glasses and random CSPs and highlight the most interesting results coming out from this analogy. We explain the replica symmetric approach to these problems and show its usefulness on the example of counting of matchings [ZDEB-1]. Then we review the survey propagation approach to constraint satisfaction on an example of 1-in-$K$ satisfiability [ZDEB-3]. Finally we summarize the main contributions of the author to the advances in the statistical physics of hard optimization problems, that are elaborated in the rest of the thesis.

Optimization is a common concept in many areas of human activities. It typically involves a large number of variables, e.g. particles, agents, cells or nodes, and a cost function depending on these variables, such as energy, measure of risk or expenses. The problem consists in finding a state of variables which minimizes the value of the cost function.

In this thesis we will concentrate on a subset of optimization problems the so-called constraint satisfaction problems (CSPs). Constraint satisfaction problems are one of the main building blocks of complex systems studied in computer science, information theory and statistical physics. Their wide range of applicability arises from their very general nature: given a set of $N$ discrete variables subject to $M$ constraints, the CSP consists in deciding whether there exists an assignment of variables which satisfies simultaneously all the constraints. And if such an assignment exists then we aim at finding it.

In computer science, CSPs are at the core of computational complexity studies: the satisfiability of boolean formulas is the canonical example of an intrinsically hard, NP-complete, problem. In information theory, error correcting codes also rely on CSPs. The transmitted information is encoded into a codeword satisfying a set of constraints, so that the information may be retrieved after transmission through a noisy channel, using the knowledge of the constraints satisfied by the codeword. Many other practical problems in scheduling a collection of tasks, in electronic design engineering or artificial intelligence are viewed as CSPs. In statistical physics the interest in CSPs stems from their close relation with the theory of spin glasses. Answering if frustration is avoidable in a system is the first, and sometimes highly nontrivial, step in understanding the low temperature behaviour.

A key point is to understand how difficult it is to solve practical instances of a constraint satisfaction problem. Everyday experience confirms that sometimes it is very hard to find a solution. Many CSPs require a combination of heuristics and combinatorial search methods to be solved in a reasonable time. A key question we address in this thesis is thus why and when are some instances of these problems intrinsically hard. Answering this question has, next to its theoretical interest, several practical motivations

Understanding where the hardness comes from helps to push the performance of CSPs solvers to its limit.

Understanding which instances are hard helps to avoid them if the nature of the given practical problem permits.

Finding the very hard problem might be interesting for cryptographic application.

A pivotal step in this direction is the understanding of the onset of hardness in random constraint satisfaction problems. In practice random constraint satisfaction problems are either regarded as extremely hard as there is no obvious structure to be explored or as extremely simple as they permit probabilistic description. Furthermore, random constraint satisfaction models are spin glasses and we shall thus borrow methods from the statistical physics of disordered systems.

2 Constraint Satisfaction Problems: Setting

Constraint Satisfaction Problem (CSP): Consider $N$ variables $s_{1}\dots,s_{N}$ taking values from the domain $\{0,\dots,q-1\}$, and a set of $M$ constraints. A constraint $a$ concerns a set of $k_{a}$ different variables which we call $\partial a$. Constraint $a$ is a function from all possible assignments of the variables $\partial a$ to $\{0,1\}$. If the constraint evaluates to $1$ we say it is satisfied, and if it evaluates to we say it is violated. The constraint satisfaction problem consists in deciding whether there exists an assignment of variables which satisfies simultaneously all the constraints. We call such an assignment a solution of the CSP.

In physics, the variables represent $q$-state Potts spins (or Ising spins if $q=2$). The constraints represent very general (non-symmetric) interactions between $k_{a}$-tuples of spins. In Boolean constraint satisfaction problems ($q=2$) a literal is a variable or its negation. A clause is then a disjunction (logical OR) of literals.

A handy representation for a CSP is the so-called factor graph, see [KFL01] for a review. Factor graph is a bipartite graph $G(V,F,E)$ where $V$ is the set of variables (variables nodes, represented by circles) and $F$ is the set of constraints (function nodes, represented by squares). An edge $(ia)\in E$ is present if the constraint $a\in F$ involves the variable $i\in V$. A constraint $a$ is connected to $k_{a}$ variables, their set is denoted $\partial a$. A variable $i$ is connected to $l_{i}$ constraints, their set is denoted $\partial i$. For clarity we specify the factor graph representation for the graph coloring and exact cover problem in fig. 1.1, both defined in the following section 1.2.2.

2.2 List of CSPs discussed in this thesis

Here we define constraint satisfaction problems which will be discussed in the following. Most of them are discussed in the classical reference book [GJ79]. The most studied constraint satisfaction problems are defined over Boolean variables, $q=2$, $s_{i}\in\{0,1\}$. Sometimes we use equivalently the notation with Ising spins $s_{i}\in\{-1,+1\}$. CSPs with Boolean variables that we shall discuss in this thesis are:

Satisfiability (SAT) problem: Constraints are clauses, that is logical disjunctions of literals (i.e., variables or their negations). Example of a satisfiable formula with 3 variables and 4 clauses (constraints) and 10 literals: $(x_{1}\vee x_{2}\vee\neg x_{3})\wedge(x_{2}\vee x_{3})\wedge(\neg x_{1}\vee\neg x_{3})\wedge(x_{1}\vee\neg x_{2}\vee x_{3})$.

$K$-SAT: Satisfiability problem where every clause involves $K$ literals, $k_{a}=K$ for all $a=1,\dots,M$.

Not-All-Equal SAT: Constraints are satisfied everytime except when all the literals they involve are TRUE or all of them are FALSE.

Bicoloring: Constraints are satisfied except when all variables they involve are equal. Bicoloring is Not-All-Equal SAT without negations.

XOR-SAT: Constraints are logical XORs of literals.

Odd (resp. Even) Parity Checks: A constraint is satisfied if the sum of variables it involves is odd (resp. even). Odd parity checks are XORs without negations.

1-in-$K$ SAT: Constraints are satisfied if exactly one of the $K$ literals they involve is TRUE.

Exact Cover, or positive 1-in-$K$ SAT: Constraints are satisfied if exactly one of the $K$ variables they involve is $1$ (occupied). Exact cover, or positive 1-in-$K$ SAT, is 1-in-$K$ SAT without negations.

Perfect matching: Nodes of the original graph become constraints, variables are on edges and determine if the edge is or is not in the matching, see fig. 1.5. Constraints are satisfied if exactly one of the $K$ variables they involve is $1$ (belongs to the matching). Note that perfect matching is just a variant of the Exact Cover

Occupation problems are defined by a binary $(K+1)$ component vector $A$. All constraints involve $K$ variables, and are satisfied if the sum of variables they involve $r=\sum_{\partial a}s_{i}$ is such that $A_{r}=1$.

Locked Occupation Problems (LOPs): If the vector $A$ is such that $A_{i}A_{i+1}=0$ for all $i=0,\dots,K-1$, and all the variables are present in at least two constraints.

We will also consider in a great detail one CSP with $q$-ary variables: The graph coloring with $q$ colors: Every constraint involves two variables and is satisfied if the two variables are not assigned the same value (color). In physics the $q$-ary variables are called Potts spins.

2.3 Random factor graphs: definition and properties

Given a constraint satisfaction problem with $N$ variables and $M$ constraints, the constraint density is defined as $\alpha=M/N$. Denote by ${\cal R}(k)$ the probability distribution of the degree of constraints (number of neighbours in the factor graph), and by ${\cal Q}(l)$ the probability distribution of the degree of variables. The average connectivity (degree) of constraints is

The constraint density is then asymptotically

A random factor graph with a given $N$ and $M$ is then created as follows: Draw a sequence $\{l_{1},\dots,l_{N}\}$ of $N$ numbers from the distribution ${\cal Q}(l)$. Subsequently, draw a sequence $\{k_{1},\dots,k_{M}\}$ of $M$ numbers from the distribution ${\cal R}(k)$, such that $\sum_{a=1}^{M}k_{i}=\sum_{i=1}^{N}l_{i}$. The random factor graph is drawn uniformly at random from all the factor graphs with $N$ variables, $M$ constraints and degree sequences $\{l_{1},\dots,l_{N}\}$ and $\{k_{1},\dots,k_{M}\}$.

Another definition leading to a Poissonian degree distribution is used often if the degree of constraints is fixed to $K$ and the number of variables is fixed to $N$. There are ${N\choose K}$ possible positions for a constraint. Each of these positions is taken with probability

The number of constraints is then a Poissonian random variable with average $M=cN/K$. The degree of variables is distributed according to a Poissonian law with average $c$

If $K=2$ these are the random Erdős-Rényi graphs [ER59]. This definition works also if constraints are changed for variables, that is if the degree of variables and the number of constraints are fixed, as in e.g. the matching problem.

The random factor graphs are called regular if both the degrees of constraints and variables are fixed, ${\cal R}(k)=\delta(k-K)$ and ${\cal Q}(l)=\delta(l-L)$. In section 4.3 we will also use the truncated Poissonian degree distribution

The average connectivity for the truncated Poissonian distribution is then

In the cavity approach, the so-called excess degree distribution is a crucial quantity. It is defined as follows: Choose an edge $(ij)$ at random and consider the probability distribution of the number of neighbours of $i$ except $j$. The variables (analogously for constraints) excess degree distribution thus reads

We will always deal with factor graphs where $K$ and $c$ are of order one, and $N\to\infty,M\to\infty$. These are called sparse random factor graphs. Concerning the physical properties of sparse random factor graphs the two definitions of a random graph with Poissonian degree distribution are equivalent. Some properties (e.g. the annealed averages) can however depend on the details of the definition.

Consider a random variable $i$ in the factor graph. We want to estimate the average length of the shortest cycle going through variable $i$. Consider a diffusion algorithm spreading into all direction but the one it came from. The probability that this diffusion will arrive back to $i$ in $d$ steps reads

where $\gamma_{l}=\overline{l^{2}}/\overline{l}-1$ and $\gamma_{k}=\overline{k^{2}}/\overline{k}-1$ are the mean values of the excess degree distribution (1.8). The probability (1.9) is almost surely zero if

An important property follows: As long as the degree distributions ${\cal R}(k)$ and ${\cal Q}(l)$ have a finite variance the sparse random factor graphs are locally trees up to a distance scaling as $\log{N}$ (1.10). We define this as the tree-like property.

In this thesis we consider only degree distributions with a finite variance. A generalization to other cases (e.g. the scale-free networks with long-tail degree distributions) is not straightforward and many of the results which are asymptotically exact on the tree-like structures would be in general only approximative. We observed, see e.g. fig. 2.2, that many of the nontrivial properties predicted asymptotically on the tree-like graphs seems to be reasonably precise even on graphs with about $N=10^{2}-10^{4}$ variables. It means that the asymptotic behaviour sets in rather early and does not, in fact, require $\log N\gg 1$.

3 Computational complexity

Theoretical computer scientists developed the computational complexity theory in order to quantify how hard problems can be in the worst possible case. The most important and discussed complexity classes are the P, NP and NP-complete.

A problem is in the P (polynomial) class if there is an algorithm which is able to solve the problem for any input instance of length $N$ in at most $cN^{k}$ steps, where $k$ and $c$ are constants independent of the input instance. The formal definitions of what is a ”problem”, its ”input instance” and an ”algorithm” was formalized in the theory of Turing machines [Pap94], where the definition would be: The complexity class P is the set of decision problems that can be solved by a deterministic Turing machine in polynomial time. A simple example of polynomial problem is sorting a list of $N$ real numbers.

A problem is in the NP class if its instance can be stored in memory of polynomial size and if the correctness of a proposed result can be checked in polynomial time. Formally, the complexity class NP is the set of decision problems that can be solved by a non-deterministic Turing machine in polynomial time [Pap94], NP stands for non-deterministic polynomial. Whereas the deterministic Turing machine is basically any of our today computers, the non-deterministic Turing machine can perform unlimited number of parallel computations. Thus, if for finite $N$ there is a finite number of possible solutions all of them can be checked simultaneously. This class contains many problems that we would like to be able to solve efficiently, including the Boolean satisfiability problem, the traveling salesman problem or the graph coloring. Problems which do not belong to the NP class are for example counting the number of solutions in Boolean satisfiability, or the random energy model [Der80, Der81].

All the polynomial problems are in the NP class. It is not known if all the NP problems are polynomial, and it is considered by many to be the most challenging problem in theoretical computer science. It is also one of the seven, and one of the six still open, Millennium Prize Problems that were stated by the Clay Mathematics Institute in 2000 (a correct solution to each of these problems results in a $1,000,000 prize for the author). A majority of computer scientists, however, believes that the negative answer is the correct one [Gas02].

The concept of NP-complete problems was introduced by Cook in 1971 [Coo71]. All the NP problems can be polynomially reduced to any NP-complete problem, thus if any NP-complete problem would be polynomial then P$=$NP. Cook proved [Coo71] that the Boolean satisfiability problem is NP-complete. Karp soon after added 21 new NP-complete problems to the list [Kar72]. Since then thousands of other problems have been shown to be NP-complete by reductions from other problems previously shown to be NP-complete; many of these are collected in the Garey and Johnson’s ”Guide to NP-Completeness” [GJ79].

Schaefer in 1978 proved a dichotomy theorem for Boolean ($q=2$) constraint satisfaction problems. He showed that if the constraint satisfaction problem has one of the following four properties then it is polynomial, otherwise it is NP-complete. (1) All constraints are such that $s_{i}=1$ for all $i$ is a solution or $s_{i}=0$ for all $i$ is a solution. (2) All constraints concern at most two variables (e.g. in 2-SAT). (3) All constraints are linear equations modulo two (e.g. in XOR-SAT). (4) All constraints are the so-called Horn clauses or all of them are the so-called dual Horn clauses. A Horn clause is a disjunction of variables such that at most one variable is not negated. A dual Horn clause is when at most one variable is negated. A similar dichotomy theorem exists for 3-state variables, $q=3$, [Bul02]. Generalization for $q>3$ is not known.

3.2 The average case hardness

Given the present knowledge, it is often said that all the polynomial problems are easy and all the NP-complete problems are very hard. But, independently if P$=$NP or not, even polynomial problems might be practically very difficult, and some (or even most) instances of the NP-complete problems might be practically very easy.

An example of a still difficult polynomial problem is the primality testing, a first polynomial algorithm was discovered by [AKS04]. But a ”proof” of remaining difficulty is the EFF prize [EFF] of $100,000 to the first individual or group who discovers the first prime number with at least 10,000,000 decimal digits.

And how hard are the NP-complete problems? One way to answer is that under restrictions on the structure an NP-complete problem might become polynomial. Maybe the most famous example is 4-coloring of maps (planar factor graphs) which is polynomial. Moreover, it was a long standing conjecture that every map is colorable with 4 colors, proven by Appel and Haken [AH77b, AH77a]. Interestingly enough 3-coloring of maps is NP-complete [GJ79].

But there are also settings under which the problem stays NP-complete and yet almost every instance can be solved in polynomial time. A historically important example is the Boolean satisfiability where each clause is generated by selecting literals with some fixed probability. Goldberg introduced this random ensemble and showed that the average running time of the Davis-Putnam algorithm [DP60, DLL62] is polynomial for almost all choices of parameter settings [Gol79, GPB82]. Thus in the eighties some computer scientist tended to think that all the NP-complete problems are in fact on average easy and it is hard to find the evil instances which makes them NP-complete.

The breakthrough came at the beginning of the nineties when Cheeseman, Kanefsky and Taylor asked ”Where the really hard problems are?” in their paper of the same name [CKT91]. Shortly after Mitchell, Selman and Levesque came up with a similar work [MSL92]. Both groups simply took a different random ensemble of the satisfiability (in the second case) and coloring (in the first case) instances: the length of clauses is fixed to be $K$ and they are drawn randomly as described in sec. 1.2.3. They observed that when the density of clauses $\alpha=M/N$ is small the existence of a solution is very likely and if $\alpha$ is large the existence of a solution is very unlikely. And the really hard instances were located nearby the critical value originally estimated to be $\alpha_{s}\approx 4.25$ in the 3-SAT [MSL92]. The hardness was judged from the median running time of the Davis-Putnam-Logemann-Loveland (DPLL) backtracking-based algorithm [DP60, DLL62], see fig. 1.2. This whipped away the thoughts that NP-complete problems might in fact be easy on average. Many other studies and observations followed.

The hard instances of random $K$-satisfiability became very fast important benchmarks for the best algorithms. Moreover, there are some indications that critically constrained instances might appear in real-world applications. One may imagine that in a real world situation the amount of constraints is given by the nature of the problem, and variables usually correspond to something costly, thus the competitive designs contain the smallest possible number of variables.

Given a random $K$-SAT formula of $N$ variables the probability that it is satisfiable, plotted in fig. 1.3 for 3-SAT, becomes more and more like a step-function as the size $N$ grows. An analogy with phase transitions in physics cannot be overlooked. The existence and sharpness of the threshold were partially proved [Fri99]. The best known probabilistic bounds of the threshold value in 3-SAT are $3.520$ for the lower bound [KKL03, HS03] and $4.506$ for the upper bound [DBM00]. Numerical estimates of the asymptotic value of the threshold are $\alpha_{s}\approx 4.17$ [KS94], $\alpha_{s}\approx 4.258$ [CA96], $\alpha_{s}\approx 4.27$ [MZK+99b, MZK+99a]. The finite size scaling of the curves in fig. 1.3 is quite involved as the crossing point is moving. That is why the early numerical estimates of the threshold were very inaccurate. The work of Wilson [Wil02], moreover, showed that the experimental sizes are too small and the asymptotic regime for the critical exponent is not reached in any of the current empirical works. The study of XOR-SAT indeed shows a crossover in the critical exponent at sizes which are not accessible for $K$-SAT [LRTZ01].

The studies of random $K$-SAT opened up the exciting possibility to connect the hardness with an algorithm-independent property, like the satisfiability phase transition. But what exactly makes the instances near to the threshold hard remained an open question.

4 Statistical physics comes to the scene

Spin glass is one of the most interesting puzzles in statistical physics. An example of a spin glass material is a piece of gold with a small fraction of iron impurities. Physicist, on contrary to the rest of the human population, are interested in the behaviour of these iron impurities and not in the piece of gold itself. A new type of a phase transition was observed from the high temperature paramagnetic phase to the low temperature spin glass phase, where the magnetization of each impurity is frozen to a non-zero value, but there is no long range ordering. More than 30 years ago Edwards and Anderson [EA75] introduced a lattice model for such magnetic disordered alloys

where $S_{i}\in\{-1,+1\}$ are Ising spins on a 3-dimensional lattice, the sum runs over all the nearest neighbours, $h$ is the external magnetic field and the interaction $J_{ij}$ is random (usually Gaussian or randomly $\pm J$). The solution of the Edwards-Anderson model stays a largely open problem even today.

The mean field version of the Edwards-Anderson model was introduced by Sherrington and Kirkpatrick [SK75], the sum in the Hamiltonian (1.11) then runs over all pairs $(ij)$ as if the underlying lattice would be fully connected. Sherrington and Kirkpatrick called their paper ”Solvable Model of a Spin-Glass”. They were indeed right, but the correct solution came only five years later by Parisi [Par80c, Par80b, Par80a]. Parisi’s replica symmetry breaking (RSB) solution of the Sherrington-Kirkpatrick model gave rise to a whole new theory of the spin glass phase and of the ideal glass transition in structural glasses. The exactness of the Parisi’s solution was, however, in doubt till 2000 when Talagrand provided its rigorous proof [Tal06]. The relevance of the RSB picture for the original Edwards-Anderson model is widely discussed but still unknown.

A different mean field version of the Edwards-Anderson model was introduced by Viana and Bray [VB85], the lattice underlying the Hamiltonian (1.11) is then a random graph of fixed average connectivity. The complete solution of the Viana-Bray model is also still an open problem.

4.2 First encounter

The Viana-Bray model of spin glasses can also be viewed as random graph bi-partitioning (or bi-coloring at a finite temperature). The peculiarity of the spin glass phase will surely have some interesting consequences for the optimization problem itself. Indeed, the close connection between optimization problems and spin glass systems brought forward a whole collection of theoretical tools to analyze the structural properties of the optimization problems.

All started in 1985 when Mézard and Parisi realized that the replica theory can be used to solve the bipartite weighted matching problem [MP85]. Let us quote from the introduction of this work: ”This being a kind of pioneering paper, we have decided to present the method {meaning the replica method} on a rather simple problem (a polynomial one) the weighted matching. In this problem one is given $2N$ points $i=1,\dots,2N$, with a matrix of distance $l_{ij}$, and one looks for a matching between the points (a set of $N$ links between two points such that at each point one and only one link arrives) of a minimal length.” Using the replica symmetric (RS) approach they computed the average minimal length, when the elements of the matrix $l_{ij}$ are random identically distributed independent variables.

Shortly after Fu and Anderson [FA86] used the replica method to treat the graph bi-partitioning problem. They were the first to suggest that, possibly, the existence of a phase transition in the average behaviour will affect the actual implementation and performance of local optimization techniques, and that this may also play an important role in the complexity theory. Only later, such a behaviour was indeed discovered empirically by computer scientists [CKT91, MSL92].

The replica method also served to compute the average minimal cost in the random traveling salesmen problem [MP86a, MP86b]. Partitioning a dense random graph into more than two groups and the coloring problem of dense random graphs were discussed in [KS87]. Later some of the early results were confirmed rigorously, mainly those concerning the matching problem [Ald01, LW04]. All these early solved models are formulated on dense or even fully connected graph. Thus the replica method and where needed the replica symmetry breaking could be used in its original form. Another example of a ”fully connected” optimization problem which was solved with a statistical physics approach is the number partitioning problem [Mer98, Mer00].

And what about our customary random $K$-satisfiability, which is defined on a sparse graph? Monasson and Zecchina worked out the replica symmetric solution in [MZ96, MZ97]. It was immediately obvious that this solution is not exact as it largely overestimates the satisfiability threshold, the replica symmetry has to be broken in random $K$-SAT.

An interesting observation was made in [MZK+99b]: They defined the backbone of a formula as the set of variables which take the same value in all the ground-state configurations In CSPs with a discrete symmetry, e.g. graph coloring, this symmetry has to be taken into account in the definition of the backbone.. No extensive backbone can exist in the satisfiable phase in the limit of large $N$. If it would, then adding an infinitesimal fraction of constraints would almost surely cause a contradiction. At the satisfiability threshold an extensive backbone may appear. The authors of [MZK+99b] suggested that the problem is computationally hard if the backbone appears discontinuously and easy if it appears continuously. They supported this by replica symmetric solution of the SAT problem with mixed 2-clauses and 3-clauses, the so-called $2+p$-SAT. Even if the replica symmetric solution is not correct in random $K$-SAT and even if it overlooks many other important phenomena the concept of backbone is fruitful and we will discuss its generalization in chapter 4.

How to deal with the replica symmetry breaking on a sparse tree-like graph was an open question since 1985, when Viana and Bray [VB85] introduced their model. The solution came only in 2000 when Mézard and Parisi published their paper ”Bethe lattice spin glass revisited” [MP01]. They showed how to treat correctly and without approximations the first step of replica symmetry breaking (1RSB) and described how, in the same way, one can in principal deal with more steps of replica symmetry breaking, this extension is however numerically very difficult. But before explaining the 1RSB method we describe the general replica symmetric solutions. And illustrate its usefulness on the problem of counting matchings in graphs [ZDEB-1]. Only then we describe the main results of the 1RSB solution and illustrate the method in the 1-in-$K$ SAT problem [ZDEB-3]. After we list several ”loose ends” which appeared in this approach. Finally we summarize the main contribution of this thesis. This will be the departure point for the following part of this thesis which contains most of the original results.

5 The replica symmetric solution

The replica symmetric (RS) solution on a locally tree-like graph consists of two steps:

Compute the partition sum and all the other quantities of interest as if the graph would be a tree.

The replica symmetric assumption: Assume that the correlations induced by long loops decay fast enough, such that this tree solution is also correct on the only locally tree-like graph.

Equivalent names used in literature for the replica symmetric solution are Bethe-Peierls approximation (in particular in the earlier physics references) or belief propagation (in computer science or when using the iterative equation as an algorithm to estimate the marginal probabilities - magnetizations in physics). Both these conveniently abbreviate to BP.

Let $\phi_{a}(\partial a)$ be the evaluating function for the constraint $a$ depending on the variables neighbourhooding with $a$ in the factor graph $G(V,F,E)$. A satisfied constraint has $\phi_{a}(\partial a)=1$ and violated constraint $\phi_{a}(\partial a)=0$. The Hamiltonian can then be written as

The energy cost is thus one for every violated constraint. The corresponding Boltzmann measure on configurations is:

where $\beta$ is the inverse temperature and $Z_{G}(\beta)$ is the partition function. The marginals (magnetizations) $\chi^{i}_{s_{i}}$ are defined as the probabilities that the variable $i$ takes value $s_{i}$

The reason for this interest is that, for reasonable graph ensembles, $F_{G}(\beta)$ is self-averaging. This means that the distribution of $F_{G}(\beta)/N$ becomes more and more sharply peaked around $f(\beta)$ when $N$ increases.

5.2 The replica symmetric solution on a single graph

First suppose that the underlying factor graph is a tree, part of this tree is depicted in fig. 1.4. We define messages $\psi_{s_{i}}^{a\to i}$ as the probability that node $i$ takes value $s_{i}$ on a modified graph where all constraints around $i$ apart $a$ were deleted, and $\chi_{s_{j}}^{j\to a}$ as the probability that variable $j$ takes value $s_{j}$ on a modified graph obtained by deleting constraint $a$. On a tree these messages can be computed recursively as

where $Z^{a\to i}$ and $Z^{j\to a}$ are normalization constants, the factor $\phi_{a}(\{s\},\beta)=1$ if the constraint $a$ is satisfied by the configuration $\{s\}$ and $\phi_{a}(\{s\},\beta)=e^{-\beta}$ if not. We denote by $\psi^{a\to i}$ the whole vector $(\psi^{a\to i}_{0},\dots,\psi^{a\to i}_{q-1})$ and analogically $\chi^{j\to a}=(\chi^{k\to a}_{0},\dots,\chi^{j\to a}_{q-1})$. This is one form of the belief propagation (BP) equations [KFL01, Pea82], sometimes called sum-product equations. The probabilities $\psi$, $\chi$ are interpreted as messages (beliefs) living on the edges of the factor graph, with the consistency rules (1.16a) and (1.16b) on the function and variable nodes. Equations (1.16) are usually solved by iteration, the name message passing is used in this context. In the following it will be simpler not to consider the ”two-levels” equations (1.16) but

where $Z^{j\to i}=Z^{a\to i}\prod_{j\in\partial a-i}Z^{j\to a}$. Notice that on simple graphs, i.e., when either $l_{i}=2$ for all $i=1,\dots,N$ or $k_{a}=2$ for all $a=1,\dots,M$, the form (1.17) simplifies further. And on constraint satisfaction problems on simple graphs (e.g. the matching or coloring problems) the ”two-levels” equations are almost never used.

Assuming that one has found the fixed point of the belief propagation equations (1.16a-1.16b), one can deduce the various marginal probabilities and the free energy, entropy etc. The marginal probability (1.14) of variable $i$ estimated by the BP equations is

To compute the free energy we first define the free energy shift $\Delta F^{a+\partial a}$ after addition of a function node $a$ and all the variables $i$ around it, and the free energy shift $\Delta F^{i}$ after addition of a variable $i$. These are given in general by:

The total free energy is then obtained by summing over all constraints and subtracting the terms counted twice [MP01, YFW03]:

This form of the free energy is variational, i.e., the derivatives $\frac{\partial(\beta F_{G}(\beta))}{\partial\chi^{i\to a}}$ and $\frac{\partial(\beta F_{G}(\beta))}{\partial\psi^{a\to i}}$ vanish if and only if the probabilities $\chi^{i\to a}$ and $\psi^{a\to i}$ satisfy (1.16a-1.16b). This allows to compute easily the internal energy as

All the equations (1.16)-(1.22) are exact if the graph $G$ is a tree. The replica symmetric approach consists in assuming that all correlations decay fast enough that application of eqs. (1.16)-(1.22) on a large tree-like graph $G$ gives asymptotically exact results. These equations can be used either on a given graph $G$ or to compute the average over the graph (and disorder) ensemble.

5.3 Average over the graph ensemble

where the functions ${\cal F}_{\psi}$ and ${\cal F}_{\chi}$ represent the BP equations (1.16a-1.16b), $q(l)$ and $r(k)$ are the excess degree distributions defined in (1.8). If there is a disorder in the interaction terms, as e.g. the negations in $K$-SAT, we average over it at the same place as over the fluctuating degree.

Solving equations (1.23a-1.23b) to obtain the distributions ${\cal P}$ and ${\cal O}$ is not straightforward. In some cases (on regular factor graphs, at zero temperature, etc.) it can be argued that the distributions ${\cal P}$, ${\cal O}$ are sums of Dirac delta functions. Then the solution of eqs. (1.23a-1.23b) can be obtained analytically. But in general distributional equations of this type are not solvable analytically. However, a numerical technique called population dynamics [MP01] is very efficient for their resolution. In appendix E we give a pseudo-code describing how the population dynamics technique works.

Once the distributions ${\cal P}$ and ${\cal O}$ are known the average of the free energy density can be computed by averaging (1.20) over ${\cal P}$. This average expression for the free energy is again in its variational form (see [MP01]), i.e., the functional derivative $\frac{\delta f(\beta)}{\delta{\cal P}(h)}$ vanishes if and only if ${\cal P}$ satisfies (1.32). The average energy and entropy density are thus expressed again via the partial derivatives.

As we mentioned, on the ensemble of random regular factor graphs (without disorder in the interactions) the solution of equations (1.23) is very simple: ${\cal P}(\psi)=\delta(\psi-\psi^{\rm reg})$, ${\cal Q}(\chi)=\delta(\chi-\chi^{\rm reg})$, where $\psi^{\rm reg}$ and $\chi^{\rm reg}$ is a self-consistent solution of (1.16). This is because in the thermodynamical limit an infinite neighbourhood of every variable is exactly identical thus also the marginal probabilities have to be identical in every physical solution.

5.4 Application for counting matchings

To demonstrate how the replica symmetric method works to compute the entropy, that is the logarithm of the number of solutions, we review the results for matching on sparse random graphs [ZDEB-1]. The reasoning why the replica symmetric solution is exact for the matching problem is done on the level of self-consistency checks in [ZDEB-1]. And [BN06] have worked out a rigorous proof for graphs with bounded degree and a large girth (length of the smallest loop).

Consider a graph $G(V,E)$ with $N$ vertices ($N=|V|$) and a set of edges $E$. A matching (dimerization) of $G$ is a subset of edges $M\subseteq E$ such that each vertex is incident with at most one edge in $M$. In other words the edges in the matching $M$ do not touch each other. The size of the matching, $|M|$, is the number of edges in $M$. Our goal is to compute the entropy of matchings of a given size on a typical large Erdős-Rényi random graph.

We describe a matching by the variables $s_{i}=s_{(ab)}\in\{0,1\}$ assigned to each edge $i=(ab)$ of $G$, with $s_{i}=1$ if $i\in M$ and $s_{i}=0$ otherwise. The constraints that two edges in a matching cannot touch impose that, on each vertex $a\in V$: $\sum_{b,(ab)\in E}s_{(ab)}\leq 1$. To complete our statistical physics description, we define for each given graph $G$ an energy (or cost) function which gives, for each matching $M=\{s\}$, the number of unmatched vertices:

where $E_{a}=1-\sum_{\partial b}s_{(ab)}$.

In the factor graph representation we transform the graph $G$ into a factor graph $F(G)$ as follows (see fig. 1.5): To each edge of $G$ corresponds a variable node (circle) in $F(G)$; to each vertex of $G$ corresponds a function node (square) in $F(G)$. We shall index the variable nodes by indices $i,j,k,\dots$ and function nodes by $a,b,c,\dots$. The variable $i$ takes value $s_{i}=1$ if the corresponding edge is in the matching, and $s_{i}=0$ if it is not. The weight of a function node $a$ is

where $\partial a$ is the set of all the variable nodes which are neighbours of the function node $a$, and the total Boltzmann weight of a configuration is $\frac{1}{Z_{G}(\beta)}\prod_{a}\phi_{a}(\{\partial a\},\beta)$.

The belief propagation equation (1.16) becomes

where $Z^{b\to a}$ is a normalization constant. In statistical physics the more common form of the BP equations uses analog of local magnetic fields instead of probabilities. For every edge between a variable $i$ and a function node $a$, we define a cavity field $h^{i\to a}$ as

The recursion relation between cavity fields is then:

The expectation value (with respect to the Boltzmann distribution) of the occupation number $s_{i}$ of a given edge $i=(ab)$ is equal to

The free energy shifts needed to compute the total free energy (1.20) are

The energy, related to the size of the matching via (1.24), is then

This is the sum of the probabilities that node $a$ is not matched.

The distributional equation (1.23) becomes

And the average free energy is explicitly

Where ${\cal R}(k)$ is the connectivity distribution of the function nodes, that is the connectivity distribution of the original graph, $c$ is the average connectivity. The distributional equations are solved via the population dynamics method, see appendix E. Fig. 1.6 then presents the resulting average entropy as a function of size of the matching.

6 Clustering and Survey propagation

As we said previously in the random $K$-SAT the replica symmetric solution is not generically correct. Mézard and Parisi [MP01] understood how to deal properly and without approximations with the replica symmetry breaking on random sparse graphs, that is how to take into account the correlations induced by long loops. More precisely in their approach only the one-step (at most two-step on the regular graphs) replica symmetry breaking solution is numerically feasible. Anyhow, such a progress opened the door to a better understanding of the optimization problems on sparse graphs. The $K$-satisfiability played again the prominent role.

To compute the ground state energy within the 1RSB approach we can restrict only to energetic considerations as described in [MP03], we call this approach the energetic zero temperature limit. Applying this method to $K$-satisfiability leads to several outstanding results [MPZ02, MZ02], we describe the three most remarkable ones. Soon after, analog results were obtained for many other optimization problems, for example graph coloring [MPWZ02, BMP+03, KPW04], vertex cover [Zho03], bicoloring of hyper-graphs [CNRTZ03], XOR-SAT [FLRTZ01, MRTZ03] or lattice glass models [BM02, RBMM04].

It was known already in the ”pre-1RSB-cavity era” that replica symmetry broken solution is needed to solve random $K$-SAT. Such a need is interpreted as the existence of many metastable well-separated states, in the case of highly degenerate ground state this leads to a clustering of solutions in the satisfiable phase [BMW00, MPZ02, MZ02]. The energetic 1RSB cavity method deals with clusters containing frozen variables (clusters with backbones), that is variables which have the same value in all the solutions in the cluster. It predicts how many of such clusters exist at a given energy, the logarithm of this number divided by the system size $N$ defines the complexity function $\Sigma(E)$. According to the energetic cavity method for 3-SAT, clusters exist, $\Sigma(0)\neq 0$, for constraint density $\alpha>\alpha_{\rm SP}=3.92$ [MPZ02, MZ02].

It was conjectured [MPZ02, MZ02] that there is a link between clustering, ergodicity breaking, existence of many metastable states and the difficulty of finding a ground state via local algorithms. The critical value $\alpha_{\rm SP}$ was called the dynamical transition and the region of $\alpha>\alpha_{\rm SP}$ the hard-SAT phase.

Clusters were viewed as a kind of pure states, however, in the view of many a good formal definition was missing. It was also often referred to some sort of geometrical separation between different clusters. A particularly popular one is the following: Clusters are connected components in the graph where solutions are the nodes and two solutions are adjacent if they differ in only $d$ variables. Depending on the model and author the value of $d$ is either one of $d$ is a finite number of $d$ is said to be any sub-extensive number. The notion of $x$-satisfiability, the existence of pairs of solutions at a distance $x$, leads to a rigorous proof of existence of exponentially many geometrically separated clusters [MMZ05, DMMZ08, ART06].

The energetic 1RSB cavity method allows to compute the ground state energy and thus also the satisfiability threshold $\alpha_{s}$. In 3-SAT its value is $\alpha_{s}=4.2667$ [MPZ02, MZ02, MMZ06]. This value is computed as a solution of a closed distributional equation. This time there is an excellent agreement with the empirical estimations. Is the one step of replica symmetry breaking sufficient to locate exactly the satisfiability threshold? The stability of the 1RSB solution was investigated in [MPRT04], the 1RSB energetic cavity was shown to describe correctly the ground state energy for $4.15<\alpha<4.39$ in 3-SAT. In particular, it yields the conjecture that the location of the satisfiability threshold is actually exact. From a rigorous point of view it was proven that the 1RSB equations give an upper bound on the satisfiability threshold [FL03, FLT03, PT04].

The most spectacular result was the development of a new message passing algorithm, the survey propagation [MZ02, BMZ05]. Before the replica and cavity analysis were used to compute the quenched averages of thermodynamical quantities. Using always the self-averaging property that the average of certain (not all) quantities is equal to their value on a large given sample. Mézard and Zecchina applied the energetic 1RSB cavity equations, later called survey propagation, on a single large graph. This resulted in an algorithm which is arguably still the best known for large instances of random 3-SAT near to the satisfiability threshold. And even more interesting than its performance is the conceptual advance this brought into applications of statistical physics to optimization problems.

7 Energetic 1RSB solution

In this section we derive the energetic zero-temperature limit of the 1RSB method. When applied to the satisfiability problem this leads, between others, to the calculation of the satisfiability threshold and to the survey propagation equations and algorithm. We illustrate this on the 1-in-3 SAT problem. Before doing so we have to introduce the warning propagation equations, on which the derivation of the survey propagation relies.

In general warning propagation (min-sum) is a zero temperature, $\beta\to\infty$, limit of the belief propagation (sum-product) equations (1.16a-1.16b). It can be used to compute the ground state energy (minimal fraction of violated constraints) at the replica symmetric level. A constraint satisfaction problem at a finite temperature gives rise to $\phi_{a}(\{\partial a\},\beta)=1$ if the constraint $a$ is satisfied by configuration $\{s_{\partial a}\}$, and $\phi_{a}(\{\partial a\},\beta)=e^{-2\beta}$ if $a$ is not satisfied by $\{s_{\partial a}\}$The factor 2 in the Hamiltonian is introduced for convenience and in agreement with the notation of [ZDEB-3].. In a general Boolean CSP, with $N$ variables $s_{i}\in\{-1,1\}$, the warning propagation can then be obtained from (1.16a-1.16b) by introducing warnings $u$ and $h$ as

This leads in the limit of zero temperature, $\beta\to\infty$, to

The warnings $u$ and $h$ can thus be interpreted in the following way

Given this interpretation the prescriptions (1.35) on how to update the warnings over the graph becomes intuitive. Variable $i$ collects the preferences from all constraints except $a$ and sends the result to $a$. Constraint $a$ then decides which value $i$ should take given the preferences of all its other neighbours.

Given the fixed point of the warning propagation (1.35) the total warning of variable $i$ is

The corresponding energy can be computed as

where $\Delta E^{a+\partial a}$ is the number of contradictions created when constraint $a$ and all its neighbours are added to the graph, $\Delta E^{i}$ is the number of contradictions created when variables $i$ is added to the graph. The energy shifts can be computed from (1.19a-1.19b) using (1.34) and taking $\beta\to\infty$ they read

To summarize, the warning propagation equations neglect every entropic information in the belief propagation (1.16a-1.16b), thus only the ground state energy can be computed. On the other hand the fact that warnings $u$ and $h$ have a discrete set of possible values simplifies considerably the average over the graph ensemble presented in sec. 1.5.3 as the distribution ${\cal P}$ is a sum of three Dirac function, and can be represented by their weights. Deeper interpretations of warning propagation and its fixed points will be given in chapter 4. Note that in the literature the value of warnings is also called $\ast$ or ”joker” [BMWZ03, BZ04].

7.2 Survey Propagation

Survey propagation (SP) [MPZ02, MZ02] is a form of belief propagation which aims to count the logarithm of the number of fixed points of warning propagation (1.35) of a given energy (1.38). For the sake of simplicity we present the most basic form of SP which aims to count the logarithm of number of fixed points of the warning propagation with zero energy.

The constraints on values of the warnings assuring that the fixed point of warning propagation corresponds to zero energy are

For all $i$ and $a\in\partial i$: the warnings $\{u^{b\to i}\}_{b\in\partial i-a}$ are all non-negative or all non-positive,

For all $a$ and $i\in\partial a$: the preferred values of all $j\in\partial a-i$ can be realized without violating the constraint $a$.

We define probabilities that warnings $u^{a\to i}$ or $h^{i\to a}$ are positive, negative or null.

where ${\cal N}_{i\to a}$ is the normalization factor. The update of surveys $q$ given the incoming $p$s depends on the details on the constraint functions. For concreteness we write the equation for the positive 1-in-3 SAT problem. The constraints assuring zero energy then forbids that both the warnings incoming to a constraint $a$ have value $+1$.

where ${\cal N}_{a\to i}=1-p^{j\to a}_{+}p^{k\to a}_{+}$ is the normalization factor, $j$ and $k$ are the other two neighbours of $a$.

The associated Shannon entropy is called complexity [Pal83] (or structural entropy in the context of glasses) and reads [MZ02]

where ${\cal N}^{a+\partial a}$ is the probability that no contradiction is created when the constraint $a$ and all its neighbours are added, ${\cal N}^{i}$ is the probability that no contradiction is created when the variable $i$ is added. Remark the exact analogy with (1.19a-1.19b). We denote ${\cal P}_{0}^{i}\equiv\prod_{a\in\partial i}q_{0}^{a\to i}$ and ${\cal P}_{\pm}^{i}\equiv\prod_{a\in\partial i}(q_{\pm}^{a\to i}+q_{0}^{a\to i})$, then

The second equation collects the contributions from all combinations of arriving surveys except the “contradictory” ones $(+,+,+)$, $(-,-,-)$, $(+,+,0)$ and $(+,+,-)$ (plus permutations of the latter).

The survey propagation equations (1.41-1.42) and the expression for the complexity function (1.43) are exact on tree graphs. In the spirit of the Bethe approximation, we will assume sufficient decay of correlations and use these equations on a random graph The fact that on a given tree with given boundary conditions the warning propagation has a unique fixed point might seem puzzling at this point. Clarification will be made in the chapter 2.. To average over the ensemble of random graphs we adopt the same equations as we did for the belief propagation in sec. 1.5.3.

7.3 Application to the exact cover (positive 1-in-3 SAT)

The 1-in-3 SAT problem (with probability of negating a variable equal to one-half) is a rare example of an NP-complete problem which is on average algorithmically easy and where the threshold can be computed rigorously [ACIM01]. In particular it was shown that for $\alpha\neq 1$ an instance of the problem can be solved in polynomial time with probability going to one as $N\to\infty$. This result was generalized into random 1-in-3 SAT where the probability of negating a variable is $p\neq 1/2$ [ZDEB-3]. In particular we showed that for all $0.273<p<0.718$ the RS solution is correct and almost every instance can be solved in polynomial time if the constraint density $\alpha\neq 1/[4p(1-p)]$. When, however, $p<0.273$ the phase diagram is more complicated, see [ZDEB-3]. For $p=0$ the solution of the positive 1-in-3 SAT (exact cover) problem becomes very similar to the one of 3-SAT [MZ02]. The result for the complexity (1.43) in the positive 1-in-3 SAT obtained from the population dynamics method is plotted in fig. 1.7. For more detailed discussion of how the phase diagram changes from the almost-always-easy to the very-hard pattern see [ZDEB-3].

Up to certain average connectivity of variables $c_{\rm SP}=1.822$ the only iterative fixed point of the population dynamics gives $q_{0}^{a\to i}=p_{0}^{i\to a}=1$ for all $(ia)$. The associated complexity function is zero. In an interval $(c_{\rm SP},c_{s})=(1.822,1.879)$ there exist a nontrivial solution giving positive complexity function. There are thus exponentially many different fixed points of the warning propagation. Asymptotically, almost every warning propagation fixed point is associated to a cluster of solutionsThere might exist fixed points of the warning propagation which are not compatible with any solution, thus do not correspond to a cluster. Such ”fake” fixed points are negligible if the 1RSB approach is correct.. Above $c_{s}=1.879$ there is a nontrivial solution to the SP equations giving a negative complexity function. There are thus almost surely no nontrivial fixed points of warning propagation at zero energy.

Before interpreting the survey propagation results, we should check that its application on tree-like random graphs is justified. The method to do this self-consistency check has been developed in [MPRT04] and is discussed in appendix D. For 1-in-3 SAT the result in that SP is stable, thus the results are believed to be correct, for $c\in(1.838,1.948)$ [ZDEB-3]. The point $c_{s}$ belongs to this interval, thus we can interpret it safely as the satisfiability threshold. However, the point $c_{\rm SP}$ has no physical meaning, and some statements that are suggested by its existence are wrong. For example it is not true that there is not exponentially many fixed points of the warning propagation, thus no clustering, for $c<c_{\rm SP}$. This has been remarked in [ZDEB-4] and a part of chapter 2 will be devoted to understanding this.

8 Loose ends

We could summarize the understanding of the subject three years ago in the following way: The 1RSB cavity method was able to compute the satisfiability threshold. The clustered phase was predicted and its existence partially proven. The conjecture that clustering is a key element in understanding of the computational hardness was accepted. The survey propagation inspired decimation algorithm was breath-taking, and the computer science community was getting gradually more and more interested in the concepts which lead to its derivation. It might have seemed that a real progress can be made only on the mathematical side of the theory, in the analytical analysis of the performance of the message passing algorithms, or in new applications. But several loose ends hanged in the air and the opinions on their resolution were diverse. I will list three of them which I consider to be the most obtruding ones.

The energetic 1RSB cavity method (survey propagation) predicts the clustering in 3-SAT at $\alpha_{\rm SP}=3.92$. But the replica symmetric solution is unstable at already $\alpha_{\rm RS}=3.86$, at this point the spin glass susceptibility diverges and equivalently the belief propagation algorithm stops to converge on a single graph, see appendix C. What is the solution in the ”no man’s land” between $\alpha_{\rm RS}$ and $\alpha_{\rm SP}$? The values are even more significant for the 3-coloring or Erdős-Rényi graphs where the corresponding average connectivities are $c_{\rm RS}=4$ and $c_{\rm SP}=4.42$.

An iterative procedure called whitening of a solution is defined as iteration of the warning propagation equations initialized from a solution. Whitening core is the corresponding fixed point. We call white those variables which are assigned the ”I do not care” state in the whitening core. A crucial asymptotic property is that if the 1RSB solution is correct then the whitening core of all solutions from one cluster is the same and the non-white variables are the frozen ones in that cluster. Consequently, knowing a solution, the whitening may be used to tell if the solution was or was not in a frozen cluster.

Survey propagation uses information only about frozen cluster. It might seem that every cluster is uniquely described by its whitening core, that is by the set and values of the frozen variables.

Yet, the solutions found by survey propagation have always a trivial, all white, whitening core. This paradox was pointed out in [MMW07] and observed also by the authors of [BZ04]. It was suggested that the concept of whitening might be meaningful only in the thermodynamical limit. But that was not a satisfactory explanation.

The clustered phase, baptized ”Hard” in [MZ02] does not seem to be that hard. There is no local algorithm which would perform well exactly up to $\alpha_{\rm SP}=3.92$. For a while it was thought that the 1RSB stability point $\alpha_{II}=4.15$, see appendix D , is a better alternative. It was argued that the full-RSB states are more ”transparent” for the dynamics than the 1RSB states which should be well defined and separated. Moreover there was at least one empirical result which suggested that the Walk-SAT algorithm stops to work in linear time at that point [AGK04]. But other version of Walk-SAT stopped before or even after, as for example the ASAT which was argued in [AA06] to work in linear time at least up to $\alpha=4.21$.

9 Summary of my contributions to the field

In my first works [ZDEB-1, ZDEB-2, ZDEB-3] I applied the replica symmetric and the energetic 1RSB method to the matching and the 1-in-$K$ SAT problems. This is why I used these two problems to illustrate the methods in sec. 1.5.4 and 1.7.

The problem of matching on graphs is a common playground for algorithmic and methodological development. I studied the problem of counting maximum matchings in a random graph in [ZDEB-1]. Finding a maximum matching is a well known polynomial problem, while their approximative counting is a much more difficult task. We showed, that the entropy of maximum matchings can be computed using the belief propagation algorithm, a result which was later on partially proved rigorously [BN06].

My interest in the 1-in-$K$ SAT problem stemmed from the work [ACIM01] where the authors computed rigorously the satisfiability threshold and showed that the NP-complete problem is in fact on average algorithmically easy. In [ZDEB-2, ZDEB-3] we studied the random 1-in-3 SAT in two-parameter space. One parameter is the classical constraint density, the other is the probability $p$ of negating a variable in a constraint ($p=1/2$ in [ACIM01]). We showed that for $0.2627<p<0.7373$ the problem is on average easy and the satisfiability threshold can be computed rigorously. On the other hand for $p<0.07$ the problem is qualitatively similar to the 3-SAT. We computed the threshold from the energetic 1RSB approach. In the intermediate region the 1RSB approach is not stable, thus it stays an open question how exactly does the problem evolve from an on average easy case to a 3-SAT like case. Qualitatively similar phase diagram was described in the $2+p$ SAT problem [MZK+99a, AKKK01]. We also found an interesting region of the parameter space in the 1-in-3 SAT where the unit clause algorithm provably finds solutions despite the replica symmetric solution being not correct (unstable).

The rest of my works [ZDEB-4, ZDEB-5, ZDEB-6, ZDEB-7, ZDEB-8, ZDEB-10, ZDEB-9] tied up the loose ends from the previous section and mainly addressed the original question of this thesis: Why are some constraint satisfaction problems intrinsically hard on average and what causes this hardness?

I used the entropic zero temperature 1RSB approach, introduced in [MPR05], to study the structure of solutions in random CSPs. In [ZDEB-4, ZDEB-5] we discovered that the true clustering (dynamical) transition does not correspond to the onset of a nontrivial solution of the survey propagation equations. We gave a proper definition of the clustering transition and formulated it in terms of extremality of the uniform measure over solutions. The clustering transition happens always before or at the same time as the replica symmetric solution ceases to be stable. This tied up the loose end (A), as in the ”no man’s land” the energetic 1RSB solution was simply incomplete.

We showed that in general there exist two distinct clustered phases below the satisfiable threshold. In the first, dynamic clustered phase, an exponentially large number of pure states is needed to cover almost all solutions. However, average properties (such as total entropy) still behave as if the splitting of the measure did not count. In particular, a simple algorithm such as belief propagation gives asymptotically correct estimates of the marginal probabilities. However, the measure over solutions is not extremal and, more importantly, the Monte Carlo equilibration time diverges, thus making the sampling of solutions a hard problem. The second kind of clustered phase is the condensed clustered phase where a finite number of pure states is sufficient to cover almost all solutions. A number of nontrivial predictions follows: for instance the total entropy has a non-analyticity at the transition to this phase, the marginal probabilities are non-self-averaging and not given anymore by the belief propagation algorithm.

In the context of the coloring problem, i.e. anti-ferromagnetic Potts glass, I also addressed related questions of what does the 1RSB solution predict for the finite temperature phase diagram and when is the 1RSB solutions correct (stable) [ZDEB-5]. We give the full phase diagram for this model and argue that in the colorable phase for at least 4 colors the 1RSB solutions is stable, and thus believed to be exact.

In order to clarify and substantiate this heuristic picture, we introduced the random subcubes model in [ZDEB-8], a generalization of the random energy model. The random subcubes model is exactly solvable and reproduces the sequence of phase transitions in the real CSPs (clustering, condensation, satisfiability threshold). Its, perhaps, most remarkable property is that it reproduces quantitatively the behaviour of random $q$-coloring and random $K$-SAT in the limit of large $q$ and $K$. We showed that the random subcubes model can also be used as a simple playground for the studies of dynamics in glassy systems.

An important and quite novel phenomena I investigated in [ZDEB-5, ZDEB-7] is the freezing of variables. A variable is frozen when in all the solutions belonging to one cluster it takes the same value. I discovered that the fraction of such frozen variables undergoes a first order phase transition when the size of states is varied. I introduced the notion of the rigidity transition as the point where almost all the dominating clusters become frozen and the freezing transition as the point where all the clusters become frozen. The solutions belonging to the frozen clusters can be recognized via the whitening procedure.

We computed the rigidity transition in the random coloring in [ZDEB-5]. And we studied the freezing transition in 3-SAT numerically [ZDEB-10], with the result $\alpha_{f}=4.254\pm 0.009$ (to be compared to the satisfiability threshold $\alpha_{s}=4.267$). This study also confirms that the notion of whitening and freezing of variables in meaningful even on relatively small systems.

This allows us to tie up the loose end (B). The survey propagation algorithm describes the most numerous frozen clusters. The range of connectivities where the SP based algorithms are able to find solutions in 3-SAT lies in the phase where most solutions are in fact unfrozen. It is thus much less surprising that the SP based algorithms always find a solution with a trivial whitening.

A very natural question cannot be avoided at this point: What happens in the frozen phase where all the solutions are frozen? We know that such a phase exists, this was shown in [ART06] and numerically in [ZDEB-10]. And we also know from several authors that the known algorithms do not seem to be able to find frozen solutions in polynomial time (that is never for sufficiently large instances). We conjectured in [ZDEB-5] that the freezing is actually a relevant concept for the algorithmical hardness. Thus the answer we suggest to tie up the loose end (C) is that the simple local algorithms stop always before the freezing transition. It is a challenging problem to design an algorithm which would be able to beat this threshold.

In the coloring and satisfiability problems (at reasonably small $q$ and $K$) the freezing transition is however very near to the satisfiability threshold, see the numbers in [ZDEB-5, ZDEB-10]. It is thus difficult to make strong empirical conclusions about the relation between hardness and freezing. Motivated by the need of problems where the freezing and satisfiability would be well separated I introduced the locked constraint satisfaction problems where the freezing transition coincides with the clustering one [ZDEB-9]. The locked CSPs are very interesting from several points of view. The clusters in locked CSPs are point-like, this is why the clustering and freezing coincide. This is also connected with a remarkable technical simplification, as these problems can be fully described on the replica symmetric level.

On the other hand the locked problems are extremely algorithmically challenging. We implemented the best known solvers and showed that they do not find solutions starting very precisely from the clustering (= freezing) transition. At the same time this transition is very well separated from the satisfiability threshold.

A remarkable point about a subclass of the locked problems which we called balanced is that the satisfiability threshold can be obtained exactly from the first and second moment calculation. This adds a huge class of constraint satisfaction problems to a handful of other NP-complete CSPs where the threshold is known rigorously. And it also brings the understanding of which properties of the problem introduce fluctuations which make the second moment method fail.

The numerical work on the 3-SAT problems [ZDEB-10] also addresses another important and almost untouched question: How much are the asymptotic results relevant for systems of practical sizes. We counted the number of clusters in random 3-SAT on instances up to size $N=150$ and compared to the analytical prediction. We saw that the comparison is strikingly good for already so small systems. This should encourage the application of statistical physics methods to the real world problems.

Chapter 2 Clustering

In this chapter we introduce the concept of clustering of solutions. First we investigate when does the replica symmetric solution fail. Then we derive the one-step replica symmetry breaking equations on trees and give their interpretation on random graphs. We discuss how several geometrical definitions of clusters might be related to the pure states and review the properties of the clustered phase. Finally, we revise how is the clustering related to the algorithmical hardness and conclude that it is considerably less than previously anticipated. The original contributions to this chapter were published in [ZDEB-4, ZDEB-5, ZDEB-10].

How to recognize when is the replica symmetric solution correct? First we have to explain what do we precisely mean by ”being correct”. We obviously require that quantities like the free energy, energy, entropy, marginal probabilities (magnetizations) are asymptotically exact when computed in the replica symmetric approach. But this is not enough, as this is also satisfied in the phase which we will call later the clustered (dynamical) 1RSB phase.

A commonly used necessary condition for the validity of the RS solution is referred to as the local stability towards 1RSB. It consists in checking that the spin glass susceptibility does not diverge, or equivalently that the belief propagation algorithm converges on a large single graph, or in the probability theory this corresponds to the Kesten-Stigum condition [KS66a, KS66b]. These and other equivalent representations for the replica symmetric stability are discussed in detail in appendix C. If the replica symmetric solution is not stable then it predicts wrong free energy, entropy, correlation functions, etc. But the contrary is far from being true: even if stable, the RS solution might be wrong, and even unphysical (predicting negative entropies in discrete models, negative energies in models with strictly non-negative Hamiltonian function, or discontinuities in functions which physically have to be Lipschitzian).

It is tempting to say: The replica symmetric solution is correct if and only if the assumptions we used when deriving it are correct. In deriving the belief propagation (1.16) and the RS free energy (1.20) we used only one assumption: The neighbours of a variable $i$ are independent random variables, under the Boltzmann measure (1.13), when conditioned on the value of $i$. As we will see, this assumption is asymptotically correct also in the dynamical 1RSB phase, and thus the RS marginal probabilities, or the free energy function remain asymptotically exact in that phase.

We thus need a different definition for the ”RS correctness” which would determine whether the Boltzmann measure (1.13) can be asymptotically described as a single pure state, and whether the equilibration time of a local dynamics is linear in the system size. At the same time we do not want this definition to refer the RSB solution, because obviously we want to justify the need of the RSB solution by the failure of the RS solution.

A definition satisfying the above requirements appeared only recently [MM08, MS05, MS06c], and it can be written in several equivalent ways. From now on we say that the replica symmetric solution is correct if and only if one of the following is true.

The point-to-set correlations decay to zero.

Reconstruction on the underlying graph in not possible.

The uniform measure over solutions satisfies the extremality condition.

The 1RSB equations at $m=1$, initialized in a completely biased configuration, converge to a trivial fixed point.

In the rest of this section we explain these four statements, and show that they are indeed equivalent, and explain how do they correspond to the existence of a nontrivial 1RSB solution. We should mention that in the so-called locked constraint satisfaction problems this definition have to be slightly changed at zero temperature, we will discuss that in sec. 4.3. The transition from a phase where the RS solution is correct to a phase where it is not is called the clustering or the dynamical transition.

Our goal is to describe the structure of the set of solutions of a constraint satisfaction problem with $N$ variables. Let $\phi_{a}(\partial a)$ be the constraint function depending on variables $s_{i}\in\partial a$ involved in the constraint $a$, $\phi_{a}(\partial a)=1$ if the constraint is satisfied, $\phi_{a}(\partial a)=0$ if not. The uniform measure over all solutions can be written as

where $Z$ is the total number of solutions. The uniform measure over solutions is the zero temperature limit, $\beta\to\infty$, of the Boltzmann measure

The above expressions are valid on any given finite factor graph. The theory of Gibbs measures [Geo88] tries to formally define and describe the limiting object to which (2.1-2.2) converge in the thermodynamical limit, $N\to\infty$. A common way to build this theory is to ask: What is the measure induced in a finite volume $\Lambda$ when the boundary conditions are fixed? Roughly speaking, the good limiting objects, called the Gibbs measures or the pure states, are such that boundaries taken from the Gibbs measure induce the same measure inside the finite large volume $\Lambda$.

The Ising model on a 2D lattice gives an excellent example of how a phase transition is seen via Gibbs measures. Whereas in the high temperature paramagnetic phase the Gibbs measure is unique, in the ferromagnetic phase there are two extremal measures, one corresponding to the positive average magnetization, the other to the negative average magnetization. Indeed, if a boundary condition is chosen from one of these two then the correct magnetization will be induced in the bulk. In general the bulk in equilibrium can be described by a linear combination of these two extremal objects.

In the disordered models the situation might be much more complicated. Indeed the proper definition of the Gibbs measure in the Edwards-Anderson model (1.11) and other glassy models is a widely discussed but still an open problem [Bov06, Tal03, NS92].

The locally tree-like lattices, we are interested in here, are also peculiar from this point of view. The main difference is that in any reasonable definition of the boundary variables, the boundary has volume comparable to the volume of the interior. Thus again the usual theory of Gibbs measure implies very little. On the other hand the tree structure makes some considerations simpler. We will try to understand what sort of long range correlations might appear on the tree-like graphs by studying the tree graphs with general boundary conditions.

1.1 Properties and equations on trees

It is a well known fact that on arbitrary tree, with arbitrary boundary conditions, the belief propagation equations and the Bethe free energy are exact (the thermodynamical limit is not even needed here) [Pea88, KFL01, YFW00].

But what if the boundary conditions are chosen from a complicated measure? Then very little (if anything) is known in general. However, there is a way how to choose the boundary conditions such that the tree is then described by the one-step replica symmetry breaking equations. This is closely linked to the problem of reconstruction on trees, studied in mathematics [EKPS00, Mos01, Mos04]. The link with 1RSB was discovered by Mézard and Montanari [MM06a]. We chose to present the 1RSB equations in this new way, because it opens the door to further mathematical developments. For the original statistical physics derivation we refer to [MP00]. Another recent computer science-like derivation, which is based on the construction of a decorated constraint satisfaction problem and writing belief propagation on such a problem, in presented in [MM08, Mor07].

We explain the concept of reconstruction on trees [Mos04]. For simplicity we consider $q$-coloring on a rooted tree with constant branching factor $\gamma$ (sometimes also called the Cayley tree). A more general situation (with disorder, in the interaction or in the branching factor) is described in appendix A.

Create a rooted tree with branching $\gamma$ and with $L$ generations. An example of $\gamma=2$ and $L=8$ is in fig. 2.1. Assign a color $s_{0}$ to the root and broadcast over the edges towards the leaves of the tree in such a way that if a parent node $i$ was assigned color $s_{i}$ then each of its ancestors is assigned random one of the remaining $q-1$ colors.

At the end of this broadcasting, every node in the tree is assigned a color, and this assignment corresponds to a proper coloring (neighbours have different colors). Now in an imaginary experiment we forget the colors everywhere but on the leaves. The problem of reconstruction consists in deciding if there is any information left in the values on the leaves (and their correlation) about the original color $s_{0}$ of the root in the limit of infinite tree $L\to\infty$. If the answer is yes then we say that the reconstruction is possible, if the answer is no then the reconstruction is not possible.

Call $\{s\}_{l}$ the assignment of colors in the $l^{\rm th}$ generation of the tree. Consider formally the probability $\psi_{s_{0}}(\{s\}_{l})$ that a broadcasting process which finished at the configuration $\{s\}_{l}$ started from the color $s_{0}$ at the root. In other words, in what fraction of assignments in the interior of the tree (compatible with the boundary conditions $\{s\}_{l}$) is the color of the root $s_{0}$? Reconstruction is possible if and only if

Intuitively when the branching $\gamma$ is small and the number of colors large the information about the root will be lost very fast. If, on the contrary, the branching is large compared to the number of colors some information remains. A simple exercise is to analyze the so-called naive reconstruction algorithm [Sem08]. The naive reconstruction is possible if the probability that the leaves determine uniquely the root does not go to zero as the number of generation goes to infinity. We compute the probability $\eta$ that the far-away boundary is compatible with only one value of the root. Denote $\eta_{l}$ the probability that a variable in the $l^{\rm th}$ generation is directly implied conditioned on the value of its parent. The probability $\eta_{l-1}$ can be computed recursively as

The terms in this telescopic sum come from probabilities that number $r$ out of the $q-1$ colors are not present in the $\gamma$ descendants. In the last generation we know the colors by definition of the problem, thus $\eta_{\infty}=1$. If the iterative fixed point of (2.4) is positive then the reconstruction is possible.

This simple upper bound on the branching $\gamma$ for which the reconstruction is possible is actually quite nontrivial and in the limit of large number of colors it coincides with the true threshold at least in the first two orders, see [ZDEB-5] and [Sem08, Sly08]. This upper bound is connected to the presence of frozen variables and will be discussed in a greater detail in chapter 4.

The iterative equations for the reconstruction problem are equivalent to the one-step replica symmetry breaking equations with Parisi parameter $m$, $m=1$ will apply to the original question of reconstructibility. This was first derived by Mézard and Montanari [MM06a] and it has some deep consequences for the understanding of the RSB solution. We now explain this derivation, still for the coloring problem with a fixed branching $\gamma$ and $q$ colors. A more general form is presented in appendix A.

For given boundary conditions $\{s\}_{l}$, constructed as described above, we compute the probability $\psi^{i\to j}_{s_{i}}$ (over all broadcasting experiments leading to these boundary conditions) that a variables $i$ had color $s_{i}$, where $j$ is the parent of $i$ and the edge $(ij)$ has been cut. Given the probabilities on the descendants of $i$, which are indexed by $k=1,\dots,\gamma$, we can write

because the descendants can take any other color but $s_{i}$. The $Z^{i\to j}$ is a normalization constant. It should be noticed that this is in fact the belief propagation equation (1.16) for the graph coloring. This equation can also be derived by counting how many assignments are consistent with the boundary conditions $\{s\}_{l}$. This gives a natural interpretation to $Z^{i\to j}$

where $Z^{(i)}$ is the total number of solutions consistent with $\{s\}_{l}$ if $i$ were the root. Thus $Z^{i\to j}$ is a change in the number of solutions compatible with the boundary conditions when the $\gamma$ branches are merged.

Now we consider the distribution over all possible boundary conditions which are achievable by the broadcasting process defined above. We have to specify the probability distribution on the boundary conditions. We consider that the probability of every boundary conditions $\{s\}_{l}$ is proportional to the power $m$ of the number of ways by which we could create $\{s\}_{l}$, denote this number $Z(\{s\}_{l})$. In other words, the probability of a given boundary condition is proportional to the power $m$ of the number of possible assignments in the bulk of the tree.

The value of $m=1$ is natural for the original question of reconstruction, because every realization of the broadcasting experiment is then counted in a equiprobable way. We, however, introduced a general power $m$. The parameter $m$ will play a role of the Legendre parameter, changing its value focuses on boundary conditions compatible with a given number of assignments inside the tree.

Denote $P^{i\to j}(\psi^{i\to j})$ the distribution of $\psi^{i\to j}$, over the measure on the boundary conditions (2.7)

Where $Z^{(i)}(\{s\}_{l})$ is the number of solutions induced on the subtree rooted in vertex $i$, ${\cal Z}^{(i)}(m)$ is the corresponding normalization. To express the probability distribution $P^{i\to j}(\psi^{i\to j})$ as a function of $P^{k\to i}(\psi^{k\to i})$ we need that $\psi^{i\to j}={\cal F}(\{\psi^{k\to i}\})$, eq. (2.5). Moreover, $Z^{i\to j}$ is the increase in the total number of solutions after merging the branches rooted at $k=1,\dots,\gamma$ into one branch rooted at $i$. The distributional equation for $P$ is then

where ${\cal F}$ and $Z^{i\to j}$ are defined in (2.5), and ${\cal Z}^{i\to j}$ is a normalization constant equal to

where ${\cal Z}^{(i)}$ is the normalization from (2.8) if $i$ were the root. Notice that if we start from boundary conditions which are not compatible with any solution then the re-weighting $Z^{i\to j}=0$ at the merging where a contradiction is unavoidable. Initially at the leaves the colors of nodes are known. Call $\delta_{r}$ the $q$-component vector $\psi_{s_{i}}^{i\to j}=\delta(s_{i},r)$, then the initial distribution is just a sum of singletons

Denote $P_{0}(\psi)$ the distribution created from (2.11) after many iteration of (2.9) with $m=1$. The reconstruction is possible if and only if $P_{0}(\psi)$ is nontrivial, that is different from singleton on $\psi_{s_{i}}=1/q\,,\,\forall s_{i}$. We define the critical branching factor $\gamma_{d}$ in such a way that for $\gamma<\gamma_{d}$ the reconstruction is not possible, and for $\gamma\geq\gamma_{d}$ the reconstruction is possible. The critical values $\gamma_{d}=c_{d}-1$ for the coloring problem are reviewed in tab. 5.2.

If the reconstruction is not possible, then almost all (with respect to (2.7) at $m=1$) boundary conditions do not contain any information about the original color of the root. However, for rare boundary conditions this might be different. Obviously as long as $\gamma\geq q-1$ one can always construct boundary conditions which determine uniquely the value of the root (by assigning every of the $q-1$ colors to the descendants of every node). If $\gamma<q-1$ then this is no longer possible. And it was proven in [Jon02] that for $\gamma<q-1$ every boundary conditions lead to an expectation $1/q$ for every color on the root. If the reconstruction is possible, then different boundary conditions may lead to different expectations on the root.

The basic idea of the definition of clusters on a tree is the same as in the classical definition of a Gibbs measure [Geo88]. However, some more work is needed to make the following considerations rigorous. Define a $d$-neighbourhood of the root as all the nodes up to $d^{\rm th}$ generation, consider $1\ll d\ll l$. Consider the set ${\cal S}$ (resp. ${\cal S}^{\prime}$) of all assignments on the $d$-neighbourhood compatible with a given boundary condition $\{s\}_{l}$ (resp. $\{s^{\prime}\}_{l}$). Define two boundary conditions $\{s\}_{l}$ and $\{s^{\prime}\}_{l}$ as equivalent if the fraction of elements in which the two sets ${\cal S}$ and ${\cal S}^{\prime}$ differ goes to zero as $l,d\to\infty$. Clusters are then the equivalence classes in the limit $l\to\infty$, $d\to\infty$, $d\ll l$. The requirement $d\ll l$ comes from the fact that in $l-d$ iterations the equation (2.8) should converge to its iterative fixed point.

As we explained, more than one cluster exists as soon as the branching factor $\gamma\geq q-1$, but as long as the iterative fixed point of eq. (2.9) at $m=1$ is trivial all but one clusters are negligible because they contain an exponentially small fraction of solutions. Indeed, if the reconstruction is not possible it means that the information about the $d$-neighbourhood is almost surely lost at the $l^{\rm th}$ generation. Thus almost every broadcasting will lead to a boundary condition from the only relevant giant cluster.

Only for $\gamma\geq\gamma_{d}$, when the reconstruction start to be possible, the total weight of all solutions will be split into many clusters. In every of them the set of expectation values (beliefs) $\psi^{i\to j}$ will be different. This is related to another derivation of the 1RSB equations where the clusters of solutions on a given graph are identified with fixed points of the belief propagation equations [MM08, Mor07]. There are exponentially many (in the total number of variables $N$) initial conditions, it is also reasonable to expect that the number of clusters will be exponentially large in $N$.

The number of solutions compatible with a boundary condition $(\{s\}_{l})$ was denoted $Z(\{s\}_{l})$ in eq. (2.7). The associated entropy is then, due to interpretation of $Z^{i\to j}$ (2.6)

where the sum is over all the vertices $i$ in the tree, if $i$ is a leaf then $Z^{i\to j}=1$, if $i$ is the root that $j$ is a imaginary parent of the root. An intuition about this formula is the following: $\log Z^{i\to j}$ is the change in the entropy when the node $i$ and all edges $(ki)$, where $k$ are descendants of $i$, are added. Summing over all $i$ then creates the whole tree.

More commonly, we introduce also messages going from the parents to the descendants and write the expression for the entropy (2.12) in the equivalent Bethe form [YFW03]

where $\partial i$ are all the neighbours (descendants and the parent) of node $i$. The first sum in (2.13) goes over all the nodes in the tree, the root included, leaves have only one allowed color, thus eq. (2.14) changes correspondingly. Again the meaning of $\log{Z^{i+\partial i}}$ is the change in the entropy when node $i$ and his neighbourhooding edges are added, each edge is then counted twice, thus the shift in the entropy when an edge $(ij)$ is added, $\log{Z^{ij}}$, have to be subtracted.

We denote $\Phi(m)\equiv\log{{\cal Z}(m)}$ the thermodynamical potential associated to the measure (2.7). To avoid confusion with the real free energy, associated to the uniform measure over solutions (2.1), we call it the replicated free entropy. If a nonzero temperature is involved then $-\Phi(m)/(\beta m)$ is called the replicated free energy. The replicated free entropy on a tree can be expressed in totally analogous way as the entropy. From (2.10) we derive

which is usually written in the equivalent way

We denote $\Sigma(m)$ the Shannon entropy corresponding to measure on the boundary conditions (2.7), and we call it the complexity function.

where $S$ is the entropy averaged with respect to $\mu(\{s\}_{l})$

Thus the complexity can also be written as a function of the internal entropy via the Legendre transform of the replicated free entropy $\Phi(m)$

The reader familiar with the cavity approach surely recognized eqs. (2.9) and (2.15–2.20) as the 1RSB equations.

In the cavity method [MP01] the exponential of the complexity function $\Sigma(m)$ (2.18) counts the number of clusters corresponding to a given value of the parameter $m$, that is of a given entropy $S$ (2.19). Complexity defined on the full tree is never negative, as it is a Shannon entropy of a discrete random variable. The same is, of course true, about the entropy (2.12).

It is more interesting to consider the complexity (or the entropy) function $\Sigma_{d}(m)$ on the $d$-neighbourhood of the root. If the total number of generations of the tree is $l$ we take $1\ll d\ll l$. And moreover we require $l-d$ to be large enough, such that the distributional iterative equation (2.9) converges to its fixed point in less than $l-d$ iterations. The average complexity function on the $d$-neighbourhood can then be computed from this fixed point. And it can be both positive or negative. Its negative value then means that the number of clusters is decreasing as we are getting nearer to the root. Two important critical connectivities can be defined

$\gamma_{c}$: at which the complexity of the ”natural” clusters $\Sigma_{d}(m=1)$ becomes negative.

$\gamma_{s}$: at which the maximum of the complexity $\Sigma_{d}(m=0)$ becomes negative.

The connectivity $\gamma_{s}$ is the tree-analog of the satisfiability threshold. The connectivity $\gamma_{c}$ is the tree-analog of the condensation transition on random graphs, see chapter 3.

Strictly speaking, it is not known how to justify the interpretation of the complexity function as the counter of clusters in the derivation we just presented. In the original cavity derivation [MP01] or in the later derivations [MM08, Mor07] this point is well justified. We, however, find the purely tree derivation appealing for further progress on the mathematical side of the theory and that is why we have chosen to present this approach despite this current incompleteness.

1.2 Back to the sparse random graphs

We stress that the equations, derived in the previous section, are all exact on a given (even finite) tree and that we have not use any approximation. We were just describing boundary conditions correlated via (2.7). These, in nature recursive, equations are solved via the population dynamics technique, see the appendix E.

To come back to the sparse random graphs, which are only locally tree-like, we can consider equations (2.9-2.26) as an approximation on arbitrary graphs, just as we did with belief propagation. This leads to the one-step replica symmetry breaking (1RSB) approach. Note that on random graphs we will always speak about densities of the entropy, complexity or free-entropy etc. Thus on random graphs: instead of the entropy $S$ defined in (2.12) we consider $s=S/N$. The replicated free entropy $\Phi$ (2.15) and complexity $\Sigma$ (2.18) are also divided by the number of variables. We, however, denote them by the same symbol, as confusion is not possible.

Let us discuss once again, now from the random graph perspective, what are the correlations which make the replica symmetric approach fail. This will finally explain the definition of the replica symmetric solution being correct given at the beginning of this section.

The concept of the point-to-set correlations is common in the theory of glassy systems. Usually it is considered in the phenomenology of the real glassy systems on finite-dimensional lattices, see for example [BB04] and references therein. Here we restrict the discussion to properties relevant for the tree-like lattices.

Call $\overline{B}_{d}(i)$ all vertices of the graph which are at distance at least $d$ from $i$, define point-to-set correlation function as

where $\mu(\cdot)$ is the uniform measure over solutions (2.1), and the total variation distance of two probability distributions is defined as $||q-p||_{\rm TV}=\sum_{x}|q(x)-p(x)|/2$. The average point-to-set correlation is

The reconstruction on graphs is then defined via the decay of this correlation function. The reconstruction on tree-like graphs is not in general equivalent to the reconstruction on trees. Roughly said, it is not equivalent in the ferromagnetic models, e.g. the ferromagnetic Ising model, which spontaneously break some of the discrete symmetries. On the other hand on most of the frustrated models they are equivalent. A general condition, which might be very nontrivial to check, is given in [GM07].

If the point-to-set correlation function decays to zero, $\lim_{d\to\infty}C_{d}=0$, then almost every variable is independent of its far away neighbours. The replica symmetric approach then has to be asymptotically correct on locally tree-like lattices.

On the other hand if the point-to-set correlations do not decay to zero, then the far-away neighbours influence the value of the variable $i$. And the replica symmetric solution fails to give the correct picture of the properties of the model. The lack of decay of the point-to-set correlations is equivalent to the reconstruction on graphs, and is also equivalent to the existence of a nontrivial solution of the 1RSB equation (2.9) at $m=1$. This is also equivalent to the extremality condition for the uniform measure (2.1), which was used in definition of [ZDEB-4] and reads

where the external average is over quenched disorder (in interactions or connectivities).

The point-to-set correlations do not decay to zero for example in the low temperature phase of the ferromagnetic Ising model on a random graph. There it is sufficient to introduce the pure state ”up” and the pure state ”down” and within these pure states the point-to-set correlations will decay to zero again. On the frustrated models the situation is more complicated but the idea of the resolution is the same: If we manage to split the set of solutions into clusters (pure states) such that within each cluster the point-to-set correlations again decay, the situation is fixed. A statistical description of the properties of clusters can be obtained using the one-step replica symmetry breaking (1RSB) equations, derived in the previous section 2.1.1 and summarized in the next section 2.1.3.

However, the correlations might be more complicated and might not be captured fully by the 1RSB approach. In particular the 1RSB approach is correct if and only if the point-to-set correlation decay to zero within clusters and if the replica symmetric statistical description of clusters is correct. In appendix D we will discuss a necessary condition for the 1RSB approach being correct. In case the 1RSB approach does not fully describe the system further steps of replica symmetry breaking might provide a better approximation (that means splitting clusters into sub-clusters or aggregation of clusters) [MP00]. However, on the tree-like lattices, the exact solutions is not known in such cases.

In glasses, the clustering transition is usually studied at finite temperature and is called the dynamical transition. The clustered phase with $\Sigma(m=1)>0$ is called the dynamical 1RSB phase. This phase, where most of the static properties do not differ from the replica symmetric (liquid) ones, was first described and discuss in [KT87a, KT87b]. The dynamical transition is associated with a critical slowing down of the dynamical properties, e.g. the equilibration time is expected to diverge at this point. Note that such a purely dynamical phase transition is typical for mean-field models. In the finite dimensional glassy systems the barriers between a metastable and an equilibrium state are finite (independent of the system size). This is because the nucleation length might be large but have to be finite. Thus instead of a sharp dynamical transition in finite dimensional systems we observe only a crossover.

However, even at the mean field level, the exact dynamical description is known only in a few toy models, e.g. the spherical $p$-spin model [CK93] or the random subcubes model [ZDEB-8]. In general, the dynamical solution is only approximative, still many very interesting results were obtained. For a review see [BCKM98]. In the models on sparse random lattices even the approximation schemes are rather poor, see e.g. [SW04]. Thus the exact general relation between dynamics and the dynamical (clustering) transition is not known.

An important contribution in establishing the link between dynamics and the static solution on random graphs is [MS05, MS06b, MS06c] where the divergence of the point-to-set correlation length is linked with divergence of the equilibration time of the Glauber dynamics. This suggests that beyond the clustering transition the Monte Carlo sampling (or maybe even sampling in general) will be a hard task.

Note also that in the mathematical literature the Glauber dynamics is often studied. Many results exist about the so-called rapid mixing of the associated Markov chain [Sin93]. But the rapid mixing questions equilibration in polynomial time, whereas in physics the relevant time scale is linear. Moreover rapid mixing is defined as convergence to the equilibrium measure from any possible initial conditions, whereas in physics of glasses the notion of a typical initial condition should be used instead.

1.3 Compendium of the 1RSB cavity equations

We review the 1RSB equations on a general CSP. The order parameter is a probability distribution of the cavity field (BP message) $\psi^{a\to i}=(\psi^{a\to i}_{0},\dots,\psi^{a\to i}_{q-1})$. The self-consistent equation for $P^{a\to i}$ reads

where the function ${\cal F}(\{\psi^{b\to j}\})$ and the term $Z^{j\to i}$ are defined by the BP equation (1.17), ${\cal Z}^{j\to i}$ is a normalization constant.

The associated thermodynamical potential (2.15) is computed as

where the terms $Z^{a+\partial a}$ and $Z^{i}$ are the partition sum contributions defined in (1.19).

The logarithm of the number of states divided by the system size defines the complexity function $\Sigma$. Inversely the number of states is $e^{N\Sigma}$. At finite temperature the complexity of states with a given internal free energy is a Legendre transformation of the potential $\Phi(m)$

Useful relations between the free energy, complexity and potential $\Phi$ are

At zero energy, $E=0$, and zero temperature, $\beta\to\infty$, the free energy becomes entropy $-\beta f\to s$. Then the complexity is a function of the internal entropy of states and (2.26) becomes

This is called the entropic zero temperature limit. The internal entropy is expressed as

where $\Delta S^{a+\partial a}$ ($\Delta S^{i}$ resp.) is an internal entropy shift when the constraint $a$ and all its neighbour (the variable $i$ resp.) are added to the graph.

In the energetic zero temperature limit, described in sec. 1.6 for zero energy, the Parisi parameter $y=\beta m$ is kept constant, thus $m\to 0$. The free energy then converges to the energy, and (2.26) becomes

where the complexity is this time a function of the energy density $e$. The survey propagation equations generalized to nonzero $y$ are called the SP-$y$ equations.

Equations (2.24-2.31) are defined on a single instance of the constraint satisfaction problem. Averages ${\cal P}$ over the graph ensemble are obtained in a similar manner as in sec. 1.5.3 for the replica symmetric solution.

where in the sum over $\{l_{i}\}$, $i\in\{1,\dots,K-1\}$, and the functional ${\cal F}_{2}$ is defined by (2.24). Analogical expression holds for the average of the complexity or internal entropy. A general method to solve the equation (2.33) is the population of populations described in appendix E.5.

2 Geometrical definitions of clusters

Up to now we were describing clusters, i.e., partitions of the space of solutions, in a very abstract way which was defined only in the thermodynamical limit. We showed how to compute the number of clusters of a given size (internal entropy) (2.28), and we argued that the description makes sense if the point-to-set correlation (2.22) decays to zero within almost every cluster of that size. In this last sense clusters are what we would call in statistical physics pure equilibrium states.

On a very intuitive level, cluster are groups of nearby solutions which are in some sense separated from each other. Several geometrical definitions are used in the literature, we want to review the most common ones and state their relation to the definition above. We want to stress that it is not know whether any of the geometric definitions is equivalent to the description given above and used usually in the statistical physics literature.

First rigorous proofs of existence of an exponential number of clusters of solutions in the random $K$-SAT were based on the concept of $x$-satisfiability. Two solutions are at distance $x$ if they differ in exactly $xN$ variables. A formula is said $x$-satisfiable if there is a pair of solutions at distance $x$, and $x$-unsatisfiable if there is not.

Mora, Mézard and Zecchina [MMZ05, DMMZ08] managed to prove that for $K\geq 8$ and a constraint density $\alpha$ near enough to the satisfiability threshold the formulas are almost surely $x$-satisfiable for $x<x_{0}$, almost surely $x$-unsatisfiable for $x_{1}<x<x_{2}$, and almost surely $x$-satisfiable at $x_{3}<x<x_{4}$, where obviously $0<x_{0}<x_{1}<x_{2}<x_{3}<x_{4}<1$. This means that at least two well separated clusters of solutions exist. Proving that there is an exponentially smaller number of pairs of solutions at distances $x<x_{1}$ than at distances $x>x_{2}$ leads to the conclusion that an exponential number of well geometrically separated clusters exists [ART06].

However, the $x$-satisfiability gives too strong conditions of separability. This is illustrated for example in the XOR-SAT problem [MM06b]. It is still an open question if there is or not a gap in the $x$-satisfiability in the random 3-SAT near to the satisfiability threshold.

Another popular choice of a geometrical definition of clusters is that clusters are connected components in a graph where every solution is a vertex and solutions which differ in $d$ or less variables are connected. The distance $d$ is often said to be any sub-extensive (in the number of variables $N$) distance, that is $d=o(N)$. However, such a rule is not very practical for numerical investigations.

In $K$-SAT, in fact, $d=1$ seems to be a more reasonable choice. There are two reasons: First, clusters defined via $d=1$ have correct ”whitening” properties as we explain in the next paragraph. Second, we numerically investigated the complexity of $d=1$ connected-components clusters, fig. 2.2 right, and the agreement with the total number of clusters computed from (2.28) at $m=0$ is strikingly good. In particular, near to the satisfiability threshold $\alpha>4.15$, where the 1RSB result for the total complexity function is believed to be correct (stable) [MPRT04].

Formally, connected-components clusters have no reason to be equivalent to the notion of pure states. They are not able to reproduce purely entropic separation between clusters, which might exist in models like 3-SAT. However, fig. 2.2 suggests that there is more in this definition than it might seem at a first glance.

We define the whitening of a solution as iterations of the warning propagation equations (1.35) initialized in the solution. The fixed point is then called the whitening core. Note, that the whitening core is well defined in the sense that the fixed point of the warning propagation initialized in a solution does not depend on the order in which the messages were updated. A whitening core is called trivial if all the warning messages are , that is ”I do not care”.

The 1RSB equations at $m=0$, which give the total complexity function, can be derived as belief propagation counting of all possible whitening cores [MZ02, BZ04, MMW07]. Thus another reasonable definition of clusters is that two solutions belong to the same cluster if and only if their whitening core is identical. In fig. 2.2 left we plot numerically computed complexity of the whitening-core clusters compared to the complexity computed from (2.28) at $m=0$. The agreement is again good, in particular near to the satisfiability threshold, $\alpha>4.15$, where the SP gives a correct result.

The $d=1$ connected-components clusters share the property that all the solution from one clusters have the same whitening core. Proof: If this would not be true then there have to exist a pair of solutions which do not have the same whitening core but differ in only one variable, this is not possible because then the whitening could be started in that variable.

Note, however, that the definition of whitening-core clusters put all the solutions with a trivial whitening core into one cluster. This is not correct as, at least near to the clustering threshold, there are many pure states with a trivial whitening core. This is closely connected to the properties of frozen variables which will be discussed in chapter 4.

In order to obtain the data in fig. 2.2 we generate instances of the random 3-SAT problem with $N$ variables and $M$ clauses, constraint density is then $\alpha=M/N$. We count number of solutions in $A=999$ random instances and choose the median one where we count the number of connected-components and whitening-core clusters ${\cal S}$. This is repeated $B=1000$ times. The average complexity is then computed as $\Sigma=\sum_{i=1}^{B}\log{\cal S}_{i}/(BN)$, if the median instance was unsatisfiable then we count zero to the average, that is if all the $B$ instances are unsatisfiable then the complexity is zero. We do such a non-traditional sampling to avoid rare instances with very many solutions, which we would not be able to cluster.

3 Physical properties of the clustered phase

Let us give a summary of the properties of the clustered phase, also called the dynamical 1RSB phase. We describe only the situation when $\Sigma(m=1)>0$ (2.28), when the opposite is true the properties are completely different as we will discuss in the next chapter 3.

The complexity function computed from (2.28) is the log-number of clusters of a given internal entropy. If a solution is chosen uniformly at random it will almost surely belong to a cluster with entropy $s^{*}$ such that $\Sigma(s)+s$ is maximized in $s^{*}$, $\partial_{s}\Sigma(s^{*})=-1$, that is $m=1$. At $m=1$ the total entropy $\Sigma(s^{*})+s^{*}=\Phi(m=1)$. The replicated free entropy at $\Phi(m=1)$ is equal to the replica symmetric entropy. Thus the total entropy in the dynamical 1RSB phase is equal to the RS entropy. Also the marginal probabilities at $m=1$ are equal to the replica symmetric ones

Thus the clustering transition is not a phase transition in the Ehrenfest sense, because the thermodynamical potential, entropy in our case, in analytical at the transition.

The overlap (or here distance) distribution, which is often used to describe the spin glass phase, is also trivial and equal to the replica symmetric one in the dynamical 1RSB phase. Indeed, if exponentially many clusters are needed to cover almost all solutions, then the probability that two solutions happen to belong to the same cluster is zero.

The correlation function between two variables at a distance (shortest path in the graph) $d$ is defined as $\langle s_{i}s_{j}\rangle_{c}=||\mu(s_{i},s_{j})-\mu(s_{i})\mu(s_{j})||_{\rm TV}$. The variance of the overlap distribution, which is negligible compared to $1$ as we explained, can be expressed as $\sum_{i,j}\langle s_{i}s_{j}\rangle^{2}_{c}/N^{2}$, and thus the two-point correlation have to decay faster with distance than the number on vertices at that distance is growing. This means in particular that two neighbours of a node $i$ are independent if we condition on the value of $i$, this is again consistent with the fact that the belief propagation equations predict correct total entropy and marginal probabilities.

So far nothing is different form the replica symmetric phase. It is thus not straightforward to recognize the dynamical 1RSB phase based on the original replica computation. Presence of this phase was discovered and discussed in [KT87a, KT87b]. Later purely static methods were developed to identify this phase. The most remarkable is perhaps the $\epsilon$-coupling and the ”potential” of [FP95, FP97].

In our setting the main difference between the replica symmetric phase and the dynamical 1RSB phase is that in the later the point-to-set correlations do not decay to zero. Consequently the equilibration time of the local Monte Carlo dynamics diverges and Monte Carlo sampling becomes difficult [MS06b].

4 Is the clustered phase algorithmically hard?

Clustering has important implications for the dynamical behaviour. It slows down the equilibration and thus uniform sampling of solutions via local single spin flip Monte Carlo is not possible, or exponentially slow, beyond the dynamical threshold. But finding one solution is a much simpler problem than sampling.

In the 3-coloring of Erdős-Rényi graphs the clustering threshold is $c_{d}=4$, as at this point the spin glass susceptibility diverges, see appendix C. In the terms of the reconstruction problem the Kesten-Stigum [KS66a, KS66b] bound is sharp. On the other hand Achlioptas and Moore [AM03] proved that a simple heuristic algorithm is able to find a solution in average polynomial time up to at least $c=4.03$. This shows that the RSB phase is not necessarily hard.

A similar observation was made in the 1-in-3 SAT problem in [ZDEB-3]. There is a region in the values of the average density of constraints and the probability of negating a variable in a clause in which the replica symmetric solution is unstable and yet the unit clause propagation algorithm with the short clause heuristics was proven to find a solution in polynomial average time.

We should mention a common contra-argument; which is that in the above mentioned regions the 1RSB approach might not be correct, and the presumably full-RSB phase [Par80c] is more ”transparent” for the dynamics of algorithms, see e.g. [MRT04]. However, at least in the 3-coloring, the 1RSB approach seems to be correct in the interval in question, as we argue in appendix D.

There is a lot of numerical evidence that relatively simple single spin flip stochastic local search algorithms are able to find solutions in linear time deep in the clustered region. Examples of works where performance of such algorithms was analyzed are [KK07, SAO05, AA06, AAA+07]. In fig. 2.3 we give an example of performance of the ASAT algorithm [AA06] in 4-coloring of Erdős-Rényi random graphs [ZDEB-5]. The algorithm is described in appendix F.2.2. In the 4-coloring ASAT is able to find solutions in linear time beyond the clustering transition $c_{d}=8.35$

There is no paradox in the observations above. Quantitative statements are, however, difficult to make. Let us describe on an intuitive level the behaviour of an algorithm (dynamics) which satisfies the detailed balance condition and thus in infinite time samples uniformly from the uniform measure (2.1). We think for example about the simulated annealing [KGV83]. Above the dynamical temperature $T_{d}$ corresponding to an energy $E_{d}$ the point-to-set correlation function (2.22) decay fast and thus simulated annealing is able to reach the equilibrium. Below temperature $T_{d}$ this is not the case anymore and the dynamics is stuck for a very long time in one of the clusters, states. But the bottom of this state $E_{\rm bottom}$ lies lower than $E_{d}$, thus when lowering the temperature the average energy seen by the simulated annealing also decreases. If $E_{\rm bottom}=0$ then the algorithm will find a solution. It is not known how to compute $E_{\rm bottom}$ in general. Sometimes, far from the clustering transition, the iso-complexity approach [MRT04] gives a lower bound on $E_{\rm bottom}$. But in general, as far as we know, there is no argument saying $E_{\rm bottom}>0$. This picture can be substantiated for several simple models as the spherical $p$-spin model [CK93] or the random subcubes model [ZDEB-8]. The connection with the optimization problems was remarked in [KK07].

For the stochastic local search algorithm, which does not satisfy the detailed balanced condition, the situation might be similar. At a point the algorithm is stuck in a cluster, but if this cluster goes down to the zero energy then it might be able to find solutions even in the clustered phase.

However, the current understanding of the dynamics of the mean field glassy systems is far from complete. More studies are needed to understand better the link between the static clustered phase and the dynamical behaviour.

Chapter 3 Condensation

In this chapter we will describe the so-called condensed clustered phase. Before turning to the models of our interest we present the random subcubes model [ZDEB-8], where the condensation of clusters can be understood on a very elementary probabilistic level. After mentioning that the condensed phase is in fact very well known in spin glasses we describe the Poisson-Dirichlet process which determines the distribution of sizes of clusters in that phase. Further, we discuss general properties of the condensed phase in random CSPs. And finally we address our original question and conclude that the condensation is not much significant for the hardness of finding a solution [ZDEB-5].

The random-subcubes model [ZDEB-8] is defined by its solution space $S\subseteq\{0,1\}^{N}$; we define $S$ as the union of $\lfloor 2^{(1-\alpha)N}\rfloor$ random clusters (where $\lfloor x\rfloor$ denotes the integer value of $x$). A random cluster $A$ being defined as:

such that for each variable $i$, $\pi^{A}_{i}=\{0\}$ with probability $p/2$, $\{1\}$ with probability $p/2$, and $\{0,1\}$ with probability $1-p$. A cluster is here a random subcube of $\{0,1\}^{N}$. If $\pi^{A}_{i}=\{0\}$ or $\{1\}$, variable $i$ is said “frozen” in $A$; otherwise it is said “free” in $A$. In this model one given configuration $\sigma$ might belong to zero, one or several clusters.

We describe the static properties of the set of solutions $S$ in the random-subcubes model in the thermodynamic limit $N\to\infty$ (the two parameters $0\leq\alpha\leq 1$ and $0\leq p\leq 1$ being fixed and independent of $N$). The internal entropy $s$ of a cluster $A$ is defined as $\frac{1}{N}\log_{2}|A|$, i.e., the fraction of free variables in $A$. The probability ${\cal P}(s)$ that a cluster has internal entropy $s$ follows the binomial distribution

Then the number of clusters of entropy $s$, denoted ${\cal N}(s)$, is with high probability

where $D(x\parallel y)\equiv x\log_{2}{\frac{x}{y}}+(1-x)\log_{2}{\frac{1-x}{1-y}}$ is the binary Kullback-Leibler divergence.

We compute the total entropy $s_{\rm tot}=\frac{1}{N}\log_{2}|S|$. First note that a random configuration belongs on average to $2^{N(1-\alpha)}(1-\frac{p}{2})^{N}$ clusters. Therefore, if

then with high probability the total entropy is $s_{\rm tot}=1$.

Now assume $\alpha>\alpha_{d}$. The total entropy is given by a saddle-point estimation:

We denote by $s^{*}={\rm argmax}_{s}{[\Sigma(s)+s\,|\,\Sigma(s)\geq 0]}$ the fraction of free variables in the clusters that dominate the sum. Note that our estimation is valid (there is no double counting) since in every cluster the fraction of solutions belonging to more than one cluster is exponentially small as long as $\alpha>\alpha_{d}$.

For $\alpha>\alpha_{c}$, the maximum in (3.8) is realized by the largest possible cluster entropy $s_{\rm max}$, which is given by the largest root of $\Sigma(s)$. Then $s_{\rm tot}=s^{*}=s_{\rm max}$. We will show in the next section that in such a case almost all solutions belong to only a finite number of largest clusters. This phase is thus called condensed, in the sense that almost all solutions are ”condensed” in a small number of clusters.

In summary, for a fixed value of the parameter $p$, and for increasing values of $\alpha$, four different phases can be distinguished:

Liquid (replica symmetric) phase, $\alpha<\alpha_{d}$: almost all configurations are solutions.

Clustered (dynamical 1RSB) phase with many states, $\alpha_{d}<\alpha<\alpha_{c}$: an exponential number of clusters is needed to cover almost all the solutions.

Condensed clustered phase, $\alpha_{c}<\alpha<1$: a finite number of the biggest clusters covers almost all the solutions.

Unsatisfiable phase, $\alpha>1$: no cluster, hence no solution, exists.

2 New in CSPs, well known in spin glasses

The complexity function $\Sigma(s)$ (2.26) in random CSPs is counting the logarithm of the number of clusters per variable which have internal entropy $s$ per variable. We define dominating clusters in the same way as in the random subcubes model, that is clusters of entropy $s^{*}$ such that

In chap. 2 we discussed properties of the dynamical 1RSB phase, that is when $\Sigma(s^{*})>0$, in other words when there are exponentially many dominating clusters.

The condensed phase with $\Sigma(s^{*})=0$, described in the random subcubes model, exists also in random CSPs. And in the context of constraint satisfaction problems it was first computed and discussed in [MPR05] and [ZDEB-4]. However, historically it was the condensed phase where the 1RSB solution was first worked out [Par80c]. A very simple example of condensation can also be found in the random energy model [Der80, Der81]. As we discussed in the previous chapter 2, the dynamical 1RSB phase is well hidden within the replica solution — the total entropy is equal to the replica symmetric entropy, the overlap distribution is trivial and the two-point correlation functions decay to zero etc. All this changes in the condensed phase.

A small digression to the physics of glasses: In structural glasses, the analog of the condensation transition is well known for a long time, its discovery goes back to Kauzmann in 1948 who studied the configurational entropy of glassy materials. Configurational entropy is the difference between the total (experimentally measured) entropy and the entropy of a solid material, this thus corresponds to the complexity function. In the so called fragile structural glasses [Ang95] the extrapolated configurational entropy becomes zero at a positive temperature, nowadays called the Kauzmann temperature. The Kauzmann temperature in the real glasses is, however, only extrapolation. The equilibration time in glasses exceeds the observation time high above the Kauzmann temperature. It is a widely discussed question if there exists a true phase transition at the Kauzmann temperature or not, for a recent discussion see [DS01].

As we said, it is the condensed phase which was originally described by Parisi and his one-step replica symmetry breaking solution [Par80c]. Let us now briefly clarify the relation to the replica solution, similar reasoning first appeared in [Mon95]. In sec. 2.1 we called the Legendre transform of the complexity function the replicated free entropy $\Phi(m)$ (2.26). In the replica approach the replicated entropy $\Omega(m)=\Phi(m)/m$ is computed. From (2.28) follows