Combinatorial optimization problems near algorithmic phase transitions represent a fundamental challenge for both classical algorithms and machine learning approaches. Among them, graph coloring is a prototypical constraint satisfaction problem characterized by phase transitions that fragment the solution space into distinct clusters, making the search for valid solutions challenging. Here we introduce a physics-inspired neural framework that learns to solve large-scale graph coloring instances by combining graph neural networks with statistical-mechanics principles. Our approach integrates a planting-based supervised signal, symmetry-breaking regularization, and iterative noise-annealed neural dynamics to navigate complex solution landscapes. When the number of iterations scales quadratically with graph size, the learned solver successfully tackles increasingly challenging instances, even in regimes where finding valid solutions becomes particularly difficult, demonstrating its ability to extract meaningful information from the graph structure. The model generalizes from small training graphs to instances up to 20 times larger, demonstrating that neural architectures can learn scalable algorithmic strategies that remain effective in hard connectivity regions. These results establish a general framework for learning neural solvers that operate near fundamental phase boundaries in combinatorial optimization and inference.
Figures & tables
Figure 1: Schematic illustration of the evolution of the solution-space landscape with increasing average connectivity c . The red cross represent the planted solution. At low connectivity, the planted solution belongs to a large connected region of solutions. Above the dynamical transition cd , the solution space fragments into an exponential number of clusters. At the condensation transition cc , the solution measure concentrates onto a subexponential number of dominant clusters, one of which contains the planted solution. The satisfiability threshold cs marks the transition beyond which typical random instances become unsatisfiable. The Kesten–Stigum threshold cKS marks the detectability of the planted structure. Blue regions denote solutions, whereas orange regions represent local minima.
Figure 2: (DTR Model) Extensive energy as a function of the number of iterations for four representative connectivity values, c=12.0,13.0,13.5,14.0 . Each data point is obtained by averaging results over 50 graphs with N= 10k nodes. Error bars indicate one standard deviation across graph instances. Circles indicate values for planted graphs, while crosses indicate values for random graphs. Error bars indicate standard error of the mean across graph instances.
Figure 3: (DTR Model) Main panel Mean intensive energy as a function of connectivity, for both random and planted graphs, and a number of nodes varying from N=2500 to N=50 k. Each data point is obtained by running the algorithm for Niter=200 k iterations and averaging the energy over 50 graphs. Inset The mean extensive energy grows with the size at each value of c and does not allow to determine the algorithmic threshold in the large N limit. Error bars indicate the standard error of the mean.
Figure 4: (DTR Model) Probability of finding a solution as a function of the graph connectivity c , for sizes N∈[1024,2048,4096,8192] and a number of iterations scaling as Niter=0.04N2 . Each point is the average over 400 graphs. Error bars indicate 95% Wilson binomial confidence intervals.
Figure 5: (DTR Model) Comparison of the probability of finding a solution as a function of the graph connectivity c for different solving algorithms. Values are obtained by averaging over 400 graphs of N=1024 nodes. The number of iterations Niter is chosen to ensure the same wall-clock execution time across all algorithms. The tested algorithms are a recurrent Physics-Inspired GNN (rPI) [ 35 ] , Belief Propagation (BP) [ 37 ] , Simulated Annealing (SA) [ 2 ] , our graph neural network (GNN), and Focused Metropolis Search (FMS) [ 36 ] . Error bars indicate 95% Wilson binomial confidence intervals.
Figure 6: (KSR Model) Main panel Probability of finding a solution as a function of the graph connectivity c , for sizes 1000≤N≤20000 and a number of iterations scaling as Niter=10−4N2 . Each point is the average over 2106/N graphs. Inset The mean overlap with the planted coloring of the configuration output by the algorithm. Error bars indicate 95% Wilson binomial confidence intervals.
Combinatorial optimization finds the best feasible choices among discrete alternatives in routing, scheduling, graph problems, and resource allocation.
Dept. of Computer Science and Operations Research, Université de Montréal · Mila – Quebec AI Institute · Dept. of Mathematics & Statistics, Université de Montréal +2