Efficient ground state preparation in variational quantum eigensolver with symmetry-breaking layers

Experimental progress in controlling quantum systems has allowed for the first quantum computational advantage claim¹ in recent years. This achievement initiated the race between noisy intermediate-scale quantum² (NISQ) computers and the state-of-the-art classical algorithms.³ Now, a NISQ device with hundreds of qubits is available,⁴ and those devices are expected to solve a practical computational problem beyond the reach of classical computers. Among possible applications, a variational quantum eigensolver (VQE)^5,6 that solves the ground state problem of quantum many-body Hamiltonian has gained much attention recently (see Ref. 7 for a recent review).

The VQE combines a parameterized quantum circuit and a classical optimization algorithm: A quantum circuit evaluates the expectation value of the Hamiltonian and its derivatives for the output quantum state, whereas a classical optimizer finds better parameters that minimize the energy. As solving the ground state problem of quantum Hamiltonians is difficult for classical computers, one may easily get a possible advantage of VQEs. Still, it is less understood which Ansatz and classical optimizer should be used to solve a given Hamiltonian efficiently.

One of the most widely studied Ansätze for solving many-body spin Hamiltonians is the Hamiltonian variational Ansatz (HVA).^8,9 Inspired by a short-depth quantum algorithm for solving combinatory optimization problems,¹⁰ the HVA is constructed using rotating gates whose generators are the terms of the Hamiltonian, which resembles the Suzuki–Trotter decomposition.¹¹ Even though the HVA solves many different Hamiltonians reliably,^12–14 however, it is not necessarily good at solving problems with a symmetry-broken nature, which prevails in many-body systems. Indeed, there are local Hamiltonians whose ground states cannot be generated by this Ansatz in a constant depth, albeit such a circuit exists, as the circuit obeys the same symmetry as the Hamiltonian.¹⁵ 

In this paper, we devise a symmetry-breaking Ansatz and explore its power for solving many-body Hamiltonians. We construct our Ansatz by adding symmetry-breaking layers to the HVA. The power of those symmetry-breaking Ansätze is numerically tested using the transverse-field Ising (TFI) and the transverse-field cluster (TFC) models, which have the $Z2$ and $Z2×Z2$ symmetries, respectively. Our results show that the symmetry-breaking Ansatz can find the ground state in a constant depth when the Hamiltonian has symmetry-broken ground states, whereas the bare HVA requires a linear depth for solving the same problem. We also show that, by adding a symmetry penalizing term to the loss function, one can choose a particular symmetry-broken state among degenerate ground states.

From the above construction, it is clear that the HVA describes the ground state when D is sufficiently large.

However, it is still unclear whether the depth provided by the HVA is optimal (up to a constant factor) for describing the ground state. In this paper, we show that the answer is negative using symmetry-broken systems.

Quantum many-body systems often break a symmetry characterized by a group. When it happens, a system has more than two ground states connected by a group generator. For example, if the $Z2$ symmetry is broken, we can find two ground states |GS₀⟩ and |GS₁⟩ such that |GS₀⟩ = G|GS₁⟩, where $G≠1$ is a generator of this symmetry (or a representation of it, more generally), which satisfies $G2=1$^16,17 (see also Ref. 18 for a general introduction).

Conventionally, the symmetry-breaking phenomena are understood using Landau’s theory.¹⁹ In this theory, symmetry-broken ground states are determined by an order parameter defined by an expectation value of a local observable O. For the $Z2$ symmetry-broken ground states, one has an order parameter O₀ = ⟨GS₀|O|GS₀⟩ and the other has O₁ = ⟨GS₁|O|GS₁⟩ = −O₀. Contrarily, symmetry-preserving ground states have ⟨GS|O|GS⟩ = 0. Thus, one can diagnose whether the system breaks the symmetry by computing the order parameter.

A phase transition occurs when the ground state breaks (or restores) a new symmetry as the control parameter of the Hamiltonian changes. For example, the TFI given in Eq. (2) has the strength of the external magnetic field h as a parameter. The Hamiltonian has the $Z2$ symmetry under the spin-flip; it does not change under the transformation Z_i ↔ −Z_i. However, the symmetry of the ground state is broken when |h| < 1, i.e., it has two symmetry-broken ground states when |h| < 1 in the thermodynamic limit. In Fig. 1, we plot the order parameter ⟨Z_i⟩ as a function of h²⁰ for one of the ground states. One finds that ⟨Z_i⟩ ≠ 0 for h < 1, whereas ⟨Z_i⟩ = 0 for all h > 1. This implies that a phase transition occurs at h = 1, and the symmetry is broken when h < 1.

In general, true symmetry breaking is only observed in the thermodynamic limit, i.e., a system often has a symmetry-preserving unique ground state for any finite N even when it has degenerated ground states in the thermodynamic limit. Thus, we need a diagnostic that can tell whether a system breaks symmetry using the ground states for finite-size systems. Correlation functions are widely used for that purpose. For example, in the 1D system, consider local observables {O_i} where each O_i supports site i. Then, the correlation function of the ground state is given by C(i, j) = ⟨O_iO_j⟩ − ⟨O_i⟩⟨O_j⟩, where ⟨⋅⟩ is the expectation value for the ground state from a system of size N. If C(i, j) converges to a constant as |j − i| → ∞ (we also need to increase N accordingly), the system has symmetry broken ground states in the thermodynamic limit.¹⁸ 

With the advent of topological physics,^21–23 however, systems with degenerate ground states related by a symmetry that an order parameter or correlation functions cannot diagnose are discovered. Such an example includes the AKLT model,²⁴ which has been extensively studied over decades. The model has two ground states in the open boundary condition, but the correlation function decays exponentially, and none of the local observable can distinguish the two states. This type of phase is named symmetry-protected topological (SPT) phase.²⁵ 

To incorporate both the conventional (symmetry-broken) and the SPT phases under a single theory, a modern definition of quantum phases has been introduced.^26,27 This theory is technically involved, but in essence, it tells that two Hamiltonians (H₀ and H₁) are in the same phase if there exists a local and gapped parameterized Hamiltonian H(τ) (0 ≤ τ ≤ 1) such that H(0) = H₀, H(1) = H₁, and H(τ) has the same symmetry for all 0 ≤ τ ≤ 1 (i.e., there exists an operator P such that [P, H(τ)] = 0). Here, we say H(τ) is gapped when the difference between the lowest and the first-excited eigenvalues of H(τ) is bounded below by a constant regardless of the system size N. This definition generalizes Landau’s theory since the gap closes in the thermodynamic limit during symmetry-breaking phase transitions.

Since the adiabatic evolution of a gapped Hamiltonian is described by a constant-depth quantum circuit,²⁸ the above definition is directly translated into a language of quantum circuits. Let us assume that Hamiltonians H₀ and H₁ have ground states |ψ₀⟩ and |ψ₁⟩, respectively. Then, |ψ₀⟩ and |ψ₁⟩ are in the same SPT phase if there exists a unitary operator U such that (1) U can be expressed by a local and constant depth quantum circuit, (2) |ψ₁⟩ = U|ψ₀⟩, and (3) [U, P] = 0, where P is the symmetry generator.

If we drop the symmetry restriction from the definition of the SPT phases, we arrive at the definition of topological phases.²⁹ Namely, two ground states are in different topological phases if we cannot find any local constant-depth circuit U that connects the two states. It is also known that true topological phases do not exist in the one-dimensional spin systems, i.e., one can always find a constant-depth circuit connecting between ground states of any gapped 1D local Hamiltonians.^26,27,30,31 In the following section, we will use such a distinction between the SPT and topological phases to design an efficient quantum circuit Ansatz.

The first model we study is the TFI model for N qubits defined in Eq. (2), where we impose the periodic boundary condition Z_N+1 = Z₁. This model has two distinct (ferromagnetic and paramagnetic) phases depending on the strength of h protected by the spin-flip symmetry P = ∏_iX_i.³² The critical point of this model h = 1 is well known.²⁰ This implies that if there is a circuit that commutes with P and connects two ground states in different phases (i.e., ground states for h > 1 and h < 1), the circuit depth must be larger than a constant.²⁶ On the other hand, a finite-depth circuit that connects two different ground states exists if we do not restrict such symmetry to a circuit since the system is gapped unless h = 1.^26,31

Symmetry-preserving/breaking Ansätze with the same D [Eqs. (6) and (9)] have different numbers of total gates. However, the number of entangling gates of both circuits is the same. When implemented on NISQ devices, the latter dominates the quality of the output states since the expected fidelity of entangling gates is much lower than single-qubit gates (see, e.g., Refs. 34 and 35). Thus, in the following discussion, we mainly compare the results from these Ansätze with the same D.

To observe the effect of symmetry-breaking layers, we simulated the VQE for the TFI using the bare HVA and our symmetry-breaking Ansatz given in Eq. (9). We optimize parameters of the Ansatz using the quantum natural gradient:^36–38 For each epoch t, we update parameters as $θt+1=θt−η(F+λt1)−1∇θ⟨ψ(θt)|H|ψ(θt)⟩$⁠, where η is the learning rate, $F=(Fij)$ is the quantum Fisher matrix, and λ_t is a (step dependent) regularization constant. We choose this optimizer as it works more reliably in solving the ground problem for both classical neural networks³⁹ and VQEs.¹³ For the quantum Fisher matrix, we use the centered one $F=Fijc=Re{⟨∂θiψ(θt)|∂θjψ(θt)⟩−⟨∂θiψ(θt)|ψ(θt)⟩⟨ψ(θt)|∂θjψ(θt)⟩}$³⁷ mostly (unless otherwise stated), but the uncentered one $Fijnc=Re{⟨∂θiψ(θt)|∂θjψ(θt)⟩}$³⁶ is also considered when it improves the performance. The difference between the two can be understood using the notion of the projected Hilbert space.³⁸ For the hyperparameters, we typically use η = 0.01 and λ_t = max(100.0 × 0.9^t, 10⁻³).

We present the optimized normalized energies $Ẽ=(EVQE−EGS)/EGS$ for different N and D in Fig. 2 when h = 0.5. We have used the Fisher matrix with the centering term $Fc$⁠, and initial values of parameters {ϑ_j, κ_j ϕ_j} are sampled from the normal distribution $N(0,σ2)$ with σ = 0.001 besides D = 3, 5 at Fig. 2(b). When D = 3 (all N) and 5 (for N ≥ 12), this initialization does not give any better energy than the Ansatz without symmetry-breaking layers. We instead found that initializing parameters for the symmetry breaking layers {ϕ_j} with samples from $N(2π/D,σ2)$ finds better optima in these cases. However, this initialization does not change the results for D = 4 and performs even worse when D ≥ 6 (see the  Appendix for detailed comparisons).

Figure 2(a) shows that the ground state is only found for D ≥ N/2 when symmetry-breaking layers are absent, which is consistent with Refs. 12 and 13. On the other hand, results with symmetry-breaking layers [Fig. 2(b)] clearly demonstrate that converged energies are significantly improved for 5 ≤ D < N/2. Most importantly, converged normalized energies are $≤10−7$ for all N when D ≥ 9, which implies that our Ansatz finds a constant-depth circuit for solving the ground state. In addition, the results show that there is a finite-size effect up to N = 18 where the accurate ground state is only obtained when D ≥ N/2.

However, the converged energies and the true ground state show a large discrepancy when D = 4. In fact, the results for D = 4 are even worse than those of D = 3, which signals that the optimizer gets stuck in local minima. This type of convergence problem is already observed in Ref. 13. To obtain a better convergence for D = 4, we employ a transfer learning technique. Instead of starting from a randomly initialized circuit, we insert a block in the middle of the converged circuit from D = 3 and perturb all parameters by adding small numbers sampled from $N(0,σ′2)$ (where we typically use σ′ = 0.01). We then optimize the full circuit using the quantum natural gradient.

We show converged energies from the transfer learning in Fig. 3. Figure 3(a) illustrates that the transfer learning succeeds in finding a better optimum for D = 4. However, the result for (N, D) = (10, 5) is worse than that of the random initialization, which implies that the transfer learning does not necessarily find the global optima. We additionally plot a learning curve for (N, D) = (18, 4) in Fig. 3(b) where the initial energy is much higher due to an added perturbation, but it eventually finds a quantum state with lower energy.

In summary, our Ansatz does not find a constant-depth circuit for solving the ground state of the TFC based on the results up to N = 20. We attribute this to the finite-size effect and expect that a much larger system size N (beyond our computational capacity) is required to take advantage of symmetry-breaking Ansatz.

The main idea is minimizing the expectation value of O = H_cluster + α₁P₁ + α₂P₂ for suitable α₁ and α₂, instead of the Hamiltonian itself. We expect the obtained state to be an eigenstate with the corresponding eigenvalue ±1 of P_i when sign(α_i) = ∓1.

For example, we study the VQE for preparing the ground state with P₁ = P₂ = −1 using (N, D) = (14, 8). We have tested VQEs with the varying learning rates η ∈ [0.01, 0.025, 0.05, 0.1] and α₁ = α₂ = α ∈ [0.1, 0.3, 0.5, 1.0, 2.0, 4.0] and found that the desired state can be prepared by the best for η = 0.025, α = 2.0. In particular, all VQE runs have converged to ⟨P₁⟩ = ⟨P₂⟩ = −1 (accuracy within 10⁻⁸) with these values of η and α. The resulting learning curve is shown in Fig. 5(a). However, the converged energy $Ẽ≈1.502×10−3$ for D = 8 is significantly worse than the results from other models we have studied above (we still note that this is comparable to results from the Adam optimizer, which gives normalized energies $≈10−2−10−4$⁠).¹⁴ We thus further study whether increasing D helps the convergence in Fig. 5(b). The results show that the converged energies are getting more accurate as we increase D. Precisely, 47 instances out of 48 independent VQEs runs have converged to $Ẽ≤10−8$ when D = 20. This contrasts the barren plateaus narrative⁴² that more expressive quantum circuits are prone to suffer the optimization problem.⁴³ We believe that the behavior of our circuits is related to the known fact in classical machine learning that over-parameterized networks converge better,⁴⁴ but we leave a detailed study for future work.

We studied the limitation of the bare HVA for solving a system with symmetry-broken ground states and proposed an Ansatz with symmetry-breaking layers to overcome this problem. Such a symmetry-breaking Ansatz could enhance the overall performance of the VQE. Especially for the transverse-field Ising model, we observed that our symmetry-breaking Ansatz finds a constant-depth circuit for the ground state, whereas the bare HVA requires the depth to be linear in the system size. We further proposed a technique to choose a specific symmetry broken ground state among possible ground states. This was possible by adding a penalizing term to the cost function.

We expect that one can observe similar behavior for continuous symmetry, e.g., U(1) or SU(2). However, given that those symmetries are not broken in 1D,⁴⁵ we should use 2D models to observe the effect of symmetry-breaking layers. Symmetry-broken 2D spin models are also physically interesting as they are related to the lattice gauge theory.⁴⁶ However, we leave the detailed study of 2D models for future work, as numerical simulation of such models is more demanding.

Our symmetry-breaking Ansatz can be implemented on real quantum hardware. In this case, one should decompose a gate into the set of gates supported by a target quantum device. For example, a gate $Lzxz(θj)$ used for the transverse-field cluster model is not naturally supported by most quantum computing architectures.^47–52 Instead, it can be decomposed into, e.g., $e−iθZi−1XiZi+1=CZ1,2CZ2,3e−iθXiCZ1,2CZ2,3$⁠, which may increase overall circuit depths. Still, as our symmetry-breaking Ansatz can provide a scaling advantage, there is a certain number of qubits that our circuit can solve the problem in a shorter depth than the bare HVA even after such a decomposition.

Incoherent noise can be a significant problem when implemented on noisy quantum hardware. While a weak incoherent noise may help to find a symmetry-broken ground state when the target state is not entangled,⁵³ we generally expect that the noise destroys the quantumness of the circuit, and the circuit outputs a decohered state that is far from a true ground state.⁵⁴ As the circuit depth is the main factor limited by this type of noise, our symmetry-breaking Ansatz may again have an advantage over the bare HVA when the noise is considered.

The author thanks David Wierichs for helpful discussions. This project was funded by the Deutsche Forschungsgemeinschaft under Germany’s Excellence Strategy—Cluster of Excellence Matter and Light for Quantum Computing (ML4Q) EXC 2004/1-390534769 and within the CRC network TR 183 (Project Grant No. 277101999) as part of project B01. Numerical simulations were performed on the JUWELS cluster at the Forschungszentrum Juelich.

The authors have no conflicts to disclose.

Chae-Yeun Park: Resources (equal); Software (equal); Supervision (equal); Writing – original draft (lead).

The source code for the current paper can be found in the Github repository.⁵⁵ 

When the Ansatz contains symmetry-breaking layers, an optimization algorithm easily gets stuck in local minima. We here compare several different setups for the transverse field Ising model we have studied in Sec. III. We show results from the transverse Ising model with different initial values and using the centered and uncentered Fisher matrix in Fig. 6. We can see that the result with {ϕ_j} initialized from $N(0,σ2)$ and using the centered Fisher matrix is the most natural. However, when D = 3, 5, the results with {ψ_j} initialized from $N(2π/D,σ2)$ show better convergence.

We have also found that using the uncentered Fisher matrix improves convergence for D = 4 where the centered Fisher matrix failed to find an appropriate optimum. Our results suggest that the learning landscape of the VQEs with symmetry-breaking layers is rugged, especially when the parameters are not sufficient to describe the ground state accurately.