The first step in performing quantum chemistry calculations is to generate the chemical Hamiltonian, which describes the energy operator of the molecular system in terms of the positions of its constituent atoms. In the Chemistry module of QREChem, we generate the Hamiltonian using the self-consistent field (SCF) methods implemented in the PySCF program (Sun et al., 2018; Sun et al., 2020).

The Hamiltonian matrix elements, which represent the contributions of the various terms in the Hamiltonian, are defined in terms of one-electron integrals, h _pq , and two-electron integrals, h _pqrs (Low et al., 2019). These integrals depend on the molecular properties such as the atomic coordinates, the charge, and the choice of basis set. SCF calculations, such as the restricted Hartree-Fock (RHF) method, are performed to obtain the one- and two-electron integrals for the chosen molecular system. QREChem requires the definition of the chemical system of interest, in terms of the atomic coordinates, the charge, and the target basis set.

It is important to note that the SCF calculations can be time-consuming, especially for larger numbers of atoms and larger basis sets, and can be tricky to properly converge. However, QREChem supports the standard fcidump file format, which stores the one- and two-electron integrals, allowing a user to generate these integrals using a different program and then interface with the other modules of QREChem. The fcidump file format is a widely used format for storing quantum chemistry Hamiltonians (Knowles and Handy, 1989) and can be produced by other quantum chemistry packages, such as Gaussian (Foresman and Frish, 1996), MolPro (Werner et al., 2012), or Psi4 (Smith et al., 2020). Once the one- and two-electron integrals are obtained, the data is available to other modules within QREChem.

There are many proposed quantum algorithms for solving for the ground state energy in chemical problems, including quantum phase estimation (QPE) (Kitaev, 1997; Abrams and Lloyd, 1999), the variational quantum eigensolver (VQE) (Peruzzo et al., 2014), combinations of the two (Otten et al., 2022; D’Cunha et al., 2023), and quantum machine learning methods (Xia and Kais, 2018). Within each of these families of algorithms, there are a substantial number of possible variations. In this work, we focus on QPE using Trotterization (Ortiz et al., 2001; Babbush et al., 2015). Here, we provide a brief overview of QPE using Trotterization.

QPE solves for the eigenvalue, λ _k , for an eigenvector |v _k ⟩ of some unitary matrix, U. Aside from its use in ground state energy estimation in quantum chemistry, it also finds use in Shor’s prime number factoring algorithm (Shor, 1999) and the Hassidim–Harrow–Lloyd algorithm for matrix inversion (Harrow et al., 2009). Given a Hamiltonian generated by the Chemistry module, the unitary matrix, U is can be written as: U | v k 〉 = e i H ̂ τ | v k 〉 = e i 2 π ϕ | v k 〉 , (1) where τ is a scale factor to map the eigenvalues of Hτ onto (0, 2π] or (-π, π]. Assuming |v _k ⟩ is the ground state, the ground state energy is then mapped to the phase acquired, E = 2πϕ/τ, where units have been chosen such that ℏ = 1. The computational complexity of QPE is dependent on how the unitary matrix of Eq. 1 is implemented. In QREChem, we focus on Trotterization (Ortiz et al., 2001; Babbush et al., 2015), but many other strategies, such as Taylorization (Berry et al., 2015) and qubitization (Low and Chuang, 2019) have been proposed. The Trotterized version of the propagator, U, is

By choosing some finite number of Trotter steps, n, U is only represented approximately. A first-order Trotter formula truncates Eq. (2) at some number of steps. Choosing a sufficient number of Trotter steps for a given accuracy is important for obtaining an accurate estimate of the total resources. Higher-order Trotter-Suzuki formulas (Suzuki, 1993) can be used to decrease the number of steps at the cost of increasing the complexity of each step. The standard fermionic quantum chemistry Hamiltonian has O(N ⁴) terms. A fermion-to-spin mapping is required to implement each fermionic term on qubits. Using the Jordan–Wigner (Jordan and Wigner, 1993) transformation introduces an O(N) overhead, leading to a total complexity of O(N ⁵) for the Trotterized evolution. The evolution of the phase is mapped to an ancilla register (introducing limits on the precision, based upon the number of ancilla) and, using the quantum Fourier transform (QFT) (Shor, 1994), the ground state energy can be read out. In realistic settings, the true eigenstate |v _k ⟩ is unknown and an approximation must be used. This introduces an additional overhead in the success probability which scales as the overlap of the approximate state, |ϕ⟩, with the true eigenstate, i.e., |⟨ϕ|v _k ⟩|².

Within QREChem, we provide resource estimates for Trotterized QPE by first estimating the resources required for a single Trotter time step and then estimating the total number of Trotter steps. The resource estimation tools within Microsoft Quantum Development Kit (QDK) (Svore et al., 2018; Low et al., 2019) are able to efficiently provide such logical gate estimates, even for very large systems. We estimate the total number of Trotter steps necessary as n o 3 / 2 , where n _o is the number of orbitals used in the Hamiltonian which is based on heuristic estimates (Poulin et al., 2015; Reiher et al., 2017). Beyond number of Trotter steps, we also need to accurately compute the number of ancilla necessary to reach the a user-defined desired precision ϵ _p (which, by default, we take to be 1 mHa). This allows us to calculate the base number of binary digits necessary, n _b = −log₂(ϵ _p /ΔE _R ), where ΔE _R is a scaling factor that estimates the spectral range (i.e., ΔE _R ≈ E _max − E _min) (D’Cunha et al., 2023) and is taken to be 1 Ha (a number chosen heuristically to cover all studied molecules). The phase resulting from a QPE is given as a binary fraction and the number of bits in this fraction (the precision) is determined by the number of ancilla qubits used in the QPE algorithm. If the eigenvalue is not exactly representable with n _b bits of precision, the returned value will, instead, be mapped into the finite precision of n _b bits, causing a chance for error. To increase the QPE success probability, additional ancilla can be used. If the eigenvector is known precisely, the total number of ancilla, n _a , is a function of the desired failure probability, p _f (Nielsen and Chuang, 2000), n a = n b + log 2 2 + 1 2 p f . (3)

It is very unlikely to know the true eigenstate a priori. More accurate formulas can be derived which take into account errors in the Trotterization, ϵ _t , as well as the true gap ΔE = E ₁ − E ₀ (Li, 2022) n a = n b + log 2 2 + ϵ t 2 2 p f Δ E 2 . (4)

Since the true gap, ΔE, is generally unknown, we instead use equation-of-motion (EOM) coupled-cluster with singles and doubles (CCSD) as implemented in PySCF (Sun et al., 2018) to estimate the gap. In cases where CCSD becomes too expensive, other methods with tunable cost and accuracy, such as selected configuration interaction (Chien et al., 2018) (which can be as cheap as Hartree-Fock) can be used. We also use a target Trotter error of chemical accuracy (ϵ _t = 1.6 mHa), rather than an observed Trotter error.

To calculate the total number of rotation gates, CNOT gates, and the total depth, we combine the estimates for a single Trotter step (using Q# and the Microsoft QDK) with the estimate of the total number of Trotter steps ( n o 3 / 2 ) and multiply that by the number of ancilla, n _a , as each ancilla will require evolution to some long time, giving, for example, the total number of rotation gates n r , t o t = n r n o 2 / 3 n a , (5) where n _r is the number of rotation gates for a single Trotter step. Similar equations are used for the total depth and total number of CNOTs.

Using the same Hamiltonians generated in the Chemistry module, we utilized TFermion (Casares et al., 2022) and OpenFermion (McClean et al., 2020) to provide comparison logical resource estimates. TFermion provides estimates of a variety of quantum algorithms, including variants of Trotterization (Campbell, 2019) and Taylorization (Babbush et al., 2016), among others. It uses properties of the computed Hamiltonians, such as the 1-norm, combined with analytic formulas derived from the literature. OpenFermion provides estimation of more advanced algorithms, such as qubitization with low rank factorization (Berry et al., 2019).

Future improvements to the Algorithm module would involve the implementation of resource estimates for other evolution algorithms with explicit circuit constructions, such as qubitization and Taylorization, as well as other algorithms, such as VQE. Furthermore, overheads relating to the preparation of sufficient overlap initial states will also be implemented.

Quantum error correction (QEC) is an essential feature of any viable quantum computing system due to the intrinsic susceptibility of quantum systems to errors. These errors can be caused by a variety of factors, including decoherence and operational imperfections. In simple terms, QEC codes encode a logical qubit into several physical qubits, and through the process of measurement and classical post-processing, corrects the inevitable errors, extending the effective processing time before errors destroy the quantum state.

There are many proposed codes within QEC. Among them, the surface code (Fowler et al., 2012) stands out due to its high error threshold, relative ease of implementation, and planar geometry, which matches many proposed quantum architectures. Like any QEC code, the implementation of surface code necessitates a large number of physical qubits to encode a single logical qubit and presents a significant overhead in terms of resource requirements. Various implementations of the surface code have a given distance, d, which refers to the minimum number of physical qubits that must be affected by errors in order to cause a logical error. Within QREChem, we provide estimates of the QEC space and time overheads of the surface code using OpenFermion (McClean et al., 2020). The surface code cost estimator takes in the number of Toffoli gates, the number of logical qubits, a physical error rate, and the estimated surface code cycle time. It then estimates the total error, including error contributions from both magic state distillation (which is necessary to produce the Toffoli gates) (Litinski, 2019a) and due to the physical error rate. The physical error rate, p _P , logical error rate, p _L , and distance are approximately related via (Fowler and Gidney, 2018; Webber et al., 2022) p L = 0.1 100 p P d + 1 / 2 . (6)

The total error and space-time (number of qubits times number of seconds) is estimated for a various distances d of the surface code. The best, in terms of space-time, distance d estimate which has total error below a threshold ϵ _sc (which we by default take to be 0.1) is returned as the optimal surface code.

QREChem allows a user to input the desired total algorithmic success probability. Along with the total depth, which is estimated in the Algorithm module, and the physical error rate and cycle time, which is provided by the Hardware module, the error correction module provides the best surface code distance d. This provides an initial estimate of the QEC overhead in terms of the number of physical qubits and total runtime. Future module development will include more accurate estimates of the overhead in numbers of gates needed for the surface code (via, for example, the methods in Ref (Litinski, 2019b).), as well as estimates for other QEC codes (such as color codes (Litinski et al., 2017)).

The underlying quantum computing architecture plays a pivotal role in the overall resource estimates. Various quantum hardware can have vastly different error rates, gate times, connectivities, native gate sets, etc (Menickelly et al., 2023). The Hardware module captures several hardware-specific factors that can significantly affect the performance and resources of quantum algorithms. Within the initial release of QREChem, the Hardware module consists of a high-level description of the underlying hardware, including gate times and physical error rates. With this simplistic Hardware module, we estimate the logical depth by assuming that all single qubit gates can be performed in parallel batches using n _o , the number of orbitals, qubits, and the CNOT gates cannot be performed in parallel ( d = n r n o + n c , where n _r is the number of rotation gates and n _c is the number of CNOTs). To estimate the number of T gates n _t , which are likely to be the most expensive gate for fault-tolerant, error-corrected quantum computers (Bravyi and König, 2013; Beverland et al., 2016), we use estimates of the circuit synthesis cost of arbitrary rotations from Clifford + T gates (Selinger, 2015) n t = n r 10 + 12 ⁡ log 2 ϵ s s − 1 , (7) where ϵ _ss is the synthesis error, which we take to be 10^–9. We choose this value of the threshold to keep the synthesis error well below the standard 1 / N g bound, where N _g is the number of gates (Hastings, 2016), for all circuits studied. This is a tunable parameter which can be varied by a user. Total runtime is computed by interfacing with the Error Correction module, which requires a surface code cycle time, physical error rate, number of Toffoli gates (which is related to the number of T gate) and number of logical qubits. In the initial release of QREChem, we use experimentally demonstrated cycle times to abstract away the hardware details.

Future development of this module will incorporate underlying connectivity, specific noise models, and compilation to native gate sets, along with gate times, to computer the surface code cycle time, as it is evident from the results of this work that a fault-tolerant quantum computer will be required to execute QPE quantum circuits.