Statistical physics extracts macroscopic laws by averaging over the many microscopic degrees of freedom of a system. Disordered systems demand a second and far harder average, one over the quenched randomness itself. The classic analytical routes, the replica and cavity methods, become uncontrolled outside mean-field or tree-like limits, and conventional numerical algorithms like parallel tempering require expensive, independent equilibration for every disorder realization. In this work, we introduce a universal neural variational framework that amortizes inference across the disorder ensemble, eliminating both the need for per-instance Markov chain equilibration and the cost of retraining instance-specific variational ansatzes. Built on an encoder-decoder Transformer architecture, after training once, it produces an explicit approximation to the Boltzmann distribution given previously unseen disorder realizations without further optimization. We validate this framework on 2D Edwards-Anderson models, and apply it to the random-bond Ising model, successfully capturing the Binder cumulant crossings near the Nishimori multicritical point. These results shift the object of variational inference from the single instance to the disorder ensemble, opening a route to frustrated many-body systems where instance-by-instance computation is prohibitive.
Statistical physics extracts macroscopic laws by averaging over the many microscopic degrees of freedom of a system. Disordered systems demand a second and far harder average, one over the quenched randomness itself. The classic analytical routes, the replica and cavity methods, become uncontrolled outside mean-field or tree-like limits, and conventional numerical algorithms like parallel tempering require expensive, independent equilibration for every disorder realization. In this work, we introduce a universal neural variational framework that amortizes inference across the disorder ensemble, eliminating both the need for per-instance Markov chain equilibration and the cost of retraining instance-specific variational ansatzes. Built on an encoder-decoder Transformer architecture, after training once, it produces an explicit approximation to the Boltzmann distribution given previously unseen disorder realizations without further optimization. We validate this framework on 2D Edwards-Anderson models, and apply it to the random-bond Ising model, successfully capturing the Binder cumulant crossings near the Nishimori multicritical point. These results shift the object of variational inference from the single instance to the disorder ensemble, opening a route to frustrated many-body systems where instance-by-instance computation is prohibitive.
When independent learners are trained on data from the same disordered system, finite and noisy samples can lead them to different descriptions of the underlying physics. We quantify the reproducibility of such descriptions through consensus, the degree of agreement among independently trained models. In a teacher-student framework on the two-dimensional random-bond Ising model, an ensemble of students each infers a single uniform effective coupling from a finite set of equilibrium configurations drawn from one quenched bond realization. Aggregating the inferred couplings across students and realizations yields the consensus distribution, whose width quantifies the agreement among learners. The variance of this distribution separates exactly into a finite-sampling contribution, bounded below by inverse of Fisher information, and a quenched-disorder contribution. Sum of these two variance is minimized near critical region, so the consensus width has a minimum near the pseudocritical temperature that sharpens with increasing disorder. The two minima arise from different mechanisms: the sampling term inherits the critical peak of the Fisher information, while the disorder term is temperature-independent at leading order and acquires its critical minimum only at fourth order in the disorder strength, through the covariance between the linear and cubic responses to the bond disorder. Because the students fit a uniform coupling to a heterogeneous lattice, the inferred coupling also carries a misspecification bias that survives in the infinite-data limit. This bias displaces the minimum of the mean-squared error above the consensus minimum, so the temperature at which independent learners agree most closely is not the temperature at which their shared answer is most faithful.
Variational neural network models have achieved remarkable success in solving ground-state problems of quantum many-body systems. However, addressing classical and quantum spin glasses remains challenging, as disorder and energy frustration give rise to an exponentially large number of local energy minima separated by high-energy barriers, hindering the efficiency of conventional Metropolis-based Monte Carlo methods. To bridge this gap, we introduce Deep Boltzmann Quantum States, a class of neural quantum states inspired by deep Boltzmann machines that inherit efficient block Gibbs sampling. We also propose two key advances in the training algorithm. Firstly, we combine natural-gradient updates with state-of-the-art stochastic optimizers. Secondly, we gradually tune the hardness of the problem Hamiltonian by interpolating from an easy to a hard regime, without the need to closely approximate the instantaneous adiabatic state at intermediate times. We match the exact solution or the best available estimate for several instances of classical and quantum Ising spin-glass models with infinite-range interactions and hundreds of spins. We also solve instances of the NP-hard Job Shop Scheduling Problem exceeding the current limitations of quantum annealing hardware. To summarize, deep neural architectures with efficient global update rules and trained within an annealing-like scheme, provide a powerful framework for solving real-world hard combinatorial optimization and for investigating disordered quantum many-body systems.
Efficient sampling of two-dimensional statistical physics systems remains a central challenge in computational statistical physics. Traditional Markov chain Monte Carlo (MCMC) methods, including cluster algorithms, provide only partial solutions, as their efficiency collapses for large systems in the presence of frustration and quenched disorder. The recently proposed Tensor Network Monte Carlo (TNMC) method offers a promising alternative, yet its original implementation suffers from inefficiencies due to the lack of scalable parallel sampling. In this work, we introduce BatchTNMC, a GPU-optimized and parallelized implementation of TNMC tailored for large-scale simulations of two-dimensional spin glasses. By leveraging batch processing and parallel sampling across multiple disorder realizations, our implementation achieves speedups of up to five orders of magnitude compared with the original serial scheme. Benchmarking on two-dimensional spin glasses demonstrates dramatic gains in efficiency: for instance, on a single GPU, BatchTNMC concurrently produces 1000 uncorrelated and unbiased samples across 1000 disorder realizations on 1024×1024 lattices in just 3.3 hours, with an acceptance probability of 37%. These results establish BatchTNMC as a scalable and powerful computational framework for the study of two-dimensional disordered spin glass systems.
We introduce a Metropolis-Hastings Markov chain for Boltzmann distributions of classical spin systems. It relies on approximate tensor network contractions to propose correlated collective updates at each step of the evolution. We present benchmark computations for a wide variety of instances of the two-dimensional Ising model, including ferromagnetic, antiferromagnetic, (fully) frustrated and Edwards-Anderson spin glass instances, and we show that, with modest computational effort, our Markov chain achieves sizeable acceptance rates, even in the vicinity of critical points. In each of the situations we have considered, the Markov chain compares well with other Monte Carlo schemes such as the Metropolis or Wolff’s algorithm: equilibration times appear to be reduced by a factor that varies between 40 40 and 2000 2000 , depending on the model and the observable being monitored. We also present an extension to three spatial dimensions, and demonstrate that it exhibits fast equilibration for finite ferro- and antiferromagnetic instances. Additionally, and although it is originally designed for a square lattice of finite degrees of freedom with open boundary conditions, the proposed scheme can be used as such, or with slight modifications, to study triangular lattices, systems with continuous degrees of freedom, matrix models, a confined gas of hard spheres, or to deal with arbitrary boundary conditions.
Thermal equilibrium states of many-body Hamiltonians are essential for probing quantum chaos, finite-temperature phases of matter, and training quantum machine learning models, yet generating large collections of such states across different Hamiltonians remains costly with existing methods. We introduce a powerful operation, the quantum thermal-drift channel, to construct a measurement-controlled sampling algorithm that autonomously generates thermal states together with their system Hamiltonians as labels for general physical models. We prove that our algorithm is efficient: the total gate count scales polynomially with system size and quadratically with inverse temperature, providing the first polynomial resource bound for random thermal state generation. We characterize the distribution of sampled Hamiltonians as a normal distribution reweighted by partition functions, which quantifies a trade-off between sampling accuracy and effective label range. Level-spacing statistics computed from sampled thermal states of a 2D transverse-field Ising model show a crossover to Wigner-Dyson universality, confirming that the sampler captures nontrivial chaotic correlations. Finally, a variational quantum classifier trained on the generated dataset achieves near-optimal accuracy in predicting Hamiltonian properties of unseen states. These results establish a scalable, quantum-native route for thermodynamic simulation and labeled quantum data generation in many-body systems.
We present a bootstrap framework that yields rigorous two-sided bounds on ground-state observables of quantum systems with quenched disorder. The key idea is to extend the operator algebra to include classical disorder variables and to encode the disorder distribution through a set of moments. The method directly bounds disorder-averaged expectation values for finite or infinite lattice systems and for arbitrary disorder distributions. In contrast to conventional approaches, we do not need to sample disorder realizations, which enables us to make thorough use of symmetries that hold only on average rather than for each disorder realization. Additionally, we show that multi-replica operators can be used to access disorder fluctuations such as the variance of an observable. We apply the method to the one-dimensional random transverse-field Ising model (RTFIM) with discrete, uniform, and Gaussian disorder, obtaining tight bounds on the disorder-averaged ground-state energy density and short-range spin-spin correlators, as well as bounds on the disorder variance of the energy density. Additionally, we show that the difference between the upper and lower bounds on the disorder-averaged ground-state energy density, which we term the bootstrap gap, can be used to map out the phase diagram of the RTFIM. The bootstrap gap is small and flat within the paramagnetic and ferromagnetic phases at weak disorder, grows upon entering the Griffiths regions, and peaks along the critical line.