Bridging physical intuition and hardware efficiency for correlated electronic states: the local unitary cluster Jastrow ansatz for electronic structure

Abstract

A prominent goal in quantum chemistry is to solve the molecular electronic structure problem for ground state energy with high accuracy. While classical quantum chemistry is a relatively mature field, the accurate and scalable prediction of strongly correlated states found, e.g., in bond breaking and polynuclear transition metal compounds remains an open problem. Within the context of a variational quantum eigensolver, we propose a new family of ansatzes which provides a more physically appropriate description of strongly correlated electrons than a unitary coupled cluster with single and double excitations (qUCCSD), with vastly reduced quantum resource requirements. Specifically, we present a set of local approximations to the unitary cluster Jastrow wavefunction motivated by Hubbard physics. As in the case of qUCCSD, exactly computing the energy scales factorially with system size on classical computers but polynomially on quantum devices. The local unitary cluster Jastrow ansatz removes the need for SWAP gates, can be tailored to arbitrary qubit topologies (e.g., square, hex, and heavy-hex), and is well-suited to take advantage of continuous sets of quantum gates recently realized on superconducting devices with tunable couplers. The proposed family of ansatzes demonstrates that hardware efficiency and physical transparency are not mutually exclusive; indeed, chemical and physical intuition regarding electron correlation can illuminate a useful path towards hardware-friendly quantum circuits.

---

1 Introduction

Obtaining the ground-state energy of a correlated electronic system has been recognized as an important field of application of quantum computers.^1–6 In recent years, various algorithms have delivered promising results in the calculation of potential energy curves, ground- and excited-state energies and correlation functions for various molecular systems.^7–10 Some prominent examples of such algorithms include the variational quantum eigensolver (VQE),^6,11 quantum subspace expansion (QSE),^12,13 imaginary-time evolution,¹⁴ adiabatic state preparation (ASE),¹⁵ and quantum phase estimation (QPE).¹⁶ The feasibility of such algorithms on near-term quantum devices depends on the ability of approximate wavefunction ansatzes to yield accurate energies and properties while requiring the smallest possible amount of quantum resources.

The most widely used ansatz for eigenstates of electronic Hamiltonians is the unitary coupled cluster form with single and double excitations (qUCCSD):¹⁷

|ψ〉 = e^T̂−T̂†|Φ₀〉 (1)

where |Φ₀〉 is a reference determinant, often the Hartree–Fock state, and As with classical (projective) coupled cluster methods,¹⁸ qUCCSD preserves several desirable properties that are advantageous for the computationally efficient modeling of chemical reactions including size consistency, size extensivity, differentiability, and invariance under occupied–occupied and virtual–virtual orbital rotations. The exponential ansatz also provides a compact, approximate description of infinite-order electronic excitations, critical for even qualitatively correct physical descriptions of, e.g., many-body screening and polarization effects.¹⁹ While less studied relative to classical coupled cluster models, the unitary variant is better suited to describe some static correlation effects beyond two-electron problems when solved variationally.

The qUCCSD ansatz is especially suited for quantum computers, because the cost of exactly (i.e., variationally) evaluating the energy E = 〈Ψ|Ĥ|Ψ〉 over a qUCCSD state Ψ is polynomial with system size on a quantum computer and factorial on a classical one.²⁰ Despite such favorable polynomial scaling, the qUCCSD ansatz contains parameters, where N_o and N_v are the numbers of occupied and virtual orbitals in a mean-field reference wavefunction, and a quantum circuit of depth where N = N_o + N_v and N_t is the number of Trotter steps used to approximately implement the qUCCSD ansatz. As a result, quantum simulations based on the qUCCSD ansatz vastly exceed the capabilities of contemporary quantum devices in the noisy intermediate-scale quantum (NISQ) era in terms of gate fidelity, qubit coherence time, and qubit connectivity.

The high computational cost of qUCCSD calculations originates from the rapidly growing number of double-excitation amplitudes. However, chemical intuition (based roughly on the locality of most electronic correlations despite growing system size) suggests that more compact ansatzes can be designed, leveraging formal properties of the cluster operator. This observation has prompted various authors to introduce approximate flavors of qUCCSD including the low-depth ansatz,²¹ the low-rank approximation,²² and the k-fold product of a unitary pair coupled cluster.²³ Yet, further quantum resource reductions would be highly desirable.

An alternative approach to formulating a compact many-body wavefunction was introduced by Jastrow in 1955,²⁴ which exponentiates a function of the relative position of two particles, i.e., e^f(r_ij), such that relatively short-ranged two-body correlation effects are efficiently captured. This wavefunction is widely used in variational Monte Carlo for continuous systems^24–26 and lattice models^27,28 in both first^29–32 and second^33,34 quantization. More recently, the Jastrow has been used in the context of quantum computation.^35,36 The starting point for this study is the unitary cluster Jastrow ansatz,³⁶ which is a unitary variant of the cluster Jastrow form introduced by Neuscamman in the context of the antisymmetrized geminal power reference state.^33,34,37 This ansatz has been demonstrated to yield promising accuracy with reduced circuit depths in comparison to qUCCSD.^36,38

A challenging frontier of many-body quantum mechanics involves the description of strongly correlated fermionic systems. By “strongly correlated” we refer to quantum states in which more than one Slater determinant has significant weight in the wavefunction, a situation which results from the presence of nearly degenerate eigenstates. In the chemistry community, this is also known (equivalently) as “multi-reference character,” “static” or “non-dynamic” correlation. Strong correlation is challenging for most available quantum chemical methods. For instance, it defeats non-degenerate perturbation and response theories based on a single reference configuration, can require often prohibitively high orders of cluster operators in coupled cluster approaches, and is difficult to treat in a balanced fashion along with dynamic correlation (which can be thought of as short-range, instantaneous electron–electron repulsion). Molecular systems which exhibit strong correlation in ground (and many low-lying excited) states include di-/poly-radical organic chromophores, polynuclear transition metal compounds, and organometallic complexes with redox-noninnocent ligands.³⁹ Even simple molecules which are stable (closed-shell) at their equilibrium geometry can be pushed into regimes of strong correlation when bonds are stretched and before they are broken.

The unifying feature of most strongly correlated electronic states is that localized, aligned spins are recoupled through anti-ferromagnetic interactions into low-spin states of open-shell character. A catalytically relevant example includes the Cu₂O₂ functional unit,^40,41 in which each Cu²⁺ d⁹ ion has a localized unpaired electron which couples to form an open-shell singlet state. A simpler example, which nevertheless preserves much of the same physics, can be found in the intermediate region along the potential energy curve of stretched H₂, i.e. a diradicaloid singlet. In both cases, restricted Hartree–Fock (RHF) produces a qualitatively incorrect closed-shell electron configuration and unrestricted HF (UHF) suffers from spin contamination. To rectify this, we draw inspiration from the Hubbard Hamiltonian:

which contains a nearest-neighbor hopping term and an on-site repulsive interaction, which penalizes double-occupancy. This model has been used successfully to describe emergent effects due to strong electron correlation such as Mott insulator transitions^42–44 and high-temperature superconductivity in cuprates.^45–48

Herein we present a variational ansatz resembling a reference RHF state time-evolved under a Hubbard Hamiltonian. Starting from the unitary cluster Jastrow wavefunction, we keep only on-site, opposite-spin and nearest-neighbor, same-spin number–number terms. The result is a single-reference ansatz that can access open-shell character via a correlator which penalizes double-occupancy, and that captures strong (and weak) correlations in a systematically improvable way. The physically transparent form of the local cluster Jastrow family of ansatzes also provides a natural path to hardware efficiency: SWAP gates are no longer necessary for all qubit topologies (square, hex, and heavy-hex), and the ansatz can be implemented via a sequence of fSim gates that is of lower depth than that of qUCCSD, native to state-of-the-art quantum hardware utilizing tunable couplers.⁴⁹

2 Methods

2.1 Canonical ansatz

The unitary cluster Jastrow (UCJ) ansatz has the form³⁶ of L-fold products (or a product of L layers)whereare one-body operators. In eqn (4), p, q = 0…N − 1 label molecular spatial orbitals and σ, τ label spin polarizations (α, β for spin-up and spin-down electrons, respectively). K_pq^μ has complex matrix elements and is anti-Hermitian, whereas J_pq,στ^μ has real matrix elements and is symmetric. The UCJ ansatz commutes with the total spin operator S_z if

J_pq,αα^μ = J_pq,ββ^μ, J_pq,αβ^μ = J_pq,βα^μ, (5)

(note that commutation with the S² operator requires J_pq,αα^μ = J_pq,αβ^μ) leading to the expressionwhich is assumed throughout this work. The UCJ ansatz can be related to a twice-factorized low-rank decomposition of the qUCCD ansatz,^22,36 and the L-product form is such that the exact full configuration interaction wavefunction can be obtained viaeqn (3).^36,50

Despite the many desirable properties of the UCJ ansatz, its implementation on contemporary quantum devices is challenging due to limited qubit connectivity, gate fidelity, and qubit coherence. To understand this observation, let us consider the widely used Jordan–Wigner (JW) representation of electrons in N spatial orbitals onto a system of 2N qubits,

where (S₊)_p = (X_p + iY_p)/2 is the raising operator for qubit p, and X, Y, Z denote the 2 × 2 Pauli matrices. In the JW representation, the operator e^−K̂_μ can be implemented⁵¹ by using a quantum circuit of Givens rotations acting on adjacent qubits on a device with linear connectivity, using e.g. the design introduced originally by Reck et al.,⁵² or the more economical one proposed by Clements et al.⁵³

The operator e^iĴ_μ is mapped onto a quantum circuit of the form

where is a single-qubit gate andis a two-qubit gate. While the circuit in eqn (8) comprises two-qubit gates, it requires all-to-all qubit connectivity or the use of a fermionic swap network.⁵¹

2.2 Local ansatzes

The key idea of this work is to introduce a “local” approximation of the UCJ ansatz, which makes the following modifications for opposite-spin and same-spin number–number terms:In other words, we introduce a set of sparse matrices J_pα,qα^μ, J_pα,qβ^μ such that the quantum circuit eqn (8) only comprises gates v_p,q acting on adjacent qubits in the topology of a certain device,where S and S′ denote sets of qubits encoding the occupation of spin-up/down orbitals. To justify eqn (12), and define S and S′, it is useful to examine Fig. 1. It is natural to arrange qubits encoding the occupation of spin-up/down orbitals (represented by red/blue circles, respectively) along horizontal segments, so the operations e^K̂_μ can be implemented acting exclusively on neighboring qubits (connected by black lines). This qubit arrangement leads naturally to the choice

S′ = {(p,p + 1), p = 0…N − 2}, (13)

i.e. to retain same-spin number–number terms acting on adjacent orbitals (circles of identical colors connected by black lines). Similarly, we restrict opposite-spin terms to pairs of adjacent qubits (circles of different colors connected by thick black lines). As an example, consider the device with hexagonal-lattice topology shown in Fig. 1b: for even integers p = 2k, k = 0…⌊(N − 1)/2⌋, where N is the number of spatial orbitals and the p-th qubit of the spin-up line is connected with the p-th qubit of the spin-down line (by thick black lines connecting qubits of different colors), so that v_p,p+N gates can be applied to qubits with p ∈ S_hexagonal = {2k, k = 0…⌊(N − 1)/2⌋}. Similarly, for the device topologies sketched in Fig. 1, one has

Fig. 1 Device topologies considered in this work. Gray and colored circles indicate qubits, and black lines connect adjacent qubits in the topology of the device. Red and blue circles indicate qubits onto which, respectively, spin-up and spin-down molecular orbitals are mapped, and green circles indicate ancillae. In the case of the heavy-hex lattice (panel (c)), we consider both a “zig-zag” pattern (left) and a loop (top right). A number inside a red/blue circle denotes the index (p in the main text) of the spatial orbital whose occupation is encoded in the corresponding qubit, and thicker black lines connect (p, α) and (p, β) spin orbitals.

The resulting “local” ansatz will be abbreviated as LUCJ.

2.3 Quantum circuits

Fig. 2 provides a schematic implementation of a single layer of the ansatz written in eqn (3), for the case of a four orbital problem (8 spin-orbitals).

Fig. 2 Implementation of a layer of the Jastrow ansatz using a change-of-basis unitary exp(K_α + K_β) (green blocks), and a cluster operator comprising a same-spin term exp(iJ_αα + iJ_ββ) (orange and blue blocks) and an opposite-spin term exp(iJ_αβ) (red block), illustrated for N = 4 spatial orbitals.

The one-body blocks (Fig. 2 green), which correspond to orbital rotations within a spin sector (α or β), are composed of single-qubit R_z gates followed by a sequence of Givens' rotations, illustrated in Fig. 3. These are unaffected by any local approximations to the UCJ wavefunction.

Fig. 3 Implementation of a change of basis exp(K_α + K_β) as a composition of 2N R_z gates and N(N − 1) XX + YY gates G(θ, β) arranged in N layers (marked as G, omitting parameters to avoid clutter), illustrated for a system of N = 4 spatial orbitals.

Currently only trapped ion architectures have all-to-all connectivity. All other devices, e.g. superconducting qubits, have limited qubit connectivity. As a consequence, qubits must be swapped to realize gates acting on non-adjacent spin-orbitals. For example, in Fig. 4 we show an example of how the same spin (αα, orange/blue) and opposite spin (αβ, red) components of the e^iĴ part can be optimally implemented on a device with linear connectivity using a SWAP network.⁵⁴

Fig. 4 Implementation of a cluster operator exp(iJ_αα + iJ_ββ + iJ_αβ) with all-to-all ansatz connectivity, illustrated for N = 4 spatial orbitals. Linear device connectivity is assumed, and a generalized SWAP network is employed so that logical operations are performed on physically adjacent qubits only. Boxes connected by a vertical line denote a number–number gate U_nn(φ), and numbers on connected boxes denote the pair of qubits acted upon by a number–number gate, numbered from 0 to 2N − 1. Crosses connected by a vertical line denote a SWAP gate. Orange, blue, and red boxes denote terms of the J_αα, J_ββ and J_αβ operators, respectively.

In contrast, implementation of the LUCJ ansatzes shown in Fig. 5a–d for square, hex, and heavy-hex qubit topologies, respectively, does not require any SWAP gates. Furthermore, the opposite-spin number–number part of the LUCJ ansatz, for all qubit topologies, can be implemented in a single layer of gates, while the same-spin parts require a constant number of layers (in the present case, three). Thus, the circuit depth of the LUCJ ansatz associated with the exponential of the Ĵ operator no longer depends on the number of spatial orbitals, N (c.f.Fig. 5).

Fig. 5 Implementation of a local cluster operator exp(iJ_αα + iJ_ββ + iJ_αβ), illustrated for a system of N = 4 spatial operators. (Left) the opposite-spin operator is limited to U_nn gates acting on adjacent qubits in the topology of a square (a), hex (b), linear (c), and heavy-hex (d) lattice. For the heavy-hex lattice, the qubit interaction is generally mediated by an ancilla. (Right) (e) the same-spin operator is limited to U_nn gates compatible with linear connectivity, the same assumed in the implementation of change-of-basis unitaries.

For the heavy-hex qubit topology, as seen in Fig. 1c, qubits associated with spin-up and spin-down orbitals (red and blue circles) are generally connected by ancillas (green circles). Correspondingly, the circuit in Fig. 5d couples two qubits (marked as 0α and 0β) with an ancilla initialized in |0〉. To implement the unitary exp(−iφn_0αn_0β), a gate U_nZn(φ) = exp(−iφn_0αZ_an_0β) is applied to the three qubits. The ancilla can be measured and reset to zero, and post-selection can be carried out on the measurement outcomes to mitigate errors.

2.4 Connection with the Hubbard model

Jastrow wavefunctions are physically motivated ansatzes. For example, the Jastrow factor used in standard variational Monte Carlo calculations can be derived by propagating a mean-field wavefunction in imaginary time using the Feynman–Kac formula and a short-time approximation.^55,56 The UCJ wavefunction can be derived by treating a twice-factorized low-rank decomposition of the qUCCD ansatz²² as a variational wavefunction.³⁶ Here, we further interpret the LUCJ ansatz in terms of adiabatic state preparation (ASP).

The ASP approach starts with a Hamiltonian of a quantum system whose ground state, Φ₀, is easy to prepare (e.g. the Fock operator F̂). The Hamiltonian is then slowly varied toward the target Hamiltonian, Ĥ, whose ground state encodes the solution to the problem. If the change is slow enough, and the Hamiltonian remains gapped along the path, the ground state of Ĥ can be prepared. More specifically, a quantum state |Ψ_t〉 evolves in time according to the Schrödinger equation from time t = 0 to time t = T starting from |Ψ_t=0〉 = |Φ₀〉. We require that Ĥt = 0 = F̂ and Ĥ_t = T = Ĥ, which holds, e.g., when Under suitable conditions,¹⁵ the state |Ψ_T〉 converges to the ground state of Ĥ in the limit of T → ∞.

A connection between the LUCJ ansatz and ASP emerges by approximating the state |Ψ_T〉 with a first-order Trotter formula using M steps of length Δt = T/M,

The electronic structure Hamiltonian can be expressed as

i.e., the sum of a one-body operator and a sum of squares of Hermitian one-body operators L̂_γ, obtained, e.g., from a density fitting or Cholesky representation of the electron–electron repulsion. Approximating the exponential, e^{−iΔtĤ_kΔt}, in 15 with a first-order Trotter formula, we obtain

The UCJ functional form is retrieved by introducing the decomposition where L_p^γ denotes the eigenvalues of the operator L̂_γ, and e^K̂_γ implements a change of basis from the molecular-orbital basis to the eigenbasis of L̂_γ. The eigendecomposition of L̂_γ leads to a factorized form for the J tensors, It is known^57,58 that the operators L̂_γ have low rank, i.e. the eigenvalues L_p^γ of L̂_γ decay rapidly. In the extreme case of a rank-1 operator, one has

which is a local UCJ form (i.e., the J_pσ,qτ^kγ only has a non-zero element corresponding to p = q = 0) with L = MN_γ terms (i.e. one for each of the MN_γ matrices J).

Finally, we remark that eqn (17) contains the exponential of a one-body operator (a linear combination of Ĥ₁ and F̂), which is absent in the UCJ ansatz eqn (3). Indeed, the UCJ ansatz derives from a twice-factorized low-rank decomposition of an operator comprising double excitations only, making it akin to classical CCD (coupled-cluster doubles). Like CCD, LUCJ is not stationary under orbital rotations in the active space, which are described by the exponential of a one-body operator. Therefore, to ensure stationarity under orbital rotations, and in accordance with eqn (17), we append the exponential e^X̂ of a one-body orbital optimization operator X to the LUCJ ansatz. We also remark that recent work⁵⁹ has indicated that various quantum computing ansatzes are not stationary under orbital rotations, and that their flexibility and accuracy improves upon addition of an orbital optimization.

In summary, ASP is designed to find the ground state of a Hamiltonian Ĥ, and will do so if the simulation time T is long enough. A Trotter approximation of ASP consists of an alternation of operators with the same functional form as the UCJ (and LUCJ assuming low-rank operators) ansatz. This observation provides a connection between ASP and LUCJ, and a further physical interpretation of this ansatz, along with the twice-factorized low-rank decomposition of the qUCCD ansatz.^22,36

2.5 Computational details

All correlated calculations use the STO-6G minimal basis set. qUCCSD calculations were performed with the Qiskit program (version 0.39.0, and Qiskit nature version 0.4.0),⁶⁰ with restricted orbitals. LUCJ calculations were performed with an in-house classical implementation (detailed in the ESI†) using the BFGS optimizer (SciPy) and PySCF version 1.7.6.^61,62 All simulations in this work are carried out via linear algebra operations which do not consider the effects of shot noise or decoherence.

The bond-stretching coordinate of ethene, preserving a planar nuclear configuration (D_2h point group), was obtained from a series of constrained geometry optimizations, performed with Q-Chem⁶³ version 5.4.2. These optimizations used the ωB97X-V functional⁶⁴ and the def2-QZVPP basis.^65,66 To consistently choose the four active orbitals along the dissociation curve, at each geometry we generated unrestricted MP2 natural orbitals and selected those with a large overlap with the target atomic orbitals, i.e. carbon 2p_z and sp² orbitals.

Unless otherwise specified, when active spaces of canonical molecular orbitals were used, orbitals of increasing (i.e., less negative) orbital energy were assigned to qubits of increasing index. When unrestricted MP2 natural orbitals were correlated, orbitals of decreasing occupation number were assigned to qubits of increasing index. The orbital to qubit assignments for cyclobutadiene, stretched ethene, and benzene are shown explicitly in Fig. S4† for each topology, corresponding to the qubit indices shown in Fig. 1.

The variational optimization of LUCJ parameters is a rather delicate operation. For instance, initializing the LUCJ ansatz with K_μ = J_μ = 0 leads to a zero-valued energy gradient vector. Indeed, expanding the unitary operator (3) to first order in K_μ, J_μ leads to

and for the energy,

One can immediately see that the gradient with respect to K_μ is zero. Furthermore, when Φ₀ is the HF state, the operators Ĵ_μ preserve it, i.e. Ĵ_μ|Φ₀〉 = j_μ|Φ₀〉. In this situation,

〈Φ₀|[Ĥ,Ĵ_μ]|Φ₀〉 = 0 (21)

and thus, the parameter configuration K_μ = J_μ = 0 is a stationary point for the LUCJ energy.

In light of the above, we initialized the parameters of an LUCJ/all-to-all calculation from a truncated doubly factorized low-rank decomposition of CCSD or MP2 t₂ amplitudes. We then used a “bootstrapping” procedure, employing the converged LUCJ/all-to-all, LUCJ/square, and LUCJ/hex parameters as initial guesses for LUCJ/square, LUCJ/hex, and LUCJ/heavy-hex calculations, respectively. For single-point calculations, when initializing e.g. a LUCJ/hex calculation from parameters of a converged LUCJ/square calculation, we permuted orbitals so that p ∈ S corresponded to the largest values of |J_pp,αβ^μ| (see the ESI†). The LUCJ energy must decrease monotonically with L for a fixed connectivity and decrease monotonically with connectivity for a fixed L. We reinitialized calculations from parameters pertaining to different values of L and/or connectivity until we no longer observed deviations from this behavior.

3 Results

3.1 Single-point calculations for square cyclobutadiene, stretched ethene, and benzene

Cyclobutadiene is a highly strained, anti-aromatic molecule which exhibits a particularly challenging open-shell singlet ground-state at square geometries.⁶⁷ We investigate a transition state geometry corresponding to the rectangular to square coordinate,⁶⁸ wherein four carbon atoms form a square with a side length of 1.456 angstrom, with carbon–hydrogen bond lengths of 1.069 angstrom. We treat only the π-space, which consists of 4 electrons in 4 orbitals. While RHF is qualitatively wrong (predicts a triplet ground-state) and qUCCSD is beyond our target accuracy of 1.6 mHa vs. the exact total energy, Fig. 6 shows that all of the versions of the LUCJ ansatz can recover the exact energy. Specifically, this is accomplished by two layers (L = 2) for all-to-all and square approximations, by three layers for the hex ansatz, and by four layers for the heavy-hex ansatz. We note that double factorization of the qUCCSD t₂ amplitudes implies L = 8, which enables us to quantify the resource savings for each species (e.g. factors of 4 fewer layers for all-to-all and square). Taken together, we find that our LUCJ approach gives better accuracy than qUCCSD, and does so with lower circuit depths. We do not expect L = 1 to yield meaningful results, given that the singular value decomposition of qUCCSD amplitudes implies two layers for every singular value. In what follows, we focus on L ≥ 2.

Fig. 6 LUCJ (colored circles), qUCCSD (purple dot-dashed line) and FCI (black dashed line) ground-state energy for cyclobutadiene in a (4e,4o) active space. The SCF and qUCCSD energies are −153.169094 and −153.337275E_h, respectively.

Another model (4e,4o) problem which exhibits strong correlation is stretched ethene. At dissociation, each methylene has two non-bonding orbitals of p and sp² characters. As two well-separated methylenes are brought together to form ethene at its equilibrium geometry, these four orbitals form two bonding molecular orbitals of σ and π characters, along with the corresponding antibonding orbitals. Fig. 7 shows various calculations in an intermediate region of the double-bond breaking coordinate, at R = 2.0 angstroms. Unrestricted qUCCSD is shown in the dotted line, deviating from the exact energy by roughly 3 mHa. The all-to-all, square, and hex LUCJ models are converged to within 1.6 mHa of FCI by an L of 2, 3, and 4, respectively. The accuracy of the heavy-hex LUCJ model exceeds that of qUCCSD by L = 4. As was the case for square cyclobutadiene, the twice-factorized qUCCSD ansatz requires an L of 8, thereby requiring significantly more quantum resources to implement on the device.

Fig. 7 LUCJ (colored circles), qUCCSD (purple dot-dashed line) and FCI (black dashed line) ground-state energy for stretched ethene in a (4e,4o) active space. The SCF and qUCCSD energies are −77.445626 and −77.626017E_h, respectively.

Next we turn to the benzene molecule at equilibrium, and consider the (6e,6o) π space representation. Unlike square cyclobutadiene and stretched ethene, this system does not exhibit strong static correlation in the ground state. Indeed, the molecule possesses exceptional chemical stability due to its aromaticity (delocalized π system with 4k + 2 electrons, k = 1), and the wavefunction is qualitatively well-described by the RHF state. Nevertheless, there is substantial dynamic correlation (e.g., the deviation of MP2 from the exact FCI energy in the active space is roughly 31 kcal mol⁻¹). Fig. 8 shows the predictions of the LUCJ models alongside qUCCSD. Due to the absence of multireference characters, qUCCSD, as expected, gets rather close to the exact value, although classical CCSD is notably more accurate with a deviation of 0.2 kcal mol⁻¹ (the difference between qUCCSD and CCSD is due to approximations, e.g., Trotter errors incurred by the quantum implementation of the former). We note that SVD of the qUCCSD t₂ amplitudes implies L = 18. The all-to-all UCJ ansatz provides better accuracy than qUCCSD at and after L = 2; the square LUCJ model by L = 5. While comparable accuracy is obtained with the hex and heavy-hex LUCJ ansatzes by roughly 6 layers, these circuits are still relatively shallow compared to that of qUCCSD.

Fig. 8 LUCJ (colored circles), qUCCSD (purple dot-dashed line) and FCI (black dashed line) ground-state energy for benzene in a (6e,6o) active space. The SCF and qUCCSD energies are −230.130155 and −230.236428E_h, respectively.

3.2 Potential energy curves

Accurate, single-point correlation energy predictions computed at representative geometries along a reaction coordinate optimized, e.g., with density functional theory, can predict most thermochemical properties (such as bond dissociation energies, redox potentials, pK_a's, spin gaps, and reaction rates) with chemically useful accuracy. However, the ability to produce smooth potential energy curves which describe both the equilibrium region and the dissociation limit is also highly desirable, e.g., in simulations of chemical dynamics. In addition, pushing a chemical system at equilibrium through bond-breaking coordinates can reveal instances of strong correlation – the description of which is the motivating goal of our LUCJ ansatzes. In this section, we test the ability of LUCJ models to dissociate the single and double bonds of H₂/LiH and ethene.

Fig. 9 illustrates the dissociation of H₂ in a minimal basis (2e,2o) comparing RHF and exact values with LUCJ ansatzes with all-to-all, square, and hex qubit connectivities. As is well known, the doubly occupied RHF electron configuration cannot properly dissociate H₂. We find that the exact dissociation curve is recovered when (1) both J_αα and X, (2) only J_αα, and (3) only X are retained in the ansatz. With neither J_αα nor X (for L = 1,2), the potential energy curve of LUCJ is between those of RHF and FCI; the errors near the dissociation limit approach 0.12 Ha, far in excess of our target accuracy of 0.001 Ha. This can be understood by noting that in the same-spin sector of the number–number interaction term, n̂_p,σn̂_p,σ = n̂_p,σ. In other words, this term of the number–number operator reduces to a one-body operator, which can be represented equivalently by the X operator. A similar behavior is found for the dissociation of LiH in a (2e,3o) space, shown in Fig. 10. In this case, while ansatzes with L = 2 yield results of FCI accuracy, ansatzes with L = 1 are inaccurate if either the J_αα or the X term is missing.

Fig. 9 Potential energy curve of H₂ in a (2e,2o) active space of molecular orbitals, using the LUCJ ansatz with all-to-all, square-lattice, and hex-lattice connectivities (top, middle, and bottom), with or without J_αα and X terms. All curves except those without J_αα and X terms and hex-lattice connectivity (green lines, right panels) agree with FCI within 10⁻⁸ hartree.

Fig. 10 Potential energy curve of LiH in a (2e,3o) active space of UMP2 natural orbitals, using the LUCJ ansatz with all-to-all, square-lattice, and hex-lattice connectivities (top, middle, and bottom), with or without J_αα and X terms.

Now we examine the breaking of the carbon–carbon double-bond in ethene (C₂H₄). At equilibrium, ethene is a closed-shell singlet; at dissociation, two methylene (CH₂) molecules each have triplet ground state multiplicity. We limit our calculations to a (4e,4o) space, keeping only the orbitals deriving from the non-bonding molecular orbitals of isolated, triplet methylene – namely 2p_z and sp². Our geometry optimizations enforced planar symmetry (D_2h point group) along the dissociation. More details can be found in Section 2.5.

Data from RHF and a variety of LUCJ ansatzes are compared with those of FCI as shown in Fig. 11 (results from even values of L are shown in the ESI†). Note that the RHF curve is not monotonic (with a “hump” between R = 2.3–2.5 Å) due to our choice to constrain the geometries to a plane. The LUCJ ansatzes corresponding to all-to-all, square, hex, and heavy-hex connectivities can produce exact energies at L > 1, 2, 3, and 4, respectively. Encouragingly, all-to-all, square, and hex LUCJ ansatzes achieve our target accuracy with L ≥ 2; the heavy-hex LUCJ ansatz with L ≥ 4. All of these LUCJ models require quantum circuits that are significantly shallower than qUCCSD.

Fig. 11 Potential energy curve of ethene in a (4e,4o) active space, using the LUCJ ansatz with all-to-all, square-lattice, hex-lattice and heavy-hex lattice (green, blue, orange, and red markers) connectivity and L = 2, 4 (crosses and circles, respectively).

A close inspection of the deviations of the LUCJ energy curves from exact values (Fig. 11, bottom) reveals small but notable cusps in the cases of square and hex ansatzes with L = 2. Albeit relatively mild (on the order of 0.1 mHa), we note that this issue has been previously reported for hardware-efficient ansatzes,⁶⁹ and could complicate calculations of interatomic forces and vibrational properties. It is likely that these non-monotonicities are related to unsuccessful parameter optimization, as discussed further in Section 3.4.

3.3 Hardware perspective

The implementation of quantum algorithms on NISQ hardware faces a number of challenges. The physical time to execute an ansatz circuit must be considerably less than the qubit coherence times, T₁ and T₂, which currently are around 100 μs on state-of-the-art superconducting devices. By contrast, the duration of a two-qubit gate is around 400 ns, so that around 250 layers of two-qubit gates can be executed before qubit decoherence. A second, and arguably more stringent, challenge arises due to various types of errors that result from applying any quantum operation, e.g. a 2 qubit gate. As a simple numerical example, even with 2 qubit gate fidelities of 99%, the application of 11 2 qubit gates will lead to fidelities of <90%.

Localization of the previously proposed UCJ ansatz enables a number of significant advantages from the perspective of quantum resources required to encode the state. Gate counts in terms of the number of spatial orbitals are shown in Table 1. All flavors of the LUCJ approximation do not require SWAP gates. Furthermore, with regard to the e^iĴ part, the circuit depth no longer depends on N. Another key point is that the circuit depth for the square, hex, and linear connectivities is the same – namely, 4 (three from the same-spin part and one from the opposite-spin part). Finally, the required number of number–number gates, U_nn, is reduced from quadratic to linear in N. Device connectivities with fewer neighboring qubits require proportionally fewer U_nn gates.

Table 1 Computational cost of the operations comprising a layer of the Jastrow ansatz, for a system of electrons in N spatial orbitals

Instruction	R z	G	SWAP	U nn	cX	Depth
e^K	2N	N(N − 1)	0	0	0	(1 + N)
e^iJ all-to-all	2N	0			0	(1 + 4N)
e^iJ square	2N	0	0	N + 2(N − 1)	0	4
e^iJ hex	2N	0	0		0	4
e^iJ heavy-hex	2N	0	0		N	8
e^iJ linear	2N	0	0	1 + 2(N − 1)	0	4

---

In the case of a heavy-hex topology, as indicated in Fig. 5, 4 CNOT gates are required to implement an R_zzz gate in the presence of an ancilla. This operation has to be repeated N/4 times since the ratio between ancillas and orbitals, as illustrated in Fig. 2, is 4. Therefore, we have a total of N CNOT gates.

In addition, another salient feature of the LUCJ ansatzes is natural and efficient implementation on superconducting quantum hardware with tunable couplers. Fig. 12 shows the decomposition of the U_nn gates and Givens rotations into sequences of single-qubit gates along with R_xx, R_yy, and R_zz 2 qubit gates (bottom three rows). These R gates are then further decomposed into native 2 qubit gates, R_zx, highlighted in yellow, which can be implemented in ∼400 ns. On such fixed-frequency devices, each number–number gate requires an R_zx native gate, while each Givens rotation requires two native gates. Furthermore, the all-to-all UCJ model requires many U_nn·SWAP gates, each of which requires 3 native gates.

Fig. 12 Implementation of CNOT, U_nn, SWAP·U_nn, and XX + YY gates (rows 1–4 from above) with single- and two-qubit gates. The R_xx, R_yy, and R_zz gates (rows 5–7) along with the CNOT gate can be implemented with single-qubit gates and the R_zx gate (yellow block), native to superconducting devices with fixed-frequency qubits.

In contrast, Fig. 13 shows the decomposition of the 2 qubit gates required in the LUCJ ansatzes onto the native gates on devices with tunable couplers,^49,70 which is denoted as fSim and is represented by the matrix

Fig. 13 Implementation of CNOT, U_nn, U_nn·SWAP, and XX + YY gates (top to bottom) with single-qubit gates and the two-qubit fSim(θ, ϕ) gate (yellow block), native to superconducting devices with tunable couplers.

The U_nn gate is itself a native gate, and can be implemented with 30 ns gate times. This is an order of magnitude faster than is possible on fixed-frequency devices. A Givens rotation can be implemented with a single fSim native gate. When accounting for the 10× speed-up for each native gate, the total gate-time for a Givens rotation is 20× faster when utilizing devices with tunable couplers.

3.4 Parameter optimization

The parameters in the LUCJ ansatzes are optimized variationally. For the majority of the systems investigated in this work, parameter optimization with the BFGS algorithm has not presented insurmountable difficulties, owing to the bootstrapping and parameter initialization protocols developed and employed. However, in Fig. 14 we show a difficult case, namely the optimization of LUCJ parameters for the N₂ molecule near equilibrium (R = 1.2 Å) and near dissociation (R = 1.9 Å). In both cases, the optimization procedure features cusps in the energy and pronounced oscillations both in the norm of the gradient and in the distance between variational parameters in consecutive iterations.

Fig. 14 Optimization of the LUCJ ansatz for N₂ in a (6e,6o) active space of UMP2 natural orbitals at bond lengths R = 1.20 and 1.90 Å (top and bottom), using the LUCJ ansatz with square-lattice connectivity and L = 4. Left, middle, and right panels show the energy, energy gradient, and parameter variations (red, blue, and green, respectively) as the optimization proceeds.

Future work will involve, first and foremost, accelerating and increasing the robustness of the parameter optimization scheme. Second-order solvers which make use of the exact Hessian will be explored. Given the success of classical variational quantum Monte Carlo, a promising step forward will be to optimize the LUCJ parameters via the linear method.^71–73 Another complementary strategy involves further reducing the number of parameters in the LUCJ ansatzes. This could be done, e.g., by discarding all same-spin J_σσ terms, and offloading the corresponding correlation energy contributions via polynomially scaling classical algorithms. This appears plausible given that U_nn-type correlations are predominately from opposite-spin electrons which can occupy the same orbital (same-spin correlation is longer-ranged due to the Pauli principle).

4 Conclusions and outlook

In this work we have presented a new ansatz – the localized unitary cluster Jastrow (LUCJ) wavefunction – which can describe both dynamic and static electron correlations accurately and is amenable to efficient implementation on superconducting quantum devices. While the computational cost of exactly measuring the energy is super-polynomial on classical devices, it is polynomial on quantum devices. We have drawn inspiration from the repulsive Hubbard model, which penalizes double-occupancy of lattice sites, and built this physics into a correlator which acts on a closed-shell, mean-field (RHF) state. The LUCJ ansatz is flexible enough to include the qUCCSD wavefunction, yet retains an avenue (via increasing the number of layers, L) to systematically approach the exact state, even in the presence of multi-reference characters. While we focus on the LUCJ model defined by a square configuration of qubits without SWAP gates – i.e., keeping only on-site opposite-spin and nearest-neighbor same-spin number–number interactions – we also demonstrate that with more sparse arrangements of qubits, e.g. the hex and heavy-hex topologies, one can also converge to very high accuracy with the addition of more layers. In this work we have investigated an antiaromatic compound (square cyclobutadiene) and strongly correlated regimes encountered when stretching single, double, and triple bonds. We are encouraged that compact and hardware-friendly LUCJ ansatzes led to the appropriate dissociation limits, which is challenging to achieve even with sophisticated coupled cluster methods.^74–76 We also made resource estimates which quantify substantial reductions in gates and gate-times, in particular for devices with tunable couplers.

We envision that the LUCJ ansatz will not only be a leading way to perform variational quantum simulations,⁶ but will be an advantageous choice to initialize quantum subspace expansion (QSE), imaginary-time evolution, adiabatic state preparation, and quantum phase estimation calculations in a fault-tolerant era. Developing more robust approaches for parameter optimization will enable us to investigate larger and more complex systems, and we are eager to test whether such a compact ansatz (small values of L) will continue to preserve near-exact accuracy. Future directions include exploring this variationally optimized LUCJ decomposition in the context of time propagation, and as a single-reference ansatz form for the recently proposed non-orthogonal quantum eigensolver approach.³⁸ In light of the LUCJ ansatzes' ability to capture strong correlations with minimal quantum resources, we will explore its use as an active space solver in the contexts of tailored and externally corrected coupled-cluster schemes.^77–81 In this work we have elected to investigate the LUCJ ansatz in the absence of shot noise and decoherence; future work will involve its characterization in more realistic simulations and on quantum hardware.

As a final comment, the research community to date has largely been pursuing electronic structure ansatzes that are either motivated completely by hardware efficiency or by physical intuition. The family of LUCJ ansatzes introduced herein demonstrates that a useful balance between the two priorities is possible, i.e. the dominant type of electron correlation in many strongly correlated systems (involving opposite-spin electrons in a spatial orbital) is faithfully retained without requiring SWAP operations which are a leading bottleneck in terms of reducing circuit depth. It is our hope that this hybrid design principle for quantum chemical ansatzes will inspire future developments in the field.

Data availability

The data that support the findings of this study are available from the corresponding authors upon reasonable request.

Author contributions

J. S. and M. M. conceived the project. M. M., K. S., and J. S. developed the quantum simulation codebase. M. M. carried out numerical simulations. J. S., M. M., M. H. G., and K. B. W. analyzed computational results. All authors participated in the composition and revision of the manuscript.

Conflicts of interest

There are no conflicts of interest to declare.

Acknowledgements

J. S. acknowledges funding from the National Institute of General Medical Sciences of the National Institutes of Health under award number F32GM142231. This research used resources of the National Energy Research Scientific Computing Center (NERSC), a U.S. Department of Energy Office of Science User Facility located at Lawrence Berkeley National Laboratory. M. H.-G., K. B. W., and J. S. acknowledge support from the NSF QCLI program through grant number QMA-2016345. We thank Guglielmo Mazzola for helpful feedback.

References

1. I. M. Georgescu, S. Ashhab and F. Nori, Quantum simulation, Rev. Mod. Phys., 2014, 86, 153.
2. Y. Cao, J. Romero, J. P. Olson, M. Degroote, P. D. Johnson, M. Kieferová, I. D. Kivlichan, T. Menke, B. Peropadre, N. P. Sawaya and et al, ., Quantum chemistry in the age of quantum computing, Chem. Rev., 2019, 119, 10856–10915.
3. B. Bauer, S. Bravyi, M. Motta and G. Kin-Lic Chan, Quantum algorithms for quantum chemistry and quantum materials science, Chem. Rev., 2020, 120, 12685–12717.
4. S. McArdle, S. Endo, A. Aspuru-Guzik, S. C. Benjamin and X. Yuan, Quantum computational chemistry, Rev. Mod. Phys., 2020, 92, 015003.
5. M. Motta and J. E. Rice, Emerging quantum computing algorithms for quantum chemistry, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2021, e1580.
6. M. Cerezo, A. Arrasmith, R. Babbush, S. C. Benjamin, S. Endo, K. Fujii, J. R. Mc-Clean, K. Mitarai, X. Yuan, L. Cincio and et al, ., Variational quantum algorithms, Nat. Rev. Phys., 2021, 1–20.
7. P. J. O'Malley, R. Babbush, I. D. Kivlichan, J. Romero, J. R. McClean, R. Barends, J. Kelly, P. Roushan, A. Tranter, N. Ding and et al, ., Scalable quantum simulation of molecular energies, Phys. Rev. X, 2016, 6, 031007.
8. A. Kandala, A. Mezzacapo, K. Temme, M. Takita, M. Brink, J. M. Chow and J. M. Gambetta, Hardware-efficient variational quantum eigensolver for small molecules and quantum magnets, Nature, 2017, 549, 242–246.
9. F. Arute and et al, ., Hartree-Fock on a superconducting qubit quantum computer, Science, 2020, 369, 1084–1089.
10. J. E. Rice, T. P. Gujarati, M. Motta, T. Y. Takeshita, E. Lee, J. A. Latone and J. M. Garcia, Quantum computation of dominant products in lithium–sulfur batteries, J. Chem. Phys., 2021, 154, 134115.
11. A. Peruzzo, J. McClean, P. Shadbolt, M.-H. Yung, X.-Q. Zhou, P. J. Love, A. Aspuru-Guzik and J. L. O'Brien, A variational eigenvalue solver on a photonic quantum processor, Nat. Commun., 2014, 5, 4213.
12. J. R. McClean, M. E. Kimchi-Schwartz, J. Carter and W. A. de Jong, Hybrid quantum-classical hierarchy for mitigation of decoherence and determination of excited states, Phys. Rev. A, 2017, 95, 042308.
13. J. I. Colless, V. V. Ramasesh, D. Dahlen, M. S. Blok, M. Kimchi-Schwartz, J. McClean, J. Carter, W. De Jong and I. Siddiqi, Computation of Molecular Spectra on a Quantum Processor with an Error-Resilient Algorithm, Phys. Rev. X, 2018, 8, 011021–011028.
14. M. Motta, C. Sun, A. T. Tan, M. J. O'Rourke, E. Ye, A. J. Minnich, F. G. Brandão and G. K.-L. Chan, Determining Eigenstates and Thermal States on a Quantum Computer Using Quantum Imaginary Time Evolution, Nat. Phys., 2020, 16, 205–210.
15. E. Farhi, J. Goldstone, S. Gutmann, J. Lapan, A. Lundgren and D. Preda, A quantum adiabatic evolution algorithm applied to random instances of an NP-complete problem, Science, 2001, 292, 472–475.
16. A. Aspuru-Guzik, A. D. Dutoi, P. J. Love and M. Head-Gordon, Simulated quantum computation of molecular energies, Science, 2005, 309, 1704–1707.
17. A. Anand, P. Schleich, S. Alperin-Lea, P. W. Jensen, S. Sim, M. Díaz-Tinoco, J. S. Kottmann, M. Degroote, A. F. Izmaylov and A. Aspuru-Guzik, A quantum computing view on unitary coupled cluster theory, Chem. Soc. Rev., 2022, 51, 1659–1684.
18. R. J. Bartlett and M. Musiał, Coupled-cluster theory in quantum chemistry, Rev. Mod. Phys., 2007, 79, 291.
19. B. D. Nguyen, G. P. Chen, M. M. Agee, A. M. Burow, M. P. Tang and F. Furche, Divergence of many-body perturbation theory for noncovalent interactions of large molecules, J. Chem. Theory Comput., 2020, 16, 2258–2273.
20. J. B. Robinson and P. J. Knowles, Approximate variational coupled cluster theory, J. Chem. Phys., 2011, 135, 044113.
21. P.-L. Dallaire-Demers, J. Romero, L. Veis, S. Sim and A. Aspuru-Guzik, Low-depth circuit ansatz for preparing correlated fermionic states on a quantum computer, Quantum Sci. Technol., 2019, 4, 045005.
22. M. Motta, E. Ye, J. R. McClean, Z. Li, A. J. Minnich, R. Babbush and G. K. Chan, Low rank representations for quantum simulation of electronic structure, npj Quantum Inf., 2021, 7, 1–7.
23. J. Lee, W. J. Huggins, M. Head-Gordon and K. B. Whaley, Generalized unitary coupled cluster wave functions for quantum computation, J. Chem. Theory Comput., 2018, 15, 311–324.
24. R. Jastrow, Many-Body Problem with Strong Forces, Phys. Rev., 1955, 98, 1479–1484.
25. W. L. McMillan, Ground State of Liquid He⁴, Phys. Rev., 1965, 138, A442–A451.
26. D. Ceperley, Ground state of the fermion one-component plasma: a Monte Carlo study in two and three dimensions, Phys. Rev. B: Condens. Matter Mater. Phys., 1978, 18, 3126.
27. E. Manousakis, The spin-1/2 Heisenberg antiferromagnet on a square lattice and its application to the cuprous oxides, Rev. Mod. Phys., 1991, 63, 1–62.
28. L. Spanu, M. Lugas, F. Becca and S. Sorella, Magnetism and superconductivity in the t-t′-J model, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 77, 024510.
29. M. Casula, C. Attaccalite and S. Sorella, Correlated geminal wave function for molecules: an efficient resonating valence bond approach, J. Chem. Phys., 2004, 121, 7110–7126.
30. A. Zen, E. Coccia, Y. Luo, S. Sorella and L. Guidoni, Static and dynamical correlation in diradical molecules by quantum Monte Carlo using the Jastrow antisymmetrized geminal power ansatz, J. Chem. Theory Comput., 2014, 10, 1048–1061.
31. G. Mazzola, R. Helled and S. Sorella, Phase Diagram of Hydrogen and a Hydrogen-Helium Mixture at Planetary Conditions by Quantum Monte Carlo Simulations, Phys. Rev. Lett., 2018, 120, 025701.
32. C. Genovese, A. Meninno and S. Sorella, Assessing the accuracy of the Jastrow antisymmetrized geminal power in the H₄ model system, J. Chem. Phys., 2019, 150, 084102.
33. E. Neuscamman, Communication: a Jastrow factor coupled cluster theory for weak and strong electron correlation, J. Chem. Phys., 2013, 139, 181101.
34. E. Neuscamman, The Jastrow antisymmetric geminal power in Hilbert space: theory, benchmarking, and application to a novel transition state, J. Chem. Phys., 2013, 139, 194105.
35. G. Mazzola, P. J. Ollitrault, P. K. Barkoutsos and I. Tavernelli, Nonunitary Operations for Ground-State Calculations in Near-Term Quantum Computers, Phys. Rev. Lett., 2019, 123, 130501.
36. Y. Matsuzawa and Y. Kurashige, Jastrow-type decomposition in quantum chemistry for low-depth quantum circuits, J. Chem. Theory Comput., 2020, 16, 944–952.
37. E. Neuscamman, Improved optimization for the cluster Jastrow antisymmetric geminal power and tests on triple-bond dissociations, J. Chem. Theory Comput., 2016, 12, 3149–3159.
38. U. Baek, D. Hait, J. Shee, O. Leimkuhler, W. J. Huggins, T. F. Stetina, M. Head-Gordon and K. B. WhaleyNO to Optimization: A Non-Orthogonal Quantum Eigensolver, arXiv, 2022, preprint, arXiv:2205.09039,  DOI:10.48550/arXiv.2205.09039.
39. J. Shee, M. Loipersberger, D. Hait, J. Lee and M. Head-Gordon, Revealing the nature of electron correlation in transition metal complexes with symmetry breaking and chemical intuition, J. Chem. Phys., 2021, 154, 194109.
40. E. I. Solomon, U. M. Sundaram and T. E. Machonkin, Multicopper oxidases and oxygenases, Chem. Rev., 1996, 96, 2563–2606.
41. E. I. Solomon, F. Tuczek, D. E. Root and C. A. Brown, Spectroscopy of binuclear dioxygen complexes, Chem. Rev., 1994, 94, 827–856.
42. N. F. Mott, The basis of the electron theory of metals, with special reference to the transition metals, Proc. Phys. Soc. A, 1949, 62, 416.
43. N. F. Mott, The transition to the metallic state, Philos. Mag., 1961, 6, 287–309.
44. N. F. Mott, Metal-insulator transition, Rev. Mod. Phys., 1968, 40, 677.
45. P. W. Anderson, The resonating valence bond state in La₂CuO₄ and superconductivity, Science, 1987, 235, 1196–1198.
46. D. Zanchi and H. Schulz, Superconducting instabilities of the non-half-filled Hubbard model in two dimensions, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 9509.
47. C. J. Halboth and W. Metzner, d-wave superconductivity and Pomeranchuk instability in the two-dimensional Hubbard model, Phys. Rev. Lett., 2000, 85, 5162.
48. S. Raghu and S. Kivelson, Superconductivity from repulsive interactions in the two-dimensional electron gas, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 83, 094518.
49. B. Foxen, C. Neill, A. Dunsworth, P. Roushan, B. Chiaro, A. Megrant, J. Kelly, Z. Chen, K. Satzinger, R. Barends and et al, ., Demonstrating a continuous set of two-qubit gates for near-term quantum algorithms, Phys. Rev. Lett., 2020, 125, 120504.
50. F. A. Evangelista, G. K.-L. Chan and G. E. Scuseria, Exact parameterization of fermionic wave functions via unitary coupled cluster theory, J. Chem. Phys., 2019, 151, 244112.
51. Z. Jiang, K. J. Sung, K. Kechedzhi, V. N. Smelyanskiy and S. Boixo, Quantum algorithms to simulate many-body physics of correlated fermions, Phys. Rev. Appl., 2018, 9, 044036.
52. M. Reck, A. Zeilinger, H. J. Bernstein and P. Bertani, Experimental realization of any discrete unitary operator, Phys. Rev. Lett., 1994, 73, 58.
53. W. R. Clements, P. C. Humphreys, B. J. Metcalf, W. S. Kolthammer and I. A. Walmsley, Optimal design for universal multiport interferometers, Optica, 2016, 3, 1460–1465.
54. B. O'Gorman, W. J. Huggins, E. G. Rieffel and K. B. Whaley, Generalized swap networks for near-term quantum computing, arXiv, 2019, preprint, arXiv:1905.05118,  DOI:10.48550/arXiv.1905.05118.
55. M. Caffarel and P. Claverie, Treatment of the Schrödinger equation through a Monte Carlo method based upon the generalized Feynman-Kac formula, J. Stat. Phys., 1986, 43, 797–801.
56. M. Caffarel and P. Claverie, Development of a pure diffusion quantum Monte Carlo method using a full generalized Feynman–Kac formula. I. Formalism, J. Chem. Phys., 1988, 88, 1088–1099.
57. B. Peng and K. Kowalski, Highly efficient and scalable compound decomposition of two-electron integral tensor and its application in coupled cluster calculations, J. Chem. Theory Comput., 2017, 13, 4179–4192.
58. M. Motta, J. Shee, S. Zhang and G. K.-L. Chan, Efficient ab initio auxiliary-field quantum Monte Carlo calculations in Gaussian bases via low-rank tensor decomposition, J. Chem. Theory Comput., 2019, 15, 3510–3521.
59. J. R. Moreno, J. Cohn, D. Sels and M. Motta, Enhancing the Expressivity of Variational Neural, and Hardware-Efficient Quantum States Through Orbital Rotations, arXiv, 2023, preprint, arXiv:2302.11588,  DOI:10.48550/arXiv.2302.11588.
60. G. Aleksandrowicz, T. Alexander, P. Barkoutsos, L. Bello, Y. Ben-Haim, D. Bucher, F. Cabrera-Hernández, J. Carballo-Franquis, A. Chen and C. Chen, et al., Qiskit: an open-source framework for quantum computing, Zenodo, 2019, vol. 16.
61. Q. Sun, T. C. Berkelbach, N. S. Blunt, G. H. Booth, S. Guo, Z. Li, J. Liu, J. D. McClain, E. R. Sayfutyarova, S. Sharma and et al, ., PySCF: the Python-based simulations of chemistry framework, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2018, 8, e1340.
62. Q. Sun, X. Zhang, S. Banerjee, P. Bao, M. Barbry, N. S. Blunt, N. A. Bogdanov, G. H. Booth, J. Chen, Z.-H. Cui and et al, ., Recent developments in the PySCF program package, J. Chem. Phys., 2020, 153, 024109.
63. E. Epifanovsky, A. T. Gilbert, X. Feng, J. Lee, Y. Mao, N. Mardirossian, P. Pokhilko, A. F. White, M. P. Coons, A. L. Dempwolff and et al, ., Software for the frontiers of quantum chemistry: an overview of developments in the Q-Chem 5 package, J. Chem. Phys., 2021, 155, 084801.
64. N. Mardirossian and M. Head-Gordon, ωB97X-V: a 10-parameter, range-separated hybrid, generalized gradient approximation density functional with nonlocal correlation, designed by a survival-of-the-fittest strategy, Phys. Chem. Chem. Phys., 2014, 16, 9904–9924.
65. F. Weigend, F. Furche and R. Ahlrichs, Gaussian basis sets of quadruple zeta valence quality for atoms H–Kr, J. Chem. Phys., 2003, 119, 12753–12762.
66. F. Weigend and R. Ahlrichs, Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: design and assessment of accuracy, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
67. A. Balková and R. J. Bartlett, A multireference coupled-cluster study of the ground state and lowest excited states of cyclobutadiene, J. Chem. Phys., 1994, 101, 8972–8987.
68. P. C. Varras and P. S. Gritzapis, The transition state of the automerization reaction of cyclobutadiene: a theoretical approach using the restricted active space self consistent field method, Chem. Phys. Lett., 2018, 711, 166–172.
69. R. D’Cunha, T. D. Crawford, M. Motta and J. E. Rice, Accuracy assessment of variational quantum computing ansatze across a database of electronic structure problems, J. Phys. Chem. A, 2023, 127(15), 3437–3448.
70. F. Yan, P. Krantz, Y. Sung, M. Kjaergaard, D. L. Campbell, T. P. Orlando, S. Gustavsson and W. D. Oliver, Tunable coupling scheme for implementing high-fidelity two-qubit gates, Phys. Rev. Appl., 2018, 10, 054062.
71. S. Sorella, Wave function optimization in the variational Monte Carlo method, Phys. Rev. B: Condens. Matter Mater. Phys., 2005, 71, 241103.
72. J. Toulouse and C. J. Umrigar, Optimization of quantum Monte Carlo wave functions by energy minimization, J. Chem. Phys., 2007, 126, 084102.
73. M. Motta, G. Bertaina, D. E. Galli and E. Vitali, Implementation of the linear method for the optimization of Jastrow–Feenberg and backflow correlations, Comput. Phys. Commun., 2015, 190, 62–71.
74. I. W. Bulik, T. M. Henderson and G. E. Scuseria, Can single-reference coupled cluster theory describe static correlation?, J. Chem. Theory Comput., 2015, 11, 3171–3179.
75. T. M. Henderson, I. W. Bulik, T. Stein and G. E. Scuseria, Seniority-based coupled cluster theory, J. Chem. Phys., 2014, 141, 244104.
76. L. Zhao, J. Goings, K. Wright, J. Nguyen, J. Kim, S. Johri, K. Shin, W. Kyoung, J. I. Fuks and J.-K. K. Rhee, et al., Orbital-optimized pair-correlated electron simulations on trapped-ion quantum computers, arXiv, 2022, preprint, arXiv:2212.02482,  DOI:10.48550/arXiv.2212.02482.
77. J. Paldus, J. Čížek and M. Takahashi, Approximate account of the connected quadruply excited clusters in the coupled-pair many-electron theory, Phys. Rev. A: At., Mol., Opt. Phys., 1984, 30, 2193.
78. J. Paldus and J. Planelles, Valence bond corrected single reference coupled cluster approach: I. general formalism, Theor. Chim. Acta, 1994, 89, 13–31.
79. O. Hino, T. Kinoshita, G. K.-L. Chan and R. J. Bartlett, Tailored coupled cluster singles and doubles method applied to calculations on molecular structure and harmonic vibrational frequencies of ozone, J. Chem. Phys., 2006, 124, 114311.
80. L. Veis, A. Antalík, J. Brabec, F. Neese, O. Legeza and J. Pittner, Coupled cluster method with single and double excitations tailored by matrix product state wave functions, J. Phys. Chem. Lett., 2016, 7, 4072–4078.
81. M. Mörchen, L. Freitag and M. Reiher, Tailored coupled cluster theory in varying correlation regimes, J. Chem. Phys., 2020, 153, 244113.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d3sc02516k

---