Quantum chemistry on quantum computers: quantum simulations of the time evolution of wave functions under the S2 operator and determination of the spin quantum number S

Quantum computers have an enormous impact on quantum chemical calculations. Approaches to calculate the energies of atoms and molecules on quantum computers by utilizing quantum phase estimation (QPE) and the variational quantum eigensolver (VQE) have been well documented, and dozens of methodological improvements to decrease computational costs and to mitigate errors have been reported until recently. However, the possible methodological implementation of observables on quantum computers such as calculating the spin quantum numbers of arbitrary wave functions, which is a crucial issue in quantum chemistry, has been discussed less. Here, we propose a quantum circuit to simulate the time evolution of wave functions under an S operator, exp( iSt)|Ci, and integrate it into the QPE circuit enabling us to determine the spin quantum number of the arbitrary wave functions. We demonstrate that the spin quantum numbers of up to three spins can be determined by only one qubit measurement in QPE.


Introduction
Quantum computing is the focus of current interest in modern society. Although currently available quantum computing devices are noisy and of intermediate scale, quantum computers consisting of more than 100 qubits are expected to appear soon. This reminds us that we are close to ''quantum supremacy'', 1 in which quantum computing devices solve problems that are practically intractable for classical computers. Among the diverse topics in quantum computing, the electronic structure calculations of atoms and molecules are one of the most studied realms. 2-28 Sophisticated quantum chemical calculations on quantum computers can bring a paradigm shift in chemistry and related fields such as nitrogen fixation processes, drug designs, elucidation of catalytic mechanisms of enzymes, and so on.
The ideas to use computers built of quantum mechanical elements to simulate quantum mechanics were proposed by Feynman in the early 1980s. 29 An approach to calculate full-CI energy was proposed by Aspuru-Guzik and coworkers in 2005 2 by using quantum phase estimation (QPE). 30 In the QPE-based full-CI calculations, the time evolutions of the wave function |Ci using a Hamiltonian are simulated and an energy eigenvalue E is extracted as a phase f, as given in eqn (1), using inverse quantum Fourier transformation. 31 exp(ÀiHt)|Ci = exp(ÀiEt)|Ci = exp(Ài2pf)|Ci (1) In 2014, a quantum-classical hybrid system known as a variational quantum eigensolver (VQE) was proposed to solve the unitary coupled cluster (UCC) ansatz. 19,20 The VQE has attracted attention as the computation model in noisy intermediatescale quantum (NISQ) 32 devices, because coherence times required for the simulations are much lower than the QPEbased approaches. [19][20][21][22][23][24][25][26] Despite the rapid progress in the methodologies to calculate molecular energies on quantum computers, approaches to calculate physical quantities other than energies have been less documented. In fact, to our knowledge, there is no theoretical proposal to determine spin quantum numbers of arbitrary wave functions on quantum computers. The determination of spin quantum numbers is an important task, especially for the study of open shell systems. Open shell molecules have many lowlying electronic states with different spin multiplicities and characterization of the spin properties of their electronic ground state is one of the most fundamental issues in quantum chemistry. Molecules and materials having open shell electronic structures play an important role in chemistry and related fields: molecules undergoing covalent bond cleavages, molecular magnets and molecular spin devices, transition metal complexes in active centers of enzymes, etc., and thus the study of open shell electronic structures is crucial in current chemistry and condensed matter physics. The current computational quantum chemistry is based on many approximations and the calculations often fail to predict chemical properties for larger systems. For systems with few electrons/spins, the Lancoz methods or exact diagonalization technique for many body interactions can afford more accurate and reliable energies, 33 but the application of such an approach on classical computers is difficult for larger systems due to the exponential scaling of the computational time against the system size. Noticeably, molecular spin quantum technologies and quantum-assisted biomolecular or organic molecular modellings are the topics of ongoing interest in this field. For example, electron spin-mediated nuclear quantum simulations with switchable interaction were studied by using a vanadyl-based two-qubit system, 34 adiabatic quantum computations using molecular spins were implemented to execute prime factorization, 35 and three-electron spin qubits were exemplified to provide a realistic path forward to solid state quantum information processing. 36 Here, we propose a quantum circuit that allows us to efficiently simulate the time evolution of wave functions under the S 2 operator to calculate the S 2 eigenvalues of arbitrary wave functions by using QPE.

Theory
The spin quantum number S of wave functions can be determined using the eigenvalue of the S 2 operator as defined in eqn (2)-(4). 37 S 2 |C(S)i = S(S + 1)|C(S)i (2) Here, the summation in eqn (3) runs over unpaired electrons.
In the quantum chemical calculations on classical computers, equations to calculate the S 2 expectation values for spinunrestricted Hartree-Fock (UHF) and electron correlated wave functions have been reported. 38,39 On the other hand, in the quantum chemical calculations on quantum computers, analogous to calculating full-CI energies using QPE, we can readout the S 2 eigenvalue as the phase j by applying QPE if we can simulate the time evolution of wave functions under the S 2 operator as given in eqn (5) on quantum computers. exp(ÀiS 2 t)|Ci = exp[Ài{S(S + 1)}t]|Ci = exp(Ài2pj)|Ci (5) It should be emphasized that the S 2 operator can be used not only for the calculation of the spin quantum numbers, but also it is applicable to accelerating quantum simulations and to implementing new algorithms. For example, we can speed up adiabatic state preparations 2 of spin-singlet states by using the S 2 operator as a penalty function, H(t) = (1 À t)H HF + tH full-CI + S 2 (0 r t r 1). The S 2 operator shifts the energy eigenvalue for spin multiplet states depending on the spin quantum number, and therefore the Hamiltonian with the S 2 operator has inherently larger energy gaps between singlet and spin-multiplet states. Importantly, the addition of the S 2 operator does not change the structure of eigenfunctions because the S 2 operator commutes with H. It should be noted that the idea to use the S 2 operator as the cost function in the quantum simulations was proposed to avoid kinks in the computations of potential energy surfaces using the VQE. 21,22 Spin annihilations based on the Löwdin projection operator 40 and adiabatic spin purifications by using a time-dependent Hamiltonian H(t) = H + tS 2 are promising applications relevant to the S 2 operator. To simulate the time evolution in eqn (5) on quantum computers, we have to define the S 2 operator in the second quantized form and transform it to strings of Pauli operators by using the Jordan-Wigner transformation (JWT) or other wave function mapping techniques. 2, [41][42][43] In JWT, each qubit represents the occupation number of a particular spin orbital: |1i if the spin orbital is occupied, otherwise |0i. The first and second terms in eqn (4) can be rewritten as in eqn (6) and (7), respectively, by using creation (a † ) and annihilation (a) operators.
a y pa a pa a y qa a qa þ a y pb a pb a y qb a qb À 1 4 a y pa a pa a y qb a qb þ a y pb a pb a y qa a qa (6) Therefore, S Z ( p)S Z (q) and [S + ( p)S À (q) + S À ( p)S + (q)] terms are expressed as the product of number operators n p (n p = a † p a p ) and two-electron excitation operators, respectively, in the second quantized form. By using the JWT and quantum circuit designs proposed by Whitfield and coworkers, 5 the quantum circuits to simulate the time evolution in eqn (5) can be constructed. For the H 2 molecule with the minimal basis set (two molecular orbitals), the quantum circuit to simulate eqn (5) can be composed of 76 one-qubit and 54 two-qubit gates in the upper bound (72 one-and 48 two-qubit gates for the [S + ( p)S À (q) + S À ( p)S + (q)] term, 4 one-and 2 two-qubit gates for the S( p) 2 term, and 4 two-qubit gates for the S Z ( p)S Z (q) term), without Trotter decompositions. For general molecules, the number of quantum gates scales to N MO 2 , where N MO is the number of molecular orbitals. The reason why more than 100 quantum gates are required for the quantum simulations of exp(ÀiS 2 t)|Ci in two molecular orbital systems is the fact that the [S + (p)S À (q) + S À (p)S + (q)] term is expressed by two-electron excitation operators in the second quantized form. However, as discussed below, the number of required quantum gates can be drastically reduced by introducing a new wave function mapping method for the spin operators S.
In our previous publications, we introduced a spin coordinate mapping (SCM) to construct spin symmetryadapted configuration state functions (CSFs) on quantum computers. 15,16 In the SCM scheme, the qubit represents the spin degree of freedom of an unpaired electron in the singly occupied molecular orbital (SOMO), namely |0i for |ai and |1i for |bi. The SCM is convenient for the manipulations of spin degree of freedom, e.g., the [S À (p) + S + (p)] operation becomes a NOT operation, although it corresponds to one-electron excitations in JWT. However, the originally proposed SCM is available only for SOMOs, and here we extend the SCM for general occupations and propose generalized spin coordinate mapping (GSCM).
The correspondence between JWT and GSCM is given in Table 1. As discussed above, in the JWT, qubits store the occupation number of spin orbitals. In Table 1, the first and second qubits correspond to a and b spin orbitals, respectively. In the GSCM, the first qubit specifies whether the molecular orbital is singly occupied or not, and the second qubit represents the occupation number of the b spin orbital. It should be noted that the GSCM has a similarity to the Bravyi-Kitaev transformation (BKT) for wave function mapping in which information on both occupation and parity is stored nonlocally. 41 Interconversions between JWT and GSCM are achieved by applying controlled-NOT (CNOT) gates; applying a NOT operation to the first qubit if the second qubit is the |1i state.
To simulate exp(ÀiS 2 t)|Ci on quantum computers using the GSCM, we rewrite the S 2 operator using permutation operators P pq and number operators of unpaired electrons n spin p = |n pa À n pb |, where n pa and n pb are number operators acting on the aand b-spin orbitals, respectively, of the p-th molecular orbital: 37 Here, P pq interchanges spin coordinates of unpaired electrons in the p-th and q-th molecular orbitals, as in eqn (9).
It is clear that n spin p can be simulated by the one-qubit operation in the GSCM, although it becomes two-qubit operations in the JWT. The time evolution operator arising from P pq operators, namely exp(ÀiP pq t), seems to be non-trivial at a glance, but it can be expanded as in eqn (10) by using the relationship (P pq ) 2 = 1 (identity operator). exp(ÀiP pq t) = cos(t) À iP pq sin(t) (10) A quantum circuit for the simulations of eqn (5) for two molecular orbital systems is illustrated in Fig. 1 Here, we introduced one ancillary qubit to simulate the P pq terms. Clearly, the number of quantum gates required for the simulations of exp(ÀiS 2 t)|Ci is dramatically reduced by using the GSCM (3, 7, and 3 for one-, two-, and three-qubit operations, respectively). It should also be noted that Trotter decomposition is required only for the P pq terms, because [n spin p , n spin q ] = 0 and [n spin p , P pq ] = 0 but [P pq , P qr ] a 0.

Results and discussion
The results of quantum circuit simulations of the time evolution exp(ÀiS 2 t)|Ci for the two-and three-spin systems starting from |Ci = |abi and |aabi are depicted in Fig. 2 and 3, respectively, and several simulations starting from other initial  states are given in the ESI. † The simulation programs were developed by using OpenFermion 44 and Cirq 45 packages. The simulations were carried out with t = 2p and 360 Trotter slices (t = 2p/360 for a single Trotter step) at the first order Trotterizations. The initial wave functions we used (|abi and |aabi) are not eigenfunctions of the S 2 operator but mixtures of two spin eigenfunctions, as given in eqn (11) and (12).
The spin eigenfunctions with different spin quantum numbers evolve in a different manner, which causes interferences to lead state transformations from |abi to |bai, and |aabi to |abai and |baai. The calculated overlaps between the wave functions obtained from the quantum circuit simulations and the wave functions under exact time evolutions are larger than 0.9999996 everywhere.
By integrating the quantum circuit for exp(ÀiS 2 t)|Ci based on the GSCM into the QPE circuit, as illustrated in Fig. 4, we can determine the spin quantum number of the wave function |Ci. For two-spin systems, the spin-triplet wave function evolves as exp(ÀiS 2 t)|Ci = exp(À2it)|Ci and therefore the measurement of the qubit for QPE always gives |1i after time evolution t = p/2, although the spin-singlet wave function does not evolve by S 2 operators and the measurement always gives |0i, because the S 2 eigenvalue of spin-singlet states is 0. For three spin systems, the spin states S = 3/2 and S = 1/2 can be discriminated deterministically, by applying time evolution exp(ÀiS 2 t)|Ci with t = p/3 and following phase shift operation Z Z with Z = p/4 to the qubit to be measured. Detailed discussions are given in the ESI. † For systems of more than three spins, measurements of two or more qubits are required to determine their spin quantum number, but estimations of the number of qubits and optimal time for evolutions are trivial, because the possible spin state is limited by the number of electrons and the eigenvalue of the S 2 operator is S(S + 1). The QPE simulation results for two-spin systems are given in Fig. 5. In the QPE simulations in Fig. 5 we gradually changed evolution time t from 0 to p, and plotted the probability to obtain the |1i state by using the spin-singlet wave function |C(S = 0)i = (|abi À |bai)/O2 and the spin-triplet wave function |C(S = 1)i = (|abi + |bai)/O2. The probability was calculated by repetitively performing the QPE quantum circuit simulations 10 000 times for each evolution time and counting how many times the |1i state was obtained. In the case of the spin-singlet state, the QPE always returns to the |0i state. By contrast, for the spin-triplet state, the probability to measure the |1i state depends on the evolution time, and after p/2 of time evolution, the probability becomes maximum as expected.   . 4 A quantum circuit for one-qubit quantum phase estimation. For two-and three-electron systems, spin quantum numbers are calculated in a deterministic manner by setting t = p/2 and Z = 0, and t = p/3 and Z = p/4, respectively (see text for details).  We also performed quantum circuit simulations of QPE using a broken-symmetry wave function |C(BS)i = |abi that corresponds to the 1 : 1 mixture of spin singlet and triplet wave functions (see eqn (11)). By running quantum circuit simulations 100 000 times for t = p/2, we obtained the |1i state 49 899 times, which is very close to the ideal probability 0.5. Similar results were obtained for the quantum circuit simulations of three-electron systems, as discussed in the ESI. †

Conclusions
In conclusion, we have developed a quantum circuit for the time evolution of wave functions under the S 2 operator, and we use it to determine the spin quantum number of wave functions by means of QPE, illustrating that the quantum circuit approach underlain by quantum algorithms affords extreme efficiency in evaluating observables, which is a seemingly intractable problem. By using the GSCM implemented for the manipulations of wave functions with spin operators, the number of quantum gates required for the simulation is drastically reduced compared with that using the JWT. The spin quantum number of wave functions can be determined from only one qubit readout in QPE for two-and three-spin systems. These results show that rapid computations of the spin quantum numbers S's on quantum computers are possible, by applying the proposed quantum circuits. The approach developed here opens the door toward theoretical characterizations of the spin properties of open shell molecules in their electronic ground state, which is one of the most fundamental issues in open shell chemistry, by using quantum computers. Again, we emphasize that the S 2 operator can be used for more than just the calculations of the spin quantum number. Applications of the S 2 operator for other quantum chemical problems on quantum computers are underway and will be discussed in a forthcoming paper.

Conflicts of interest
There are no conflicts to declare.