RL-GQE: Reinforcement Learning towards Generative Quantum Eigensolving

Application of Group Relative Policy Optimization and Self-Imitation Learning enable an autoregressive transformer to construct quantum circuits from scratch, reaching 96% of the Heisenberg ground state energy without a critic network.

Cover image for RL-GQE: Reinforcement Learning towards Generative Quantum Eigensolving

Motivation

Bypassing the exponential data-generation bottleneck inherent to supervised learning requires a generative model capable of exploring the Hilbert space autonomously through pure Reinforcement Learning. However, current policy gradient methods fail in frustrated spin landscapes due to extreme gradient instability. The core hypothesis of this research is that this instability is structurally tied to the reliance on parameterized value networks (critics) and an algorithmic inability to efficiently exploit sparse reward signals within highly chaotic energy surfaces.

To resolve this, this study introduces a novel architectural framework: a pure RL Generative Quantum Eigensolver driven by Group Relative Policy Optimization (GRPO) (Shao et al., 2024) and Self-Imitation Learning (SIL) (Oh et al., 2018). By applying GRPO, this approach completely bypasses the need for a global value network. Instead, it estimates advantages by scoring discrete circuit trajectories relative to one another within the same generation batch, relying directly on immediate physical energy evaluations to normalize the extreme variance of the quantum landscape. Ultimately, validating this gradient stability and convergence capability on an unconditioned model serves as the mandatory first step toward realizing a zero-shot, prompt-conditioned eigensolver capable of generalizing across diverse physical systems.

Methodology

1. Quantum environment and Hamiltonian formulation

The target physical system is a 1D Isotropic Heisenberg XXX model with open boundary conditions, consisting of N=4N=4 qubits. The system is subjected to an external magnetic field along the Z-axis.

The Hamiltonian is defined mathematically as:

H=Ji=1N1(σixσi+1x+σiyσi+1y+σizσi+1z)+hi=1NσizH = J \sum_{i=1}^{N-1} (\sigma_i^x \sigma_{i+1}^x + \sigma_i^y \sigma_{i+1}^y + \sigma_i^z \sigma_{i+1}^z) + h \sum_{i=1}^{N} \sigma_i^z

Based on the physical configurations tested, the spin-spin coupling constant is set to J=1J=1, and the external magnetic field strength is set to h=2h=2. Ground state energy expectation evaluations ψHψ\langle \psi | H | \psi \rangle serve as the absolute metric for policy performance, with the theoretical minimum established at E0=7.8284E_0 = -7.8284.

This specific configuration was selected as a targeted proof-of-concept benchmark. Unlike diagonal Ising models, the Heisenberg XXX model features non-commuting interactions (XXX \otimes X and YYY \otimes Y) that require the generative agent to map quantum entanglement. It is notable that the relatively strong magnetic field (h=2h=2) partially polarizes the ground state, ensuring the landscape remains tractable for a blind, unconditioned RL agent. However, the system remains mathematically complex enough to trigger the gradient collapse observed in standard policy gradient baselines. Thus, it provides an ideal environment to isolate and validate the stabilizing mechanics of the proposed architecture prior to scaling.

2. Discrete action space and ansatz construction

Unlike traditional Variational Quantum Eigensolvers (VQEs) that rely on parameterized continuous optimization, this study evaluates the ability of the Reinforcement Learning (RL) agent to strictly define circuit topology using fixed-angle operations.

The discrete action space consists of Pauli rotation gates. The allowed rotation angles are defined as:

θ{±π,±π2,±π4,±π8,±π16,±π32}\theta \in \left\{ \pm\pi, \pm\frac{\pi}{2}, \pm\frac{\pi}{4}, \pm\frac{\pi}{8}, \pm\frac{\pi}{16}, \pm\frac{\pi}{32} \right\}

The agent can select from two distinct classes of operations:

  • Two-qubit interactions: RXX(θ)R_{XX}(\theta), RYY(θ)R_{YY}(\theta), RZZ(θ)R_{ZZ}(\theta) applied to adjacent qubits (ii, i+1i+1).
  • Single-qubit rotations: RX(θ)R_X(\theta), RY(θ)R_Y(\theta), RZ(θ)R_Z(\theta) applied to any individual qubit ii.

To artificially prune nonsensical operations and constrain the search space, an Action Mask is dynamically applied during generation. This mask strictly forbids the placement of consecutive identical Pauli operations on the same subset of wires (e.g. preventing RX(θ/2)R_X(\theta/2) followed immediately by another RXR_X on qubit 0).

3. Policy architecture

The quantum circuit generator is parameterized by πθ\pi_\theta, modeled as an autoregressive transformer (GPT). The architecture features 4 transformer blocks, an embedding dimension of 384, and 8 attention heads. Generation is constrained to a maximum sequence length of 10 gates. Because the model acts unconditionally, it implicitly learns the physical energy landscape of the target Hamiltonian through trial and error.

4. Reward shaping and discounted returns

The quantum circuit assembly is modeled as a Markov Decision Process (MDP). At each time step tt, the addition of a quantum gate yields an immediate reward based on the step-wise change in energy, penalized by a small gate cost (c=0.005c=0.005) to encourage shallower circuit depths:

rt=(Et1Et)cr_t = (E_{t-1} - E_t) - c

To ensure that early gate placements are optimized relative to the final sequence state they construct, the pipeline utilizes a discounted cumulative return with a discount factor γ=0.8\gamma=0.8:

Rt=rt+γRt+1R_t = r_t + \gamma R_{t+1}

5. Policy optimization mechanisms

The training pipeline achieves a final evaluation energy of 7.5173-7.5173 (approx. 96% accuracy) using a dual-loss optimization strategy: Group Relative Policy Optimization (GRPO) anchored by Self-Imitation Learning (SIL).

Group Relative Policy Optimization (GRPO) eliminates the need for a critic or value network by estimating advantages purely through group normalization. For a generated batch of sequences, the advantages AtA_t are standard-scored against the batch mean μR\mu_R and standard deviation σR\sigma_R. The surrogate objective is clamped via a hyperparameter ϵ\epsilon, and a Kullback-Leibler (KL) divergence penalty ensures stable policy updates relative to a frozen reference model πref\pi_{\mathrm{ref}}:

LGRPO=1N[min(πθ(atst)πref(atst)At,clip(πθ(atst)πref(atst),1ϵ,1+ϵ)At)βDKL(πrefπθ)]L_{\mathrm{GRPO}} = - \frac{1}{N} \sum \left[ \min\left( \frac{\pi_\theta(a_t|s_t)}{\pi_{\mathrm{ref}}(a_t|s_t)} A_t, \text{clip}\left(\frac{\pi_\theta(a_t|s_t)}{\pi_{\mathrm{ref}}(a_t|s_t)}, 1-\epsilon, 1+\epsilon\right) A_t \right) - \beta D_{\mathrm{KL}}(\pi_{\mathrm{ref}} \| \pi_\theta) \right]

To mitigate catastrophic forgetting of optimal topologies, Self-Imitation Learning (SIL) is implemented with a Replay Buffer that stores sequences that outperform a dynamically calculated energy floor. The buffer's capacity decays exponentially during the final 75% of training, forcing the agent to exploit increasingly narrower bands of strong sequences. A linearized SIL loss maximizes the log-likelihood of recreating these trajectories:

ASIL=2.0max(Rseqfloor,0.1)A_{\mathrm{SIL}} = 2.0 \cdot \max(R_{\mathrm{seq}} - \text{floor}, 0.1) LSIL=Es,aBuffer[logπθ(as)ASIL]L_{\mathrm{SIL}} = - \mathbb{E}_{s,a \sim \mathrm{Buffer}} \left[ \log \pi_\theta(a|s) \cdot A_{\mathrm{SIL}} \right]

6. Architectural ablation studies

To isolate the driving mechanics behind the proposed architecture's convergence, several ablation configurations were systematically tested against the proposed model. The proposed model had the following architectural decisions: GRPO with Linear SIL loss, with an action mask, static temperature spiking, a shrinking replay buffer, and discounted future rewards. The following variants were evaluated:

  • Algorithm baselines: The proposed GRPO framework was compared against a pure Proximal Policy Optimization (PPO) implementation, as well as a pure GRPO implementation with the SIL loss completely removed.
  • Advantage scaling: The linear scaling of the SIL advantage function was compared against an exponential scaling variant to observe gradient stability (measured via KL divergence).
  • Search constraints: The model was trained with and without the dynamic action mask to evaluate its impact on pruning redundant gate matrices from the search space.
  • Reward dynamics: The discount factor was ablated (γ=1\gamma=1) to test the necessity of evaluating sequences based on global vs. immediate energy drops.
  • Exploration dynamics: The static temperature spiking mechanism was compared against a global temperature decay schedule. Additionally, the dynamically shrinking SIL buffer was compared against a static-capacity buffer to evaluate final-stage exploitation mechanics.

Results and ablations

The proposed Reinforcement Learning Generative Quantum Eigensolver (RL-GQE) achieved a minimum evaluation energy of E=7.5173E = -7.5173 against the theoretical ground state of E0=7.8284E_0 = -7.8284. While this represents an approximate 96% relative accuracy, the difference of 0.31 indicates that while the agent successfully navigated the macroscopic energy landscape to find an optimal discrete circuit topology, it halted in a local minimum. To validate the architectural decisions of this model, a series of ablation studies were conducted, isolating the impact of the policy algorithm, advantage scaling, search constraints, and exploration temperature dynamics.

Generative distribution convergence

The plots below track the probability distribution of the generated circuit energies across the training lifecycle at epochs 20, 60, and 260. During early epochs, the autoregressive transformer generates sequences with a wide variance, approximating random walk exploration across the Hilbert space. As the Group Relative Policy Optimization (GRPO) and Self-Imitation Learning (SIL) gradients take effect, the probability mass visibly shifts leftward, focusing on single topology evaluations. By the final epochs, the virtually 0 variance and resulting energy evaluation indicates that the transformer has successfully mapped the optimal topologies and assigns near-zero probability to highly suboptimal gate sequences. This distribution shift confirms that the agent is not relying on isolated random discoveries, but has structurally learned the target energy landscape.

Generator energy distribution at epoch 20, showing wide variance across the Hilbert space
Generator energy distribution at epoch 60, with probability mass collapsing toward lower energies
Generator energy distribution at epoch 260, concentrated on a single low-energy topology

Baseline algorithm performance

The proposed architecture was evaluated against two baselines: a Proximal Policy Optimization (PPO + SIL) agent, and a pure GRPO agent with the SIL buffer completely ablated.

Evaluation energy over training for the proposed GRPO + SIL model against PPO + SIL and pure GRPO baselines

As demonstrated in the above plot, the PPO formulation proved insufficiently sample-efficient for deep exploration in the difficult energy landscape, resulting in early stagnation. Conversely, the pure GRPO agent (Ablated SIL) suffered from severe catastrophic forgetting. While GRPO effectively navigated the landscape to find transient minimums, it lacked a mechanism to permanently anchor the policy to these discoveries, causing the evaluation energy to regress during subsequent generation batches. The combination of GRPO and SIL was strictly necessary to lock in performant topologies and force the policy to converge downward monotonically.

Architectural components and search constraints

Convergence efficiency was measured via the log percentage relative error between the predicted minimum and the theoretical ground state.

Log percentage relative error for the proposed model against the ablated action mask, ablated discount factor and static SIL buffer variants

The plot above illustrates the divergence caused by these architectural constraints. Ablating the dynamic action mask resulted in immediate performance degradation since redundant operations (e.g. placing consecutive RXR_X gates on the same qubit) expanded the Hilbert search space exponentially. This trapped the agent in suboptimal, high-depth circuits.

Furthermore, the temporal discounting of future returns was tested by setting γ=1.0\gamma=1.0 (Ablated Discount Factor). By removing the temporal discount, the model weighed immediate gate placements and distant future gate placements equally. This impaired the advantage estimator's ability to perform accurate credit assignment, failing to prioritize high-impact, early-sequence optimization and stalling convergence. Finally, replacing the dynamically shrinking SIL buffer with a static buffer prevented the agent from converging on the most performant trajectories. The shrinking capacity was vital for the final exploitation phase, artificially raising the standard for trajectory retention as training progressed.

Reward shaping and gradient stability

The mechanics of the SIL anchor rely heavily on how the advantage of past performant trajectories is scaled. The proposed model utilizes a linearized scaling function, ASILRseqfloorA_{\mathrm{SIL}} \propto R_{\mathrm{seq}} - \text{floor}. This was tested against an Exponential SIL scaling variant, ASILexp(Rseq)A_{\mathrm{SIL}} \propto \exp(R_{\mathrm{seq}}), to observe how aggressively the policy should be pulled toward the buffer.

KL divergence between the active and reference policies for linear versus exponential SIL advantage scaling

The plot above shows the gradient stability of these two approaches by measuring the Kullback-Leibler (KL) divergence between the active policy πθ\pi_\theta and the reference policy πref\pi_{\mathrm{ref}}. While both SIL scaling methods exhibited oscillating spikes, the KL divergence of the Exponential variant was more erratic and elevated throughout training. This aggressive scaling over-weighted outlier trajectories, leading to fractured policy updates and caused destructive weight modifications with the resulting energy being 7.408-7.408 compared to 7.517-7.517 for the proposed model. The linearized SIL scaling, though erratic, maintained generally lower KL divergence, ensuring more stable and reliable policy updates throughout the training.

Temperature exploration and exploitation

Since the energy landscape of spin Hamiltonians contains numerous regions of vanishing gradients and local minima, the agent requires a mechanism to periodically force exploration. The proposed model utilizes a dynamic simulated annealing approach, injecting temperature "spikes" during generation when the minimum energy stagnates for a set number of epochs (stagnation_epochs = 10).

Evaluation energy and sampling temperature over training for constant-magnitude spikes versus a linearly decaying schedule

As shown in the plot above, the proposed method utilizes constant-magnitude temperature spikes, which was ablated against a linearly decaying temperature schedule. Both methods incorporated a secondary temperature stepping mechanism: if an initial spike failed to discover a lower energy state, the temperature incrementally stepped upward to force continued exploration. The linear decay schedule was tested under the hypothesis that large, sudden temperature spikes late in long training sessions might disrupt the model's ability to exploit and refine a specific local minimum, and that decreasing the spike magnitude would yield lower energy circuits. However, this decay stunted exploration prematurely. As the agent encountered increasingly complex local minima mid-training, the weakened spikes became insufficient to eject the model from suboptimal valleys. Furthermore, the model's evaluation confidence deteriorated under the decay schedule, evidenced in the upper plot by a widening mismatch between the mean and minimum evaluation energies. This gap indicates a fractured deterministic policy that generated widely varying circuit topologies. Conversely, maintaining constant-magnitude spikes proved critical for enabling the agent to continuously break out of deep local minima. The robustness of this proposed approach is reflected in the agent's high evaluation confidence. The close alignment of the mean and minimum energies confirms that the model reliably and deterministically produced consistent, low-energy circuits.

Discussion

By achieving approximately 96% accuracy relative to the theoretical ground state, the model proved capable of solving the complex combinatorial problem of circuit topology. Because this study explicitly constrained the agent to a discrete action space of fixed-angle Pauli rotations, the convergence to this deep local minimum confirms that an optimal macroscopic structure was indeed found.

Furthermore, the baseline algorithm comparisons reveal a critical insight regarding policy optimization in quantum environments. The results demonstrate that Group Relative Policy Optimization (GRPO) is mathematically superior to standard value-based algorithms like Proximal Policy Optimization (PPO) or DQN in this domain. Quantum optimization landscapes are notoriously plagued by barren plateaus and highly non-convex local minima. In such environments, absolute value estimation, relying on a parameterized Critic network to predict expected returns, struggles to accurately model the chaotic reward surface. This leads to the early stagnation observed in the PPO baseline. GRPO circumvents this limitation entirely by eliminating the Critic network. Instead, it evaluates advantages through the relative comparison of trajectories generated within the same batch. This relative trajectory comparison intrinsically normalizes the severe variance inherent to quantum circuits, providing highly stable, directional gradient signals without the instability of a diverging value function.

While GRPO provides superior directional gradients, it must be coupled with Self-Imitation Learning (SIL) to prevent catastrophic forgetting. Because RL is known for sparse rewards, having this replay buffer mechanism was crucial to the success of the proposed architecture. The ablation studies demonstrated that without a mechanism to store and imitate previous performant trajectories, the policy rapidly unlearned optimal configurations. However, integrating this memory safely requires strict mathematical bounds. The KL divergence analysis proved that a linearized SIL scaling is necessary since the exponential scaling over-weighted outlier trajectories and violently fractured the policy's trust region, ultimately destroying the learned weights. Similarly, the structural ablations confirmed that the exponentially large Hilbert space must be physically constrained. This becomes a concern especially when scaling to higher qubit counts. Implementing a dynamic action mask to forbid mathematically redundant operations, and utilizing a temporal discount factor (γ=0.8\gamma=0.8) to optimize long-term sequence building, proved necessary to keep the RL search space tractable.

Finally, the analysis of exploration temperature dynamics highlights the unique challenges of navigating spin Hamiltonians. The failure of the decaying temperature schedule reveals that the agent encounters increasingly deep local minima even in the late stages of training. Maintaining constant-magnitude temperature spikes proved essential for providing the sudden energetic momentum required to escape these valleys, confirming that aggressive exploration mechanics must be maintained until a hard exploitation cutoff is explicitly triggered.

Ultimately, this study validates that an autoregressive policy can autonomously construct highly optimized quantum circuits via RL, utilizing GRPO and SIL to navigate complex energy landscapes without classical supervision. However, the current model acts as a static solver because it has only learned the topological rules for a single, specific 1D Isotropic Heisenberg XXX Hamiltonian. To realize the full potential of the RL-GQE as a universal combinatorial optimization tool, the framework must be expanded to understand and contextualize varying physical systems. This understanding naturally dictates the required architectural advancements, forming the basis for future research into generalized, context-aware quantum solvers.

Future work

The immediate next step is to quantify the true chemical accuracy of the currently generated topologies. Because this study explicitly omitted continuous variable tuning to isolate the RL agent's topological performance, the generated ansatz structures must be subjected to a standard VQE continuous parameter optimization phase, for example using gradient descent or Adam (Kingma & Ba, 2014). This will likely confirm that the remaining 4% energy gap is simply an artifact of utilizing discrete, fixed-angle Pauli rotations.

Following this validation, the primary objective is to implement a conditional encoder. To transition the autoregressive transformer into a context-aware solver, the Hamiltonian must be tokenized and embedded as conditioning information prior to generation. A promising approach involves encoding the target Hamiltonian matrix via a Graph Neural Network (GNN) or a flattened tensor representation, allowing the policy to dynamically adjust its generation strategy based on the specific physical parameters (e.g. coupling constants and magnetic field strengths). Once the encoding mechanism is integrated, the training pipeline will be expanded to evaluate in-distribution generalization. The agent will be trained across a parameterized distribution of Hamiltonians, such as a generalized 1D Heisenberg or Transverse Field Ising Model, where the parameters JJ and hh are dynamically sampled at each epoch. The goal is to observe whether the agent can learn the underlying physics of the family of Hamiltonians, adjusting its generated topology on-the-fly in response to the encoded inputs.

This contextual capability represents the critical frontier in generative quantum computing. Current literature demonstrates significant limitations: conditional models like GQCO are heavily specialized for classical combinatorial optimization on Ising models, while pure-RL architectures have largely been confined to the generation of simple stabilizer states. The robust gradient stability proven by the GRPO + SIL architecture in this study provides the theoretical foundation required to bridge this gap. If the RL-GQE can successfully condition its policy on embedded physical parameters, it has the potential to achieve zero-shot generation of highly optimized ansätze for heavily frustrated spin landscapes, eliminating the need for classical pre-computation entirely. Crucially, by bypassing the exponential data-generation bottleneck inherent to supervised learning, this framework is structurally positioned to scale to higher qubit counts. Ultimately, this paves the way for tackling complex eigenvalue problems within classically intractable Hilbert spaces, advancing toward practical quantum advantage.

The implementation is available at github.com/Mindbeam-AI/Generalized-GQE.

Accelerate LLM pre-training with Litespark.

Unlock higher throughput, lower energy use, and seamless integration with your existing stack.

Book a Demo