Exact tunneling splittings of rotationally excited states from symmetrized path-integral molecular dynamics

This paper presents a numerically exact, symmetrized path-integral molecular dynamics method that utilizes an Eckart spring to rigorously project molecular systems onto specific rotational manifolds, enabling the efficient, simultaneous calculation of rotationally resolved tunneling splittings for multiple angular momentum states with results that closely match experimental trends and variational benchmarks.

The Gist

The Big Picture: The Quantum "Ghost" in the Machine

Imagine a molecule, like a tiny spinning top made of atoms. Sometimes, this top can flip upside down. In the classical world (the world of baseballs and cars), if the top doesn't have enough energy to climb over a hill to flip, it just sits there.

But in the quantum world, atoms are like ghosts. They can tunnel through the hill instead of going over it. When a molecule tunnels back and forth between two shapes, it creates a tiny "split" in its energy levels. Scientists call this a tunneling splitting. Measuring this split is like taking a fingerprint of the molecule's internal landscape; it tells us exactly how high and steep the "hill" is.

The problem? Calculating these splits for molecules that are spinning fast (rotationally excited) is incredibly hard. It's like trying to predict the path of a spinning, tumbling acrobat while they are also trying to walk through a wall.

This paper introduces a new, super-precise method to calculate these splits for spinning molecules, and it does it with surprising efficiency.

---

The Old Way vs. The New Way

The Old Problem: The "Blind" Simulation

Imagine you want to count how many people are in a specific room of a giant, chaotic party.

• Traditional methods were like sending a camera into the party. But the camera was "blind" to who was spinning. It saw everyone: the dancers, the spinners, the jumpers. To get the count for just the spinners, you had to run the simulation, then run it again for a different group, then again for another. It was slow, expensive, and often got the numbers wrong because the "noise" of the spinning messed up the count.

The New Solution: The "Eckart Spring" Filter

The authors (Léa Zupan, Yu-Chen Wang, and Jeremy Richardson) built a special filter. Think of it as a magic bouncer at the door of the party.

1. The Ring Polymer: In their simulation, a molecule isn't just one point; it's a string of beads (a ring polymer) representing its quantum nature.
2. The Magic Spring: They attached a special "Eckart spring" between the start and end of this string. This spring doesn't just pull the ends together; it pulls them together only if the molecule is oriented in a specific way and spinning at a specific speed.
3. The Result: If the molecule is spinning too fast or in the wrong direction, the spring pulls it back. If it's in the right "rotational state," the spring lets it pass.

This allows them to isolate only the molecules spinning at a specific speed (called quantum number $J$) without having to run separate, expensive simulations for each speed.

The "One-Stop Shop" Trick

Here is the most clever part of their discovery.

Usually, to find out how much energy a molecule has at spin-speed 1, then spin-speed 2, then spin-speed 3, you would need to run three different computer simulations.

This new method is like a "One-Stop Shop."
They run the simulation once. Because of the math behind their "magic spring," they can take the data from that single run and, through a bit of post-processing (like sorting a deck of cards), extract the answers for all the different spin speeds ($J=0, 1, 2, 3...$) at the same time.

It's like baking one giant cake, but instead of cutting it into slices, you magically get a separate, perfect cake for every flavor you wanted, all from the same batter.

Testing the Theory: Water and Ammonia

To prove their method works, they tested it on two famous molecules:

1. Water ($H_2O$): Water doesn't really tunnel in a way that splits its energy levels, but it spins. They used water to check if their "magic spring" could correctly calculate the energy differences between different spinning states.

   • The Result: It worked perfectly. Their numbers matched the most precise experiments and other super-complex math methods. It proved the "bouncer" was doing its job correctly.

2. Ammonia ($NH_3$): This is the classic "tunneling" molecule. It flips like an umbrella in the wind.

   • The Result: They calculated the tunneling splits for ammonia spinning at different speeds. Their results matched the "gold standard" (exact quantum math) almost perfectly.
   • The Trend: They confirmed a long-standing observation: as ammonia spins faster, the tunneling split gets smaller. Their method captured this trend beautifully.

Why Does This Matter?

• Speed and Efficiency: Before this, studying spinning molecules was a computational nightmare. Now, scientists can get answers for many different spinning states from a single simulation.
• Accuracy: It is "numerically exact," meaning the only errors come from the computer's statistical noise, not from bad approximations.
• Future Applications: While they tested it on small molecules (water and ammonia), this method scales up. It can be used for huge, floppy molecules (like proteins or complex clusters) that are too big for other methods to handle.

The Bottom Line

The authors have invented a quantum filter that lets them listen to the specific "song" of a molecule spinning at a certain speed, even while it's tunneling through energy barriers. They proved that you don't need a separate concert for every instrument; you can record the whole orchestra once and extract the solo for any instrument you want. This opens the door to understanding the quantum behavior of much larger and more complex molecules than ever before.

Technical Summary

1. Problem Statement

Quantum tunneling is a fundamental phenomenon where particles traverse classically forbidden potential energy barriers, leading to the lifting of degeneracy in symmetry-related states (tunneling splittings). While these splittings are crucial for understanding molecular dynamics and reaction rates, calculating them for rotationally excited states ($J > 0$) remains a significant challenge:

• Variational Wavefunction Methods: Highly accurate but suffer from the "curse of dimensionality," making them computationally intractable for larger systems or highly fluxional molecules.
• Diffusion Monte Carlo (DMC): Requires separate calculations for ground and excited states and the construction of nodal surfaces, introducing systematic errors that are difficult to estimate, especially for excited rotational states.
• Instanton Theory: A semiclassical approximation that is efficient but struggles to achieve quantitative accuracy for weakly bound dimers or low barriers, and extensions to rotationally excited states are not yet fully validated.
• Standard Path-Integral Molecular Dynamics (PIMD): Previous PIMD approaches calculated tunneling splittings without explicitly projecting onto specific rotational manifolds. Consequently, the results were often contaminated by contributions from rotational excitations, preventing the isolation of pure tunneling effects for specific $J$ states.

The core problem addressed is the lack of a numerically exact, computationally efficient method to calculate tunneling splittings for rotationally excited states ($J > 0$) that scales favorably with system size and handles strong rovibrational coupling.

2. Methodology

The authors extend the Symmetrized Path-Integral Molecular Dynamics (PIMD) framework (previously developed for $J=0$) to rigorously project the system onto selected rotational manifolds ($J$) and discrete symmetry subspaces.

A. Theoretical Framework

• Rotational Projection: The canonical partition function $Z$ is projected onto a specific total angular momentum quantum number $J$ by averaging over the degenerate magnetic quantum numbers $M$. This is achieved using the orthogonality of Wigner $D$-matrices:
  $$Z^{(J)} = \frac{1}{8\pi^2} \int d\Omega \, \text{Tr}[D^{(J)}(\Omega)] \, Z(\Omega)$$
  where $Z(\Omega)$ is the partition function with an overall rotation $\Omega$.
• Symmetrization: To isolate specific symmetry states (irreducible representations), the method incorporates permutation-inversion operations ($\hat{P}$) into the partition function, creating symmetrized partition functions $Z_P$.
• Combined Projection: The final quantity of interest is the $J$-projected, symmetrized partition function $Z_P^{(J)}$, which isolates the lowest-lying states of a specific symmetry within a specific rotational manifold.

B. The Eckart Spring and PIMD Implementation

• Ring Polymer Representation: The quantum system is mapped to a classical ring polymer. The boundary condition connects the first bead ($r^{(1)}$) and the last bead ($r^{(N)}$) via a symmetry operation $\hat{P}$ and a rotation $\Omega$.
• The Eckart Spring: A harmonic spring connects $\hat{P}^{-1}r^{(N)}$ to the rotated configuration $\hat{R}(\Omega)r^{(1)}$. This "Eckart spring" enforces the optimal rotation that minimizes the distance between the beads, effectively satisfying Eckart conditions for a rotating molecule.
• Steepest-Descent Approximation: The integral over the rotation angle $\Omega$ is evaluated using a steepest-descent approximation. As the number of beads $N \to \infty$, this becomes exact. The rotation $\Omega$ is replaced by the optimal rotation $\tilde{\Omega}$ (found via quaternion algebra), and the integral yields a prefactor involving the trace of the Wigner $D$-matrix, $\text{Tr}[D^{(J)}(\tilde{\Omega})]$.
• Key Insight: The effective Hamiltonian used for the PIMD propagation is independent of $J$. The dependence on $J$ appears only in a post-processing prefactor. This allows a single set of PIMD simulations to generate data for multiple $J$ values simultaneously.

C. Free Energy and Splitting Calculation

• Thermodynamic Integration (TI): To compute the ratios of partition functions required for energy differences, the authors use TI to calculate the free energy difference ($\Delta F$) between different effective Hamiltonians (e.g., connecting identity to a permutation operation).
• Tunneling Splitting Formula: The energy splitting between two states is derived from the ratio of their projected partition functions:
  $$\Delta E = \frac{1}{\beta} \ln \left( \frac{\sum g_\alpha \chi^* Z^{(J)}_\alpha}{\sum g_\alpha \chi^* Z^{(J)}_\alpha} \right)$$
  This reduces the problem to evaluating ratios of symmetrized partition functions, which are factorized into free-energy contributions and $J$-dependent prefactors.

3. Key Contributions

1. Extension to $J > 0$: The first rigorous derivation and implementation of PIMD that projects onto rotationally excited states ($J > 0$) while maintaining numerical exactness (within statistical uncertainty).
2. Computational Efficiency: The method allows the extraction of tunneling splittings for multiple $J$ values from a single simulation trajectory, as the $J$-dependence is confined to a post-processing prefactor. This eliminates the need for separate, costly simulations for each rotational state.
3. Rigorous Treatment of Rovibrational Coupling: Unlike the rigid-rotor approximation, this method naturally includes all rovibrational coupling terms, which are critical for accurate predictions in excited states.
4. Validation on Benchmark Systems: The formalism was successfully validated on water ($H_2O$) and ammonia ($NH_3$), showing excellent agreement with exact variational benchmarks and experimental data.

4. Results

The method was applied to two molecular systems:

A. Water ($H_2O$)

• Goal: Validate the rotational projection formalism by computing energy differences between $J=0$ and $J=1$ manifolds.
• Findings:
  • Calculated energy differences between the three $J=1$ states ($1_{01}, 1_{10}, 1_{11}$) agreed with exact variational benchmarks and experimental values within statistical error.
  • The method correctly captured rovibrational coupling effects, outperforming the rigid-rotor approximation (which underestimated energy differences by 3–4%).
  • Convergence with respect to the number of beads ($N$) was achieved, with $N=256$ providing stable results.

B. Ammonia ($NH_3$)

• Goal: Calculate rotationally resolved tunneling splittings for the umbrella inversion mode for $J = 0, \dots, 4$ and $J=1, K=1$.
• Findings:
  • Accuracy: The PIMD results for tunneling splittings ($\Delta$) matched numerically exact variational wavefunction calculations (on the same Potential Energy Surface) within statistical uncertainty.
  • Trend: The results reproduced the experimentally observed trend that tunneling splitting decreases as $J$ increases.
  • PES Limitation: While the theoretical values were numerically exact for the given Potential Energy Surface (PES), they slightly underestimated absolute experimental values. This discrepancy was attributed to imperfections in the underlying PES rather than the methodology.
  • Statistical Uncertainty: Uncertainties increased with $J$ due to the oscillatory nature of the Wigner $D$-matrix trace and the slow sampling of relative end-to-end bead orientations, but remained manageable.

5. Significance

• Scalability: Unlike variational methods, this PIMD approach scales favorably with system size, making it applicable to larger, highly fluxional molecules (e.g., $CH_5^+$, formic acid dimers) where traditional methods fail.
• Exactness: It provides a route to "numerically exact" tunneling splittings (limited only by statistical error and the quality of the PES) for rotationally excited states, a regime previously inaccessible to rigorous quantum dynamics.
• Versatility: The method is formulated in Cartesian coordinates, making it readily applicable to a wide range of molecular systems without the need for complex internal coordinate definitions.
• Future Impact: This work establishes a foundation for studying non-adiabatic effects, complex reaction rates involving tunneling, and spectroscopic properties of excited rotational states in complex molecular environments. It also suggests potential improvements for semiclassical instanton theory by providing exact benchmarks for rotationally excited states.

In conclusion, this paper presents a breakthrough in computational quantum chemistry, offering a robust, exact, and efficient tool to probe the interplay between rotation, vibration, and tunneling in molecular systems.