Precise computation of rovibronic resonances of molecular hydrogen: $EF\ ^1\Sigma_\mathrm{g}^+$ inner-well rotational states

TL;DR

This paper accurately computes the rovibronic resonances of molecular hydrogen's $EF$ state using advanced variational methods, achieving near-exact energies and refining theoretical predictions to match experimental data.

Contribution

It introduces an energy-tracking optimization for basis representation, significantly improving the precision of resonance calculations in the four-body hydrogen molecule.

Findings

Non-relativistic energies converged within a few nano Hartree

Relativistic and QED corrections align theory with experiments

Predissociative widths are negligible at high accuracy

Abstract

Selected states of the $EF\ ^1\Sigma_\mathrm{g}^+$ electronic manifold of the hydrogen molecule are computed as resonances of the four-body problem. Systematic improvement of the basis representation for the variational treatment is achieved through an energy-tracking optimization procedure. The resulting non-relativistic energy is converged within a few nano Hartree, while the predissociative width is found to be negligible at this level of accuracy. The four-particle non-relativistic energies are appended with relativistic and quantum electrodynamics corrections which close the gap between the experimental observations and earlier theoretical work.

Tables

Table 1. Table 1: Variational, non-relativistic four-particle energy, in E h subscript E h \text{E}_{\text{h}} , corresponding to the N = 0 , 1 , … , 5 𝑁 0 1 … 5 N=0,1,...,5 rotational states of the ground vibrational state in the inner well, ‘ E ​ 0 ​ N 𝐸 0 𝑁 E0N ’, of the E ​ F 𝐸 𝐹 EF   Σ g + 1 superscript superscript subscript Σ g 1 {}^{1}\Sigma_{\mathrm{g}}^{+} electronic state of H = 2 { e − , e − , p + , p + } {}_{2}=\{\text{e}^{-},\text{e}^{-},\text{p}^{+},\text{p}^{+}\} . The term values, in cm -1 , are given with respect to the rovibronic ground state. To obtain these singlet ( S e = 0 subscript 𝑆 e 0 S_{\text{e}}=0 ) states, the parity and the proton spin were chosen to be p = ( − 1 ) N 𝑝 superscript 1 𝑁 p=(-1)^{N} and S p = ( 1 − p ) / 2 subscript 𝑆 p 1 𝑝 2 S_{\text{p}}=(1-p)/2 , respectively. The E ( 2 ) superscript 𝐸 2 E^{(2)} non-relativistic energy is estimated to be converged within a few nano Hartree.

$N$	$E^{\left(2\right)}$	$T^{\left(2\right)}$ ${}^{\text{a}}$	$\Delta T_{\text{o-c}}^{\left(2\right)}$ ${}^{\text{b}}$	$\Delta T_{\text{o-c}}^{\left(2\right)}$ Yu and Dressler (1994)
0	$-$0.712 197 577	99164.664	0.123	0.320
1	$-$0.711 908 569	99228.094	0.124	0.321
2	$-$0.711 332 945	99354.429	0.128	0.304
3	$-$0.710 475 421	99542.633	0.133	0.312
4	$-$0.709 342 932	99791.186	0.138	0.32
5	$-$0.707 944 454	100098.116	0.145	0.33

Table 2. Table 2: Perturbative relativistic and QED corrections up to estimates for m ​ α 7 𝑚 superscript 𝛼 7 m\alpha^{7} (see text), in cm -1 , to the E ​ 0 ​ N 𝐸 0 𝑁 E0N – X ​ 00 𝑋 00 X00 term values of H 2 reported in Table  1 . The relativistic and QED corrections are estimated to be accurate within 10 − 3 superscript 10 3 10^{-3}  cm -1 , which results an overall uncertainty estimate ± 0.005 plus-or-minus 0.005 \pm 0.005  cm -1  for T ( 2 ​ … ​ 7 ) superscript 𝑇 2 … 7 T^{(2\ldots 7)} .

$N$	$\delta T^{\left(4\right)}$ ${}^{\text{a}}$	$\delta T^{\left(5\right)}$ ${}^{\text{b}}$	$\delta T^{\left(\mathrm{6..7}\right)}$ ${}^{\text{c}}$	$\delta T^{\left(\mathrm{4..7}\right)}$ ${}^{\text{d}}$	$T^{\left(\mathrm{2..7}\right)}$ ${}^{\text{e}}$	$\Delta T_{\text{o-c}}^{\left(\mathrm{2..7}\right)}$ ${}^{\text{f}}$
0	0.475	$-0.351$	$-0.0027$	0.122	99164.786	$0.001$
1	0.478	$-0.351$	$-0.0027$	0.124	99228.217	$0.001$
2	0.482	$-0.352$	$-0.0027$	0.127	99354.557	$0.001$
3	0.488	$-0.353$	$-0.0027$	0.132	99542.764	$0.000$
4	0.496	$-0.355$	$-0.0027$	0.138	99791.326	$0.000$
5	0.506	$-0.357$	$-0.0028$	0.146	100098.265	$-0.001$

Equations

$$E^{(\mathrm{2..7})}=E^{(2)}+\sum _{k=2}^5{\alpha }^k⟨\psi |{\mathcal{H}}^{(k+2)}|\psi ⟩$$

$${\mathcal{H}}^{(4)}=$$

$$+\pi \delta ({\mathbf{r}}_{12})-\frac{1}{2}\left[{\mathbf{p}}_1\frac{1}{r_{12}}{\mathbf{p}}_2+{\mathbf{p}}_1\cdot {\mathbf{r}}_{12}\frac{1}{r_{12}^3}{\mathbf{r}}_{12}\cdot {\mathbf{p}}_2\right],$$

$${\mathcal{H}}^{(5)}=$$

$$+\left[\frac{164}{15}+\frac{14}{3}\ln \alpha \right]\delta ({\mathbf{r}}_{12})-\frac{7}{6\pi }P(1/r_{12}^3)$$

$${\mathcal{H}}_{\text{est}}^{(6)}$$

$${\mathcal{H}}_{\text{est}}^{(7)}$$

Full text

The lowest-energy rotational and vibrational states of the ground electronic state of the hydrogen molecule have received much attention over the past decade. We have witnessed several orders of magnitude improvement in terms of accuracy and precision both on this experimental Hölsch et al. (2019); Liu et al. (2009); Altmann et al. (2018); Cheng et al. (2018) and theoretical Wang and Yan (2018); Puchalski et al. (2019) frontier of molecular physics.

The hydrogen molecule has several electronically excited states, many of them are very interesting by their own, e.g., the famous double-well features caused by avoided crossings. The rovibronic level structure is dominated by non-adiabatic interactions among the states, not yet fully, quantitatively understood. In the meanwhile, many of these electronically excited states have been measured experimentally to high precision Altmann et al. (2018); Hölsch et al. (2018), and rovibronic states corresponding to electronic excitations, e.g., the $EF$ and $GK~{}^{1}\Sigma_{\mathrm{g}}^{+}$ manifolds, have been used in the excitation sequences resulting in the ultra-precise dissociation energy of the lowest rovibrational levels of the ground electronic state Liu et al. (2009); Altmann et al. (2018); Cheng et al. (2018). Therefore, it is not only for purely theoretical interest to aim for a better and more complete theoretical description of electronically excited states of the hydrogen molecule. In particular, the entire dynamical range of molecular hydrogen, which has already been experimentally studied, spans an about 130 000 cm*-1* broad energetic and a 15 bohr broad structural (proton-proton separation) range, and includes a very large number of sharp spectral transitions which can be measured to high precision. For this reason, we think that the computation of a variety of these rovibronic states would offer an excellent testing ground for the numerous small ‘effects’ which have been identified during the study of the ground state of H2, see for example, Ref. Puchalski et al. (2019). The present work cannot aim for including all these effects at once, but we wish to provide a good starting point by significantly improving upon earlier theory for the selected states.

H2 has several, challenging excited states, some of them are bound even in the four-body treatment, e.g., the rovibronic levels of $B\ ^{1}\Sigma_{\mathrm{u}}^{+}$. In the present work, we will look at the first electronically excited state beyond the ground state, which is the $EF\ ^{1}\Sigma_{\text{g}}^{+}$ state. The Born–Oppenheimer (BO) potential energy curve of $EF\ ^{1}\Sigma_{\text{g}}^{+}$ shows a double-well feature due to an avoided crossing with the nearby $GK\ ^{1}\Sigma_{\text{g}}^{+}$ state (Figure 1). In a pre-Born–Oppenheimer (pre-BO) description Mátyus and Reiher (2012); Mátyus (2013, 2019), ‘all’ non-adiabatic couplings and ‘effects’ are automatically included, so we will not use potential energy curves, nor coupling vectors in the computations, but the curves are useful to look at and we will continue to use the electronic state labels to have a short description and reference for the computed four-particle states. Since all non-adiabatic couplings are included, also the coupling with the $X\ ^{1}\Sigma_{\text{g}}^{+}$ (ground) state continuum is present (we cannot separate it), and thus the $EF$ states can only be obtained as resonances within the four-body problem Mátyus (2013). The lower-energy vibrations of $EF$ have been estimated to have a very long predissociative lifetime (much longer than their radiative lifetime) due to their very weak coupling to the dissociation continuum of the $X\ ^{1}\Sigma_{\text{g}}^{+}$ ground electronic state Quadrelli et al. (1990); Mátyus (2013).

The non-adiabatic manifold, which includes also the $EF$ state, was computed by Yu and Dressler Yu and Dressler (1994) by explicitly coupling of nine electronic states. Yu and Dressler used accurate potential energy curves and non-adiabatic coupling vectors, and their nine-state computation resulted in rovibronic term values within 0.1–20 cm*-1* of experiment. This nine-state computation was a significant improvement upon the earlier five-state study of Quadrelli, Dressler, and Wolniewicz Quadrelli et al. (1990), which showed a larger, 1.3–120 cm*-1* deviation from the experimental results. As it was pointed out by Hölsch, Beyer, and Merkt recently Hölsch et al. (2018), performing a non-adiabatic computation with more than nine fully coupled electronic states for this system is not obvious (one would need to include many more electronic states and probably also the interaction with the H${}_{2}^{+}\ +\text{e}^{-}$ ionization continuum) but extension of non-adiabatic perturbation theory could be possible. The effective Hamiltonian for the quantum nuclear motion over coupled electronic states which perturbatively accounts for the effect of the distant electronic states (not included in the fully coupled electronic band) has been recently formulated MaT following Refs. Teufel (2003); Panati et al. (2007) and its numerical application, by generalizing the computational approach developed for a single-state non-adiabatic Hamiltonian Mátyus (2018a, b); Przybytek et al. (2017); Pachucki and Komasa (2009), should soon follow.

However, benchmark energies and wave functions are best obtained from the direct solution of the four-particle Schrödinger equation. Some time ago Mátyus (2013), one of us has identified rovibronic resonances in pre-BO computations, also from the $EF$ manifold, by using symmetry-adapted, explicitly correlated Gaussian (ECG) basis functions optimized for bound states. The global-vector representation Suzuki and Varga (1998); Suzuki et al. (1998) of ECGs made it possible to have a general basis set for $N\geq 0$ total angular momentum quantum number. At that time, a systematic improvement scheme of the basis representation for resonances was not available, which made it difficult to assess the accuracy of the results. Although a complex analogue of the real variational principle exists for the complex-coordinate rotated (CCR) Hamiltonian Kukulin et al. (1988), the practical utilization of this complex variational principle for basis function optimization is not straightforward (of course, its utility for computing the energy and the width for a given basis set is well established). In this work, we propose a practical approach inspired by the stabilization method Hazi and Taylor (1970); Kukulin et al. (1988); Mandelshtam et al. (1994), which will be useful for the systematic improvement of the basis set for long-lived resonance states that are only weakly coupled to the continuum.

It is well known that basis functions can be optimized for the $n$th excited bound state by minimizing the $n$th eigenvalue of the Hamiltonian matrix Suzuki and Varga (1998); Aznabaev et al. (2018). If we used this procedure to minimize the energy of a resonance (by picking the $n$th eigenstate in a starting basis set, which is beyond the dissociation threshold and resembles most the resonance to be computed) then, after a few basis refinement cycles, discrete representations of the continuum would start to accumulate in our energy list and we would end up minimizing the energy converging to the dissociation threshold. In order to avoid this minimization ‘collapse’, we do not focus on the $n$th state but we define an energy threshold, $\varepsilon$, the energy of which is slightly smaller than the exact energy of the state we are looking for. (This value can be estimated and does not need to be a very tight estimate. In practice, we use the experimental term to set a loose lower bound.) We use this ‘energy-tracking’ procedure to optimize the basis representation for the selected rovibronic state by minimizing the energy of the first state above the $\varepsilon$ energy threshold. It is important to add that the energy tracking (in its simplest form) works for long-lived resonances, which couple to the continuum only weakly, and the energy of which in increasing but finite basis sets stabilizes over many iteration cycles. During the refinement cycles, discrete representations of the continuum may drop in our energy list, but we track and minimize the energy of that state which is the lowest above our $\varepsilon$ threshold. We use this procedure only to optimize the basis functions. Following the basis optimization, we studied resonance features using the complex-coordinate rotation (CCR) technique Kukulin et al. (1988) and its pre-BO implementation of Ref. Mátyus (2013). The point at which the CCR trajectories stabilize in the complex energy plane defines the full (complex) molecular wave function of the resonance (the complex energy and the CCR angle). This complex wave function can then be used to compute relativistic and quantum electrodynamics (QED) corrections, using the CCR form of the correction operators, to obtain the relativistic and QED corrections for both the energy and the width (lifetime) I.Korobov (2004).

Now, we focus our attention to the computation of the rotational excitations, with $N=0,1,2,3,4,$ and $5$ total angular momentum quantum numbers, of the lowest-energy vibration in the inner well of the $EF\ ^{1}\Sigma_{\text{g}}^{+}$ electronic state. The inner-well ground vibrational state was labelled with $E0$ in Ref. Yu and Dressler (1994), so we will use $E0N$ for the $N$th rotational excitation of this state. We can access these states by choosing the appropriate non-relativistic quantum numbers, which are listed in caption to Table 1, and by defining the $\varepsilon$ energy threshold for the energy-tracking procedure, which we set to ca. 10–20 cm*-1* lower than the energy estimated from the experimental term value. Using the QUANTEN computer program Mátyus and Reiher (2012); Mátyus (2013, 2019) (see also the Supplementary Material som ), we optimized, in repeated refinement cycles, a starting basis set which was compiled from the extensive bound-state optimization work of Ref. Mátyus (2013). In retrospect, we can confirm that the largest basis set of Ref. Mátyus (2013) gives $EF$ energies accurate within ca. 10 nE${}_{\text{h}}$, which we can estimate now from the convergence pattern and behaviour of the states upon the systematic refinement of the basis set.

In Table 1, we list the non-relativistic energy values optimized in the present work ($E^{(2)}$ column), which are estimated to be within a few times $1\ \text{nE}_{\text{h}}=10^{-9}\ \text{E}_{\text{h}}$ of the exact, non-relativistic value. The non-relativistic term values, the difference of the $E0N$ energies and the non-relativistic energy of the ground state ($X00$), reduce the 0.3 cm*-1* deviation of the nine-state non-adiabatic computation of Yu and Dressler Yu and Dressler (1994) to 0.1 cm*-1*, and confirm their electronic energy error estimate. Yu and Dressler also estimated the relativistic and QED corrections to be ca. 0.08 cm*-1* for the $E0N$–$X00$ terms. This is an average value for the states with different rotational quantum numbers, and was compiled from the expectation value of the Breit–Pauli Hamiltonian with the electronic wave function at $R=1.9$ bohr Wolniewicz and Dressler (1994) and the QED correction of H${}_{2}^{+}$ Bukowski et al. (1992) at the $R=1.9$ bohr proton-proton distance, which is near the effective structure of the $E00$ state. If we correct our non-relativistic term values with these estimates, then the deviation from experiment reduces to $(0.035,0.036,0.040,0.045,0.050,0.057)$ cm*-1* for $N=0,1,2,\ldots,5$, respectively.

We would like to have a more complete account for the relativistic and QED effects and also to possibly know the $N$ dependence of the correction. For this purpose, we compiled data from the literature Wolniewicz (1998) and carried out additional computations using the non-relativistic, pre-BO wave functions computed in the present work. The relativistic, leading and higher-order QED effects (we have explicitly considered estimates up to the so-called $m\alpha^{7}$ terms) are calculated as perturbative corrections using the non-relativistic energy and wave function, $E^{(2)}$ and $\psi$, in terms of increasing orders of the fine-structure constant, $\alpha$,

[TABLE]

where the $\mathcal{H}^{(k+2)}$ operators are reproduced from the literature in the following paragraphs (with the usual meaning of the symbols and operators, the details can be found in the corresponding references). The underlying integrals, to be described in the following paragraphs, are evaluated so that the uncertainty of each correction term is better than $0.001$ cm*-1*.

To calculate the (spin-independent) relativistic correction, we have started out from the expectation value of the Breit–Pauli Hamiltonian (of the electrons) Bethe and Salpeter (1975); Puchalski et al. (2017)

[TABLE]

where the electrons are labelled with ‘1’ and ‘2’, while the protons are labelled with ‘3’ and ‘4’. Wolniewicz already calculated the expectation value of $\mathcal{H}^{(4)}$ with the electronic wave function along a series of nuclear separations Wolniewicz (1998); Wof . We obtained the relativistic correction to each $E0N$ state by evaluating the expectation value of the BO relativistic correction curve (represented with polynomial fits) with the pre-BO wave functions. To obtain the term corrections, $\delta T^{(4)}$ (Table 2), we used the relativistic correction value, $-1.652$ cm*-1*, of the $X00$ ground state derived from a similar level of theory Puchalski et al. (2017). Note that this value is 0.002 cm*-1* smaller than the correction calculated directly with the Breit–Pauli Hamiltonian of the electrons and protons and the four-particle wave function of the $X00$ state Wang and Yan (2018), which will have to be accounted for when we estimate the uncertainties of the present results.

The leading QED contribution (to the electronic part of the problem) Araki (1957); Sucher (1958); Bethe and Salpeter (1975); Komasa et al. (2011) is

[TABLE]

which we evaluated with accurate electronic wave functions along a series of nuclear configurations. The $\ln K$ non-relativistic Bethe-logarithm was also treated within the BO approximation similarly to Ref. Piszczatowski et al. (2009). The $\ln K(R)$ values for the $EF$ electronic state were approximated with the $\ln K(R)$ function of the ion core of $EF$, so we could use the accurate $\ln K(R)$ values of (the lowest electronic state of) H${}_{2}^{+}$ computed by Korobov Korobov et al. (2013); Aznabaev et al. (2018). The Dirac delta terms containing the electron-proton and electron-electron displacement vectors were obtained similarly to the Darwin terms of the relativistic correction, i.e., by computing the expectation value of the $R$-dependent correction curves with the pre-BO wave function for each $N$. The last term in Eq. (3) is the Araki–Sucher (AS) correction, which we computed for the $EF$ state in the present work using accurate electronic wave functions (obtained within the BO module of QUANTEN using floating ECGs Mátyus (2018b, a)) and the integral transformation technique Pachucki et al. (2005). The AS correction to each $E0N$ state was obtained as the expectation value of the correction curve with the four-particle wave function, and was found to be an order of magnitude smaller, $-0.001$ cm*-1*, than the correction for the ground state, $-0.013$ cm*-1* Komasa et al. (2011); Puchalski et al. (2019). By summing up all these contributions, we obtain the leading QED correction to the energy, which changes from $0.385$ to $0.379$ cm*-1* as $N$ increases from 0 to 5. To calculate the $\delta T^{(5)}$ leading QED term corrections listed in Table 2 we used the $0.736$ cm*-1* value of the ground state compiled from Refs. Puchalski et al. (2017, 2019).

Higher-order QED corrections were estimated at the one-loop level by retaining only those terms which give the dominant corrections at these orders Piszczatowski et al. (2009); Puchalski et al. (2016, 2019)

[TABLE]

While the $\mathcal{H}^{(6)}_{\text{est}}$ contribution to the $E0N$–$X00$ term values is $-0.003$ cm*-1*, the $\mathcal{H}^{(7)}_{\text{est}}$ changes the terms by as little as $2\cdot 10^{-4}$ cm*-1*, which is negligible given the uncertainty of the current evaluation of the relativistic and QED integrals.

The overall $\delta T^{(4..7)}$ contribution (Table 2) of the relativistic, leading and higher-order QED effects increases from 0.122 to 0.146 cm*-1* upon the increase of $N=0$ to 5. The resulting $T^{(2..7)}$ term values for $N=0,\ldots,5$ show $\pm 0.001$ cm*-1* deviations from the experimental values of Ref. Dickenson et al. (2012), which is better (probably fortuitous) than the uncertainty of the present theoretical values, which are thought to be accurate within about $\pm 0.005$ cm*-1*.

The experimental values of Ref. Dickenson et al. (2012) are more precise than our theoretical results, and the additional significant digits surely hide interesting physics, so theory should aim for further improvements. In order to help future work, we close this article with commenting on the possible sources of uncertainties in our work.

First of all, the non-relativistic energy was obtained in a variational procedure (stabilized for long-lived resonances), and systematic improvement of this value is rather straightforward. The convergence pattern observed in repeated rounds of refinement cycles suggests that the non-relativistic energy, $E^{(2)}$, is converged within $0.000\,5$ cm*-1*. Assessment of the uncertainty of the relativistic and QED corrections is more delicate. In the light of the developments of recent years for the rovibronic ground state Puchalski et al. (2017); Wang and Yan (2018); Puchalski et al. (2016, 2019), we think that the largest source of error in our work must be due to the relativistic ‘recoil’ effect, on the order of a few $10^{-3}$ cm*-1*, which in simple terms means that (at least) the relativistic corrections should be computed using the full electron-nucleus Breit–Pauli Hamiltonian and the four-particle wave functions Wang and Yan (2018); Puchalski et al. (2019). Then, in order to pinpoint one-two more digits in the calculations, the current approximations used for the non-relativistic Bethe-logarithm term will have to be checked and the contribution of the neglected terms in the higher-order QED corrections (in particular, the full $m\alpha^{6}$ contribution of $\mathcal{H}^{(6)}$) will have to be elaborated.

In the usual perturbative manner, relativistic quantum electrodynamics (and possibly beyond) is adapted to molecular computations, it is necessary to evaluate and sum several, small (and often not so small) contributions (of different signs) on top of a direct, variational non-relativistic computation. We think that the extremely rich excited state, rovibronic level structure of the hydrogen molecule (Figure 1) offers an excellent opportunity to challenge and cross-check the theoretical and computational procedures both in terms of the completeness of the physical description as well as regarding the error balance of possible uncertainties and inaccuracies. The present work demonstrates that electronically excited states of H2 can be theoretically described to high precision and, with further improvements, they will provide equally useful and at the same time complementary information to the study of the ground electronic state.

Acknowledgement

We acknowledge financial support from a PROMYS Grant (no. IZ11Z0_166525) of the Swiss National Science Foundation. DF thanks a doctoral scholarship from the New National Excellence Program of the Ministry of Human Capacities of Hungary (ÚNKP-18-3-II-ELTE-133). EM is thankful to ETH Zürich for supporting a stay as visiting professor during 2019 and the Laboratory of Physical Chemistry for their hospitality, where part of this work has been completed.