Molecular Dynamics with Nuclear Effects on Quantum Computers

AI summary

Overview

Research area: Quantum computing for quantum chemistry — specifically ab-initio molecular dynamics (AIMD) with nuclear quantum effects treated via the nuclear-electronic orbital (NEO) framework.

Technical level: Advanced. The paper assumes familiarity with second quantization, Hartree-Fock, active spaces, VQE ansätze, and vibrational spectroscopy.

Scope (one sentence): The authors propose and numerically validate a hybrid quantum-classical algorithm that computes NEO ground-state energies and nuclear forces on a simulated quantum computer via VQE while propagating molecular geometries classically, demonstrated on H₂, H₂O and the Zundel cation H₅O₂⁺.

What This Paper Is About

Nuclear quantum effects (proton tunneling, zero-point motion, hydrogen bonding) strongly influence proton transfer and biological chemistry, but including them in classical quantum-chemistry calculations blows up the full-configuration-interaction space exponentially — exact NEO calculations are infeasible even for tiny molecules. The paper's goal is to show that a quantum computer can act as the energy-and-force engine inside a molecular dynamics loop, keeping qubit requirements linear in system size, and to check whether the resulting simulated vibrational spectra match experiment.

Key Contributions

1. A hybrid quantum-classical AIMD algorithm with nuclear quantum effects. The method extends the Fedorov-style quantum Born-Oppenheimer dynamics to the NEO framework: at each timestep a VQE run on the quantum device produces the instantaneous ground state and energy, Hellmann-Feynman forces are measured as expectation values of ∇_I H, and the geometry is updated classically with the velocity Verlet integrator.

2. A NEO active-space embedding with a matching gradient formulation. The inactive/active partitioning is generalized to include nuclear orbitals, and the inactive Fock operators are written in a more general (non-symmetry-exploiting) form so that differentiating them correctly captures the molecular-orbital response of both active and inactive orbitals to nuclear motion.

3. A practical active-space selection protocol. Because standard automated selectors (e.g. autoCAS) produce spaces too large for NISQ devices, the authors combine single-orbital entropies, orbital "variability" along the expected reaction path (accumulated 1 − overlap over trajectory images), and a symmetry requirement that mirror-image orbitals be included together. They also introduce a reference-geometry overlap/linear-assignment scheme ("corresponding orbitals") so orbitals do not swap in and out of the active space as geometries change.

4. A hardware-motivated qubit-reduction and benchmarking study. Parity mapping with particle-number tapering removes two qubits, with two further qubits tapered for NEO (nuclear number conservation) and the nuclear beta-spin sector projected out. Compact hardware-efficient ansätze (real amplitudes, tUPS) are benchmarked against UCCSD/NEO-UCCSD.

Main Findings

• Hydrogen molecule (Table 1, relative errors in %): UCCSD gives 4.6% (purely electronic) and −5.0% (NEO); real amplitudes 10.8% and 0.6%; tUPS 5.4% and 0.3%; harmonic Hessian at the FCI level 6.7% and −4.6%; harmonic Hessian at the DFT level (PBE0/6-31G) 1.8%. NEO tends to underestimate frequencies while the purely electronic treatment overestimates them.

• Hardware-efficient ansätze perform well under NEO. For H₂ the NEO errors of the real amplitudes and tUPS ansätze are below one percent, which the authors attribute to a "fortunate cancelling of errors" — the general overestimation by simpler ansätze offsetting the NEO tendency to underestimate.

• tUPS beats real amplitudes in the electronic case by roughly a factor of two for H₂ (5.4% vs 10.8% purely electronic; 0.3% vs 0.6% with NEO), while remaining much shallower than UCCSD.

• Water (Table 2, relative errors in %): NEO improves on the purely electronic results for all three modes. UCCSD: mode 1 6.8 → −0.9, mode 2 5.3 → −1.0, mode 3 2.8 → −1.9. Real amplitudes: 7.6 → 0.0, 8.7 → 5.8, 10.6 → 5.9. tUPS: 8.0 → 0.1, 8.4 → 4.8, 9.7 → 5.8. Harmonic Hessian at the CASCI level: 7.1 → 0.3, 1.8 → −3.4, 4.2 → −0.7. The DFT-level Hessian errors (PBE0/6-31G) are listed as −2.9, −4.0 and −2.3 percent.

• Water error pattern is not fully explained. The authors note that the NEO error is smaller for water than for hydrogen while the purely electronic error is worse for modes 1 and 2, and offer three untested hypotheses: that direct quantum-nucleus–quantum-nucleus interactions are hard to simulate, that high-frequency modes are hard to simulate, and that the small basis sets and active spaces are responsible.

• Zundel ion: nuclear quantum effects make proton shuttling effectively barrierless. The abstract reports that upon inclusion of nuclear quantum effects the proton shuttling motion in H₅O₂⁺ becomes effectively barrierless and that errors in the simulated frequencies improve significantly. The section reporting the corresponding numbers is cut off in the supplied text, so those numerical values are not available here.

• Larger nuclear basis sets degrade the dynamics. As shown in the paper's Appendix G (referenced, content not included), simulations with larger nuclear basis sets severely underestimate frequencies, because additional angular momentum functions let quantum nuclei move relative to the basis-set centers, weakening the force response the algorithm relies on.

Methodology in Plain English

The idea is a loop: the quantum computer solves the electronic-plus-nuclear electronic-structure problem for a fixed arrangement of atoms, and a classical computer then moves the atoms a small step according to the forces it derives.

• Treating protons as quantum particles. Instead of assuming nuclei are classical points, the NEO approach places selected light nuclei (here, protons) in their own orbitals and solves coupled nuclear-electronic self-consistent-field equations. The Hamiltonian then contains electronic one- and two-body terms, nuclear one- and two-body terms, and an electron-nucleus attraction term, with the proton mass set to 1836 a.u. and its charge to −1 a.u.

• Forces from the Hellmann-Feynman theorem. Because the ground state is variationally optimized at every step, the force on a nucleus is simply the expectation value of the derivative of the Hamiltonian with respect to that nuclear coordinate, plus the classical nuclear repulsion gradient.

• Keeping the circuits small. An active-space transformation freezes the always-occupied orbitals at the Hartree-Fock level and keeps only a handful of active orbitals for the quantum device. Parity mapping with particle-number tapering removes redundant qubits; for NEO, two more qubits are tapered and the nuclear beta-spin block is projected out.

• Ground states by VQE. A parameterized circuit is optimized by a classical optimizer (COBYLA) to minimize the energy expectation value, then the same optimized circuit is used to measure force operators.

• Nuclear basis centers move with the nuclei. The semi-classical approximation treats nuclear basis-set centers as the mean positions of the quantum nuclei, requiring a one-to-one mapping between quantum nuclei and basis centers and excluding higher-angular-momentum nuclear functions.

• Validation by spectroscopy. For each molecule, geometries were optimized on (NEO)FCI or (NEO)CASCI surfaces with the FIRE2 optimizer, Hessians computed, and selected normal modes perturbed by 50 meV. Trajectories of 1 ps with a 121 as timestep (33 cm⁻¹ frequency resolution) were run, and velocity power spectra and IR spectra were compared with experimental frequencies; peak positions were extracted by Gaussian fits.

• Computational stack. Electronic 6-31G basis; two nuclear basis sets (DZSNB and DZSPNB) derived from DZSPDN; active spaces of 5 orbitals/6 electrons for water and 5 orbitals/4 electrons for the Zundel ion (H₂ used all orbitals); Qiskit 1.3.2 statevector simulator, PySCF 2.7.0, NumPy 1.26.4, and the Atomic Simulation Environment Calculator interface.

Why This Matters

The work shows that a NISQ-era device, in principle, can serve as the energy-and-force engine for dynamics that include nuclear quantum effects — a regime where classical exact treatment scales exponentially. It also shows that cheap hardware-efficient circuits can rival chemically motivated UCCSD circuits once NEO is used, which matters for near-term hardware feasibility.

Real-world applications:

• Proton-transfer and proton-coupled electron transfer chemistry, including enzyme catalysis and bioenergetics pathways such as photosynthesis and respiration, which the paper cites as motivation.
• Hydrogen bonding in water and ice, relevant to electrochemistry and aqueous/environmental chemistry.
• Hydrogen technologies — the authors are at the German Aerospace Center (DLR) and Helmholtz Institute Ulm, institutions strongly associated with energy storage research, where proton transport in electrolytes and hydrogen-bonded networks is central.
• Benchmarking quantum algorithms against experimental vibrational spectra, providing a reproducible validation protocol for future quantum-chemistry hardware.

Industry relevance: The results are relevant to quantum-computing hardware and software vendors (circuit-depth-friendly ansätze, qubit tapering, active-space tooling), to the energy-storage and catalysis sectors interested in proton conduction, and to pharmaceutical/biotech R&D where proton tunneling affects reaction rates.

Future Directions

• Explaining the hydrogen-versus-water accuracy gap. The three hypotheses offered (direct quantum-nuclei interaction, high-frequency modes, small basis sets/active spaces) remain untested; the authors explicitly state this requires a broader molecule study with larger basis sets and active spaces.
• Larger nuclear basis sets. Since including p-type nuclear functions currently causes severe frequency underestimation, a scheme is needed that keeps basis centers pinned to nuclear expectation values while retaining higher-angular-momentum flexibility — for example the Lagrange-multiplier or penalty-term approach used in CNEO.
• Exploiting point-group symmetry. Z2 tapering was deliberately deferred because molecular symmetries change during dynamics and could cause inconsistent qubit counts between timesteps; a consistent scheme is an open problem.
• Running on real hardware. All results use a noiseless statevector simulator; demonstrating the algorithm on noisy devices, and extending it beyond the small systems where second quantization is favorable toward first-quantization approaches, remain future work.

Target Audience

Quantum-computing researchers working on quantum chemistry algorithms; electronic-structure theorists interested in NEO and nuclear quantum effects; method developers concerned with active spaces, qubit tapering and ansatz design for NISQ devices; and computational chemists who want a benchmarked, experimentally validated protocol for vibrational spectra from quantum simulations. Readers need a solid grasp of second quantization and VQE to follow the methods, though the comparison tables and the Zundel-ion result are accessible to a broader physical-chemistry audience.

Authors’ abstract

Nuclear quantum effects are critical for describing proton transfer and hydrogen bonding, but their incorporation into quantum chemistry calculations is often computationally prohibitive on classical hardware. A promising alternative are quantum computers due to their linear scaling in the space requirements with system size. We introduce a novel hybrid quantum-classical algorithm for ab-initio molecular dynamics that incorporates nuclear quantum effects via the nuclear-electronic orbital method. The proposed algorithm evaluates ground state energies and forces on the quantum computer using a variational quantum eigensolver, while the molecular geometries are updated classically. We validate our approach through simulations of $\text{H}_2$, $\text{H}_2\text{O}$ and the Zundel ion $\text{H}_5\text{O}_2^+$, comparing the simulated vibrational spectra with experimental data. Upon inclusion of nuclear quantum effects, the proton shuttling movement in the Zundel ion becomes effectively barrierless, and errors in the simulated frequencies improve significantly. Employing compact hardware-efficient ansatz circuits we achieve results comparable to the more accurate UCCSD ansatzes, which hints towards the feasibility of executing our algorithm on near-term quantum devices.