Variational Neural Annealing

• Variational Neural Annealing is a framework that trains autoregressive models using a temperature schedule to progressively concentrate probability mass on low-energy configurations.
• It replaces slow Markov-chain sampling with exact autoregressive sampling while minimizing a free-energy objective equivalent to a forward Kullback–Leibler divergence to the Boltzmann target.
• The method has been effectively applied to spin systems, combinatorial optimization, quantum chemistry, and portfolio optimization by mitigating mode collapse and improving convergence.

Searching arXiv for the papers on arXiv and recent VNA-related work to ground the article. arxiv_search(query="Variational Neural Annealing autoregressive annealing variational free energy", max_results=10) arxiv_search(query="Message Passing Variational Autoregressive Network annealing Ising (Ma et al., 2024)", max_results=5) Variational Neural Annealing (VNA) is a variational annealing framework in which a tractable parametric distribution is trained to follow an annealed target distribution as temperature is lowered, so that probability mass progressively concentrates on low-energy or otherwise high-quality configurations. In its canonical form, introduced as a variational alternative to simulated annealing for classical and quantum spin systems, VNA replaces slow Markov-chain sampling by exact sampling from autoregressive models and optimizes a free-energy objective equivalent to a forward Kullback–Leibler divergence to the Boltzmann target (Hibat-Allah et al., 2021). Subsequent work has extended the same organizing idea to combinatorial optimization, dense spin glasses, lattice protein folding, autoregressive neural quantum states, multimodal black-box optimization, flow-based variational inference, and large-scale portfolio optimization (Khandoker et al., 2022).

1. Origin, scope, and terminology

The name VNA is rooted in the 2021 formulation of variational classical annealing and variational quantum annealing, where annealing is recast as variational optimization over an autoregressive model rather than direct sampling from the exact Boltzmann distribution (Hibat-Allah et al., 2021). In that formulation, the learned distribution $q_\theta(x)$ is annealed by minimizing a variational free energy while the temperature decreases from a high-temperature regime, which favors broad support, to a low-temperature regime, which favors concentration near ground states.

In later literature, closely related terminology is used with domain-specific variants. The paper on recurrent neural networks for combinatorial optimization uses the name Variational Classical Annealing (VCA) but explicitly identifies it with “Variational Neural Annealing” as a variational neural approach that trains an autoregressive neural distribution to track annealed Boltzmann targets (Khandoker et al., 2022). The message-passing Ising work presents VNA as a training paradigm for a parametric sampler $p_\theta(s)$ that approximates the Boltzmann distribution $p_\beta(s)$ of an Ising Hamiltonian along an annealing schedule (Ma et al., 2024). In neural quantum states for quantum chemistry, VNA is used in a different but related sense: entropy regularization is added to the variational Monte Carlo objective and its weight is annealed to zero during training, thereby promoting exploration early and energy minimization late (Knitter et al., 2024).

A broader generalization appears in Natural Variational Annealing, where the search distribution is a variational posterior such as a mixture of Gaussians, the objective is entropy regularized, and optimization is carried out by natural-gradient learning (Minh et al., 8 Jan 2025). Flow-based variational inference with tempered posteriors, as implemented in LINFA, does not use the acronym VNA explicitly, but it realizes the same annealing pattern by optimizing a neural variational family against a sequence of tempered targets $p^t(z,x)$ with $t \uparrow 1$ (Wang et al., 2023). Across these uses, the common structure is a variational distribution, an annealing or tempering schedule, and a training rule that converts exploration at high temperature into concentration at low temperature.

2. Variational objectives and annealing mechanisms

The classical formulation of VNA starts from a Hamiltonian or cost function and the associated Boltzmann distribution. For Ising systems, one writes

$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$

A neural autoregressive model

$$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$

is then trained to minimize the forward divergence $D_{KL}(p_\theta\Vert p_\beta)$, or equivalently the variational free energy. In the notation used for fully connected Ising models,

$$J(\theta;\beta)=\mathbb{E}_{p_\theta}[\beta H(\mathbf{s})+\log p_\theta(\mathbf{s})],$$

while in temperature notation the same trade-off is written as

$$F_\theta(T)=\mathbb{E}_{p_\theta}[E(x)]-T\,S(p_\theta) =\mathbb{E}_{p_\theta}[E(x)+T\log p_\theta(x)].$$

Both forms express the same balance between energy minimization and entropy maximization (Ma et al., 2024).

The score-function identity provides the basic stochastic gradient estimator. In the classical VNA/VCA literature the gradient is written as an expectation of $p_\theta(s)$0 weighted by the local free-energy term, and samples are obtained exactly by ancestral sampling from the autoregressive factorization. This removes the dependence on Metropolis acceptance, autocorrelated chains, and equilibration of an external sampler (Hibat-Allah et al., 2021). Warm starts across the annealing schedule are a defining operational feature: parameters at the next inverse temperature are initialized from the previous one, so the learned distribution tracks a sequence of increasingly sharp targets rather than being retrained from scratch.

A related but distinct formulation is used in autoregressive neural quantum states for chemistry. There the standard variational Monte Carlo energy

$p_\theta(s)$1

is replaced by an entropy-regularized objective

$p_\theta(s)$2

with an annealing schedule

$p_\theta(s)$3

Here $p_\theta(s)$4 is the weight of the entropy regularizer rather than the inverse physical temperature. The underlying curriculum is nevertheless analogous: early training encourages broad support over the physically allowed sector, while later training removes the entropy term and recovers the pure energy objective (Knitter et al., 2024).

In multimodal optimization, Natural Variational Annealing uses

$p_\theta(s)$5

whose maximizer is the Gibbs measure

$p_\theta(s)$6

As $p_\theta(s)$7, the variational search distribution concentrates on optima, and annealing controls the exploration–exploitation transition (Minh et al., 8 Jan 2025). The same mathematical pattern also underlies tempered posterior annealing in normalizing-flow variational inference, where the target is replaced by a power posterior $p_\theta(s)$8 and the temperature parameter $p_\theta(s)$9 is increased until the true target is recovered (Wang et al., 2023).

3. Model classes and algorithmic structure

The original VNA work used autoregressive recurrent architectures because they are normalized by construction, can be sampled exactly in $p_\beta(s)$0 time, and provide exact log-probabilities needed for low-variance score-function gradients (Hibat-Allah et al., 2021). The paper develops vanilla, tensorized, two-dimensional, and dilated RNN parameterizations, with the dilated variant reducing the effective path length between distant variables to $p_\beta(s)$1 and improving performance on fully connected systems such as the Sherrington–Kirkpatrick and Wishart Planted Ensemble models (Hibat-Allah et al., 2021).

Later work greatly broadened the model family. VCA for Max-Cut, nurse scheduling, and the traveling salesman problem used vanilla and dilated RNNs with site-dependent or shared parameters, together with masking in the TSP decoder to enforce permutation constraints (Khandoker et al., 2022). The MPVAN architecture augmented a variational autoregressive network with Hamiltonian-aware message passing that explicitly uses coupling strengths and signs,

$p_\beta(s)$2

thereby encoding local-field information directly in the conditional generator rather than asking the autoregressive network to infer it implicitly (Ma et al., 2024).

Autoregressive neural quantum states introduced MADE, transformer, and RetNet ansätze whose Born distribution is exactly samplable. In that setting, exact autoregressive sampling is combined with hard masks enforcing the electron-number sector and with reverse Jordan–Wigner ordering as a heuristic for reducing the correlation burden on the conditional chain (Knitter et al., 2024). A notable computational result is the RetNet-versus-transformer inference threshold

$p_\beta(s)$3

under $p_\beta(s)$4, beyond which recurrent RetNet inference uses fewer FLOPs than transformer inference (Knitter et al., 2024).

Constraint handling has also become a central design axis. In lattice protein folding, a dilated RNN over move sequences is combined with strict masking of invalid actions by setting invalid logits to $p_\beta(s)$5 before log-softmax. Sampling uses the masked distribution, while gradients are evaluated with the unmasked probabilities in an upper-bound training scheme,

$p_\beta(s)$6

which was introduced to avoid constraining gradients by the masking projection (Khandoker et al., 28 Feb 2025). In normalizing-flow variational inference, the neural family is not autoregressive in the same sense, but annealing is inserted by tempering the target density rather than the decoder, so the ELBO becomes an annealed free-energy objective over MAF or RealNVP transformations (Wang et al., 2023).

4. Empirical domains and representative results

VNA has been applied to a wide range of optimization and inference problems. The empirical record is not uniform across all formulations, but the recurring observation is that annealed variational training can improve low-temperature sampling, residual energies, or multimodal coverage relative to non-annealed training and several classical baselines.

Domain	Formulation	Representative result
Spin glasses	VCA and MPVAN for Ising Hamiltonians	VCA outperformed SA and SQA asymptotically on random chains, 2D EA, SK, and WPE; MPVAN achieved much smaller KL and better residual energies on dense frustrated Ising models (Hibat-Allah et al., 2021)
Combinatorial optimization	RNN-based VCA	On Max-Cut, NSP, and TSP, VCA outperformed SA on average in the asymptotic limit by one or more orders of magnitude in terms of relative error, reaching up to $p_\beta(s)$7 cities for TSP (Khandoker et al., 2022)
Quantum chemistry	Autoregressive NQS with entropy-annealed training	VNA-augmented MADE, transformer, and RetNet matched or nearly matched CCSD and FCI across seven small molecules (Knitter et al., 2024)
Lattice protein folding	Masked dilated RNN with variational annealing	The method matched the best-known energies on Istrail benchmark sequences up to $p_\beta(s)$8 beads (Khandoker et al., 28 Feb 2025)
Portfolio optimization	Classical VNA with autoregressive RNN	VNA handled more than $p_\beta(s)$9 assets and achieved performance comparable to Mosek while converging faster on hard instances (Ranabhat et al., 9 Jul 2025)

The 2021 spin-glass study established the baseline empirical case for VNA. On $p^t(z,x)$0D random chains, $p^t(z,x)$1D Edwards–Anderson systems up to $p^t(z,x)$2, the Sherrington–Kirkpatrick model with $p^t(z,x)$3, and the Wishart Planted Ensemble with $p^t(z,x)$4, variational classical annealing produced residual energies that decreased faster with the number of annealing steps than simulated annealing or simulated quantum annealing, and the long-time gaps could reach orders of magnitude on glassy instances (Hibat-Allah et al., 2021). The same work also showed that “CQO (VCA with $p^t(z,x)$5, i.e., no annealing) plateaus,” which directly tied the gains to annealing rather than to the network class alone (Hibat-Allah et al., 2021).

The message-passing extension sharpened that picture for dense, frustrated Ising models. MPVAN was evaluated on the Wishart Planted Ensemble for $p^t(z,x)$6 and on the Sherrington–Kirkpatrick model for $p^t(z,x)$7. It achieved smaller $p^t(z,x)$8 on enumerable WPE instances, lower free-energy relative error, and the best residual energy per site across tested baselines, with advantages that increased with graph density and system size (Ma et al., 2024). On WPE, MPVAN often found ground states up to $p^t(z,x)$9 and attained the lowest residual energy at $t \uparrow 1$0, where none of the methods reached the true ground state (Ma et al., 2024).

Outside spin glasses, VCA was shown to generalize to structured discrete optimization. For unweighted Max-Cut, nurse scheduling with quadratic penalties, and Euclidean TSP, the recurrent VCA framework used the same entropy-regularized free-energy objective and linear temperature schedule. The paper reports that average residual error at slow annealing is one or more orders of magnitude lower than SA, while best-case performance on TSP reached the exact ground state for $t \uparrow 1$1 at $t \uparrow 1$2 and near-exact solutions for $t \uparrow 1$3 and $t \uparrow 1$4 (Khandoker et al., 2022).

The protein-folding instantiation pushed VNA into a constrained generative setting. Using masked sampling and upper-bound training on two-dimensional HP lattice folding, the method matched the best-known energies for 20merA, 20merB, 24mer, 25mer, 36mer, 48mer, 50mer, and 60mer sequences. For the $t \uparrow 1$5-mer, the reported energy was $t \uparrow 1$6, and prior works did not report results for this size (Khandoker et al., 28 Feb 2025).

In autoregressive neural quantum states for ab initio chemistry, VNA functioned as a training strategy rather than as a thermal sampler. On H$t \uparrow 1$7O, N$t \uparrow 1$8, O$t \uparrow 1$9, H$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$0S, PH$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$1, LiCl, and Li$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$2O, quartic-schedule entropy annealing upgraded MADE, transformer, and RetNet models so that all three matched or nearly matched CCSD and FCI values across the benchmark slate. The ablation without VNA showed large degradations for several RetNet and transformer cases, particularly for N$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$3, O$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$4, LiCl, and Li$$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$5O (Knitter et al., 2024).

Large-scale portfolio optimization provides a different empirical regime: dense quadratic structure, constraints, and industrial instance sizes. After mapping the optimization problem to an Ising-like Hamiltonian, classical VNA identified near-optimal solutions on the S&P 500, Russell 1000, and Russell 3000, including portfolios with $$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$6 assets. The reported hard-instance wall-clock times were $$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$7 s versus $$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$8 s on S&P 500, $$p_\beta(\mathbf{s})=\frac{1}{Z(\beta)}\exp(-\beta H(\mathbf{s})), \qquad Z(\beta)=\sum_{\mathbf{s}}\exp(-\beta H(\mathbf{s})).$$9 s versus $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$0 s on Russell 1000, and $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$1 s versus $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$2 s on Russell 3000 for VNA and Mosek respectively under the paper’s comparison settings (Ranabhat et al., 9 Jul 2025).

5. Exploration, mode collapse, and theoretical interpretation

A central theme in the VNA literature is that annealing is not merely a scheduling convenience; it is a mechanism for preserving support, delaying premature concentration, and improving low-temperature optimization. The fully connected Ising work makes this explicit: warm-starting across inverse temperatures “stabilizes training, prevents catastrophic mode collapse, and moves probability mass toward lower energy configurations,” while the Hamiltonian-aware message passing further mitigates collapse by encoding couplings directly in the network (Ma et al., 2024). In that study, negative entropy was used to monitor collapse, and MPVAN delayed the onset of mode collapse relative to VAN and VCA, with collapse occurring after at least $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$3 training steps in MPVAN versus approximately $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$4 in the baselines (Ma et al., 2024).

The chemistry formulation states the same point in a different language. There VNA “regularizes the NQS objective with an entropy term and anneals its weight to zero during training,” enlarging the effective support of the Born distribution during the early and middle stages. The paper explicitly links this to entropy-regularized reinforcement learning and to the practical need to prevent the sampler from prematurely collapsing onto a small subset of Fock configurations (Knitter et al., 2024). The same paper notes that too rapid decay of the annealing schedule can harm accuracy, whereas a quartic schedule with warm-start is effective across MADE, transformer, and RetNet (Knitter et al., 2024).

The 2026 theoretical study of annealed variational inference on Gaussian mixtures provides the sharpest formal analysis of this phenomenon. It studies annealed reverse-KL minimization to the tempered density $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$5 and derives a closed-form expression for the probability of mode collapse under an exponential schedule (Fogliani et al., 13 Feb 2026). In the high-temperature regime, entropy dominates and drives anti-alignment of the learned means; in the low-temperature regime, the means decouple and are pulled toward the target basins. The derived collapse probability

$$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$6

with

$$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$7

makes explicit the schedule dependence of collapse prevention (Fogliani et al., 13 Feb 2026). The theoretical conclusion is that successful annealing requires enough time in the high-temperature regime to separate modes, but not so high an initial temperature that separation becomes too slow.

This perspective also clarifies a recurring misconception: annealing is not equivalent to training the same variational family directly at zero temperature. In the original VNA paper, classical optimization without annealing plateaus on the two-dimensional Edwards–Anderson model, whereas the annealed version continues to improve (Hibat-Allah et al., 2021). In the protein-folding work, starting at nonzero $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$8 consistently produces lower average energies than $$p_\theta(\mathbf{s})=\prod_{i=1}^N p_\theta(s_i\mid s_{<i})$$9, again isolating entropy regularization as the operative ingredient rather than the architecture alone (Khandoker et al., 28 Feb 2025).

6. Relation to adjacent frameworks, practical limits, and open directions

VNA sits at the intersection of simulated annealing, variational inference, autoregressive generative modeling, and, in some formulations, variational Monte Carlo. It differs from simulated annealing because it does not rely on local Markov moves on the target landscape; instead it learns a reusable sampler and then draws independent samples from it (Khandoker et al., 2022). It differs from parallel tempering and simulated quantum annealing in the same respect: exact ancestral sampling replaces the need to equilibrate chains or Trotterized replicas (Hibat-Allah et al., 2021). It differs from standard deterministic annealing variational inference because much of the VNA literature emphasizes exact sampling from neural autoregressive models, explicit free-energy minimization, and extraction of low-energy configurations after the anneal.

At the same time, several nearby methods are structurally compatible with VNA. LINFA implements power-posterior annealing with normalizing flows by replacing $D_{KL}(p_\theta\Vert p_\beta)$0 with $D_{KL}(p_\theta\Vert p_\beta)$1 and increasing $D_{KL}(p_\theta\Vert p_\beta)$2 either linearly or adaptively through AdaAnn,

$D_{KL}(p_\theta\Vert p_\beta)$3

which is a deterministic annealing VI mechanism over neural variational families (Wang et al., 2023). The finite-temperature quantum many-body work based on a Rényi-2 free-energy functional does not require an annealing schedule, but it explicitly notes that warm-starting across $D_{KL}(p_\theta\Vert p_\beta)$4 values is straightforward and can improve convergence, making it naturally compatible with a VNA-style continuation strategy (Lu et al., 2024). Natural Variational Annealing extends the principle to multimodal black-box optimization with mixtures of Gaussians and natural gradients, showing that annealed entropy regularization need not be tied to Boltzmann energies over discrete states (Minh et al., 8 Jan 2025).

The main limitations reported across the literature are schedule sensitivity, model-capacity dependence, and constraint handling. In dense Ising systems, extremely large $D_{KL}(p_\theta\Vert p_\beta)$5 can still induce residual mode collapse, and the reported remedies are gradual annealing, adequate batch sizes, and Hamiltonian-aware message passing (Ma et al., 2024). In neural quantum states, overly rapid decay of the entropy regularizer harms accuracy, while aggressive caps on unique samples or discarding singleton samples may not generalize beyond the tested system sizes (Knitter et al., 2024). In lattice protein folding, removing masking or replacing the upper-bound objective with masked-log training degrades convergence, and when all moves become invalid a forced invalid action is taken and assigned $D_{KL}(p_\theta\Vert p_\beta)$6 (Khandoker et al., 28 Feb 2025). In portfolio optimization, feasibility depends strongly on penalty weights, and the paper identifies augmented Lagrangian methods as a promising direction for robust enforcement of multiple constraints (Ranabhat et al., 9 Jul 2025).

A plausible implication is that VNA is better understood as a design pattern than as a single fixed algorithm. The pattern consists of a variational family, an entropy-aware or temperature-aware objective, a continuation schedule, and an exact or otherwise tractable sampling mechanism. Within that pattern, the literature now spans tensorized and dilated RNNs, message-passing autoregressive networks, MADE, transformers, RetNet, normalizing flows, and mixture models, with applications ranging from disordered Ising Hamiltonians to quantum chemistry and financial optimization. The common claim supported across these papers is narrower but technically precise: annealed variational training can materially improve exploration, reduce mode collapse, and strengthen low-temperature or low-entropy optimization when the underlying model family is sufficiently expressive and the continuation schedule is chosen carefully (Hibat-Allah et al., 2021).