Molecular anions have recently been detected in the interstellar and circumstellar media. Accurate modeling of their abundance requires calculations of collisional data with the most abundant species that are usually He atoms and H2 molecules. In this paper, we focus on the collisional excitation of the first observed molecular anion, C6H, by He and H2. Theoretical calculations of collisional cross sections rely generally on ab initio interaction potential energy surfaces (PESs). Hence, we present here the first PESs for the C6H–H2 and C6H–He van der Waals systems. The ab initio energy data for the surfaces were computed at the explicitly correlated coupled cluster with single, double, and scaled perturbative triple excitations level of theory. The method of interpolating moving least squares was used to construct 4D and 2D analytical PESs from these data. Both surfaces are characterized by deep wells and large anisotropies. Analytical models of the PESs were used in scattering calculations to obtain cross sections for low-lying rotational transitions. As could have been anticipated, important differences exist between the He and H2 cross sections. Conversely, no significant differences exist between the collisions of C6H with the two species of H2 (para- and ortho-H2). We expect that these new data will help in accurately determining the abundance of the C6H anions in space.

The interstellar medium (ISM) is host to a number of molecules including both neutral and ionic species. Although linear carbon-chain anions have been included in interstellar chemistry models for some time, their detection has only recently occurred due to the previous lack of fundamental laboratory data.

The first observed molecular anion was C6H, detected in the circumstellar envelope of the evolved carbon-rich star IRC +10216 and within the dense cold molecular cloud TMC-1.1 Since then, the identification of five other molecular anions, C4H, C8H, CN, C3N, and C5N, has been made possible mainly through laboratory spectroscopic data.2–6 The detection of these anions and subsequent abundance analysis is not only important in constraining the chemical network of the ISM, the most modern of which includes over 4000 chemical reactions and 400 species,7 but their presence also directly impacts the free electron density and therefore affects the rates of cloud collapse and star formation.

However, to model the physical and chemical conditions in astrophysical environments containing anions, collisional rate coefficients with the dominant colliders are needed. Indeed, emission spectra are interpreted through detailed radiative transfer calculations, which requires the knowledge of both collisional and radiative rates. If theoretically calculated or experimentally generated collisional rate coefficients are unavailable, the quantum level populations are usually estimated by assuming local thermodynamic equilibrium (LTE), the approximation that the level populations follow a Boltzmann distribution. Generally, though, this is a poor approximation in cool, low-density regions where anions are detected.

Rate coefficients for the rotational excitation of C2H and OH molecules due to collision with He have been computed recently by Dumouchel et al.8 and Hauser et al.,9 respectively. However, these molecules, even if suspected to be in the ISM, have not been detected and the only fully relevant astrophysical data are those for the CN–H2 collisional system computed by Kłos and Lique.10 The CN–H2 rate coefficients were obtained from a new reliable ab initio potential energy surface of the CN–H2 complex and the quantum scattering calculations were performed to investigate rotational energy transfer in collisions of the CN molecule with both para-H2 (j = 0) and ortho-H2 (j = 1) molecules. However, no experimental scattering data for anionic systems of interstellar interest have yet been obtained.

Hence, there is a real lack of collisional data for (interstellar) molecular anions that prevent accurate analysis of the anionic emission spectra that can be recorded with highly resolved telescopes like the Atacama Large Millimeter/submillimeter Array (ALMA) interferometer. Among the interstellar molecular anions, C6H is of particular interest. C6H has been found in several star-forming regions1,11–13 and is one of the most abundant anionic species. He and H2 are generally the dominant colliders in interstellar and circumstellar media, the predominant one being H2.

Our aim is then to compute collisional rate coefficients for the C6H–He and C6H–H2 systems in order to provide the astrophysical community with data for a second detected anion besides CN. Within the Born-Oppenheimer approximation, the study of inelastic collisions requires two steps: (i) the calculation of an ab initio PES describing the interactions between the particles in collision (ii) the study of the dynamics of nuclei on this surface. Accurate rate coefficient calculations require a high-quality PES. Interaction PESs between long carbon chains and He or H2 are complex and difficult to represent because of the size of the target.14 Indeed, for short intermolecular distances, the interaction is typically moderate to possibly weakly attractive for a T-shape approach and is often extremely repulsive upon linear approach. This may lead to singularities in the angular expansion and severe oscillations in the numerical fit of the PES over the usual Legendre expansions that are used in quantum scattering calculations.

In this paper, we present the first PES for the C6H–H2 and C6H–He collisional systems. Special care has been put to their computation in order to have analytical PESs that can be used with confidence in quantum scattering calculations. Details of the electronic structure calculations are given in Section II A, while Section II B describes the surface fitting procedure. As an application, the first dynamics calculations are presented in Section III and the conclusions and prospects for the future are given in Section IV.

For the C6H–H2 system, both monomers were held rigid in collinear arrangements with the internuclear bond distance for H2 fixed at rHH = 0.766 65 Å (the vibrationally averaged bond distance for para- hydrogen, j = 0). The geometry of C6H was optimized at the explicitly correlated coupled cluster with single, double, and scaled perturbative triple excitations level of theory using an explicit correlation consistent valence triple zeta basis set [CCSD(T*)-F12b/VTZ-F12].15,16 The alternating pattern of bond distance parameters and Jacobi coordinates is shown in Figure 1.

FIG. 1.

The four intermolecular coordinates are illustrated. The optimized bond distance parameters for C6H are listed at left. The rHH distance at right is the empirical vibrationally averaged distance for para- hydrogen.

FIG. 1.

The four intermolecular coordinates are illustrated. The optimized bond distance parameters for C6H are listed at left. The rHH distance at right is the empirical vibrationally averaged distance for para- hydrogen.

Close modal

Since the PES represents the van der Waals (vdW) interaction of two closed shell species, single reference coupled-cluster based methods were chosen. The Molpro electronic structure package was used for all of the calculations reported here.17 In previous studies we have compared the performance of standard and explicitly correlated (F12) coupled cluster for a range of basis sets.18–20 Table I compares the interaction energy at the global minimum (collinear structure with intermonomer separation of 6.426 Å) obtained by several high-level methods. The higher symmetry for the collinear structure enables benchmarking with large basis sets and testing the effect of core-correlation. Explicitly correlated CCSD(T)-F12b calculations were performed using the VTZ-F12 and VQZ-F12 bases, correlating the valence electrons and the CVTZ-F12 and CVQZ-F12 bases, correlating all electrons. Standard CCSD(T) results for valence correlation only (since the effect of core-correlation was found to be negligible) were produced using the AVTZ and AVQZ bases for comparison. The explicitly correlated results are already remarkably well converged at the triple zeta level and indicate negligible contribution from core-correlation. Well depths of −710.5 cm−1 and −710.7 cm−1 were obtained at the VTZ-F12 and CVQZ-F12 levels, respectively. These are both close to the value of −712.1 cm−1 on the fitted surface. In contrast, the standard CCSD(T) interaction energy (without counterpoise correction or mid-bond functions) is −765.6 cm−1 for AVTZ and −707.0 cm−1 for AVQZ. Given the large number of electrons in this system (40), an explicitly correlated method was chosen in order to obtain well converged energies with a moderate basis set at relatively low cost. We and others21 have noted some slight numerical issues with F12 methods when considering the finest details of the long range part of the potential for neutral systems. No discrepancies in our data were noted for this ion-neutral system ranging out to 22 Å separation. The perturbative triples (T) contribution is not directly included in the F12b explicit correlation formalism and so in some cases we have used basis extrapolation in order to better converge that contribution. Here we scaled the (T) contribution based on the MP2-F12/MP2 correlation-energy ratio from second-order Møller-Plesset perturbation theory, implemented in Molpro as (T*). A lower-level guide surface using an explicit correlation consistent valence double zeta basis set [CCSD(T*)-F12a/VDZ-F12] was used to avoid computing expensive high-quality data in high-energy regions. The T1-diagnostic was monitored and found to be roughly 0.017 for all geometries in the high-level data set. For the C6H–He system (with the same number of electrons, 40) the same level of theory was used, but without a low-level guide surface. Since default auxiliary basis sets necessary for F12 calculations are not defined for helium with the VTZ-F12 basis in Molpro, the augmented correlation consistent valence quintuple zeta basis from second-order Møller-Plesset perturbation theory (av5z/mp2fit)22 and segmented contracted highly polarized quadruple zeta valence quality (def2-qzvpp/jkfit)23 basis sets were specified for the density fitting (DF) and resolution of the identity (RI) bases, respectively. In contrast to the good numerical behavior of the CCSD(T*)-F12b calculations for the C6H–H2 system out to 22 Å separation, some small numerical discrepancies were noted beyond 17 Å for C6H–He. Tests with various other auxiliary bases did not improve the behavior. Thus, to determine an analytic long-range for the 2D C6H–He system, standard CCSD(T)/AVTZ and CCSD(T)/AVQZ energies were extrapolated to the complete basis set (CBS) limit. The CCSD(T)/CBS data matched the CCSD(T*)-F12b/VTZ-F12 interaction energies very closely in the range of 9–15 Å, allowing a smooth switch to be applied (see below).

TABLE I.

Collinear arrangement of rigid fragments with center-of-mass separation of 6.426 Å.

MethodWell depth (cm−1)
CCSD(T)-F12b/VTZ-F12 −710.5 
CCSD(T)-F12b/VQZ-F12 −710.7 
(AE)-CCSD(T)-F12b/CVTZ-F12 −709.4 
(AE)-CCSD(T)-F12b/CVQZ-F12 −711.9 
CCSD(T)/AVTZ −765.6 
CCSD(T)/AVQZ −707.0 
MethodWell depth (cm−1)
CCSD(T)-F12b/VTZ-F12 −710.5 
CCSD(T)-F12b/VQZ-F12 −710.7 
(AE)-CCSD(T)-F12b/CVTZ-F12 −709.4 
(AE)-CCSD(T)-F12b/CVQZ-F12 −711.9 
CCSD(T)/AVTZ −765.6 
CCSD(T)/AVQZ −707.0 

The four dimensional (4D) intermolecular C6H–H2 potential includes 4691 symmetry-unique high-level ab initio energy data and is represented analytically by the interpolating moving least squares (IMLS) method19,24,25 using a weight function to interpolate between local fitting basis expansions. The C6H–H2 system is similar to the previously studied CO2–CS2,26 (NNO)2,19,27 (OCS)2,28 and (CO)220,29 systems from a fitting standpoint (weakly interacting rigid linear monomers) but has some differences with respect to symmetry. With two different monomers, the system lacks monomer exchange symmetry, but the H2 monomer is symmetric with respect to exchange of the two end-atom nuclei (180° rotation). The fitting basis (mentioned below) can be adapted to treat this symmetry by placing a simple constraint on the basis indices. As discussed previously,26,28 this is more complicated when the basis is used interpolatively (with changing weights) so here we employ the simple procedure of adding the symmetry partner for each symmetry-unique ab initio data point to the fitting set (flipping the H2 fragment). There is no additional cost in terms of electronic structure calculations and for cases of relatively low permutation symmetry (a factor of only two here) the fitting set does not become too unwieldy. For systems with very high permutation symmetry, the development of a permutation invariant basis would be preferred.30 The automated procedure that was developed to construct 4D PESs for CO2–CS2, (NNO)2, (OCS)2, and (CO)2, has been described in detail previously,19,20,25,26 and was employed here. The same inter-monomer coordinates and a fitting basis of 301 functions composed of products of radial functions with associated Legendre bend functions were used. The same distance metric, interpolative weight function, and SVD-based dynamic conditioning procedure were also used. For C6H–H2 the range of inter-monomer center-of-mass distances was R = [2.3, 22.0] Å, while the fitted energy range included all stable geometries (the global minimum has a well depth of −712.1 cm−1), but was restricted to about 2100 cm−1 above the separated monomers asymptote. As was done previously, to avoid computing and discarding costly high-level ab initio data in highly repulsive regions, an initial lower-level guide surface was constructed (at the CCSD(T*)-F12a/VDZ-F12 level). For the low-level surface, a set of 3000 symmetry unique points was distributed according to a Sobol sequence31 subject to an exponential R-dependent bias that favors points at R = 2.3 Å over points at R = 22.0 Å by a factor of about 20 (making the short-range repulsive region much more densely sampled). As before, the guide surface was fit using the same IMLS scheme as the final high-level PES, but with a smaller fitting basis of only 40 functions per local expansion. For the high-level PES (at the CCSD(T*)-F12b/VTZ-F12 level), 2000 initial seed points were distributed the same way according to the exponentially radially biased Sobol sequence, but with high-energy regions excluded by the lower level guide surface. Starting from the 2000 seed points, sets of 48 automatically determined points were added in each of a series of iterations until the estimated RMS fitting error was reduced to below 1.0 cm−1. The accuracy of the final PES was tested using a random set of 288 points, confirming the estimated sub-wavenumber accuracy. The PES generation algorithm and fitting error estimate method have been described previously and were applied here with the entire coordinate and energy range fit without bias, in an automated fashion. The PES generation algorithm was terminated and finalized for use in this case with a total of 4691 symmetry unique points. This is considerably more than the roughly 2000 points used to fit other 4D systems.19,20,26 The necessity for more points is likely due to the extreme anisotropy in the PES plotted in Figures 2 and 3.

FIG. 2.

Plot of C6H–H2 PES with θ2 = 0 (end-on approach).

FIG. 2.

Plot of C6H–H2 PES with θ2 = 0 (end-on approach).

Close modal
FIG. 3.

Plot of C6H–H2 PES with θ2 = 90, and ϕ = 0 (side-on, planar approach).

FIG. 3.

Plot of C6H–H2 PES with θ2 = 90, and ϕ = 0 (side-on, planar approach).

Close modal

The long length of the C6H fragment introduces anisotropy with respect to approach of the H2 fragment either towards the ends or side of C6H. In addition, there is a remarkably strong anisotropy with respect to the orientation of the H2 fragment (see Figures 2 and 3). End-on approach (θ2 fixed at 0 or π) of H2 includes geometries such as the completely collinear global minimum (−712.1 cm−1) with H2 interacting with the C-end of the C6H fragment. The side-on approach of H2 (θ2 fixed at π/2), however, contrasts with the end-on approach since the completely collinear minimum disappears.

To fit the 2D PES representing the interaction potential of the C6H–He system, a similar, but simplified procedure, was used. A pruned fitting basis of 39 functions (radial and Legendre, each with maximum order of 10) was used for the interpolation. Similar to the 4D PES described above, an initial radially biased set of points was determined, again using a Sobol sequence. Beginning with 400 points in the range R = [2.3, 17] Å, generations of automatically placed points were added, stopping at 953 points when the mean fitting error was below 0.3 cm−1. As mentioned above, a long range analytic form was determined using CCSD(T)/CBS data. An R-dependent switch of the form tanh(3.5∗(R − 11)) was used to smoothly switch between the IMLS and analytic forms. Thus, the PES is a 50:50 mixture of the two contributions at 11 Å, and the switch is numerically complete before 15 Å. The C6H–He PES is plotted in Figure 4. The global minimum energy found on the fitted PES is −68.76 cm−1 at R = 3.5055 Å and θ1 = 74.23°.

FIG. 4.

Plot of 2D C6H–He PES. The global minimum energy (−68.76 cm−1) is found at R = 3.5055 Å and θ1 = 74.23°.

FIG. 4.

Plot of 2D C6H–He PES. The global minimum energy (−68.76 cm−1) is found at R = 3.5055 Å and θ1 = 74.23°.

Close modal

We consider collisions of C6H with para-H2 and ortho-H2 such as

C6H(j1)+o-;p-H2(j2)C6H(j1)+o-;p-H2(j2)
(1)

where j1 and j2 designate the rotational levels of C6H and H2, respectively.

The calculations are restricted to low/moderate collisional energies (Ec < 500 cm−1). The rotational energy levels of the C6H and H2 molecules were computed using the experimental spectroscopic constants of BeC6H=0.045927, αC6H = 3.3356 × 10−5, DeC6H=1.079×109, and BeH2=60.853, αH2 = 3.062, DeH2=0.0471.1,32,33

We also consider collisions of the C6H anion with helium,

C6H(j1)+HeC6H(j1)+He,
(2)

where j is the orbital angular momentum quantum number.

Since the C6H–H2 PES is extremely anisotropic, the surface was further adapted from the analytical fits. Chapman and Green34 found that extremely anisotropic interactions, specifically HC3N–He, lead to expansions that may be slowly convergent and subject to numerical instability. This was the case in our present study where the PES rapidly varies, so a regularized potential was introduced to eliminate the numerical issues. The procedure detailed below, adapted from Wernli et al.14 describes the modification to the C6H–H2 PES. This modification was not necessary for the less expensive calculations on the 2D C6H–He PES.

For low energies the new regularized potential vreg(V(R, θ1, θ2, ϕ)) remained identical to the original potential, but at a threshold energy Va it begins to smoothly saturate along the repulsive curve until reaching a limiting value based on the parameter Vb. A function Sf was used in the saturation area to smoothly switch the original PES to a constant value. This regularized potential allows the Legendre expansion to be performed with computational ease and numerical accuracy and is given by

vreg(V)={V,VVa,Va+(VbVa)Sf(VVaVbVa),V<Vb,Va+(VbVa)(2π)2,VVb,
(3)

where the switching function is defined as

Sf(u)=(2π)2sin[π2sin(πu2)].
(4)

The threshold energy chosen for Va was 300 cm−1, while Vb = 2000 cm−1. These values yield a regularized potential that is extremely accurate up to 500 cm−1, wherein it only diverges from the original PES by ∼4.0 cm−1. The value for Vb prompts the PES to terminate at ∼990 cm−1. A comparison between the original PES and the regularized PES is shown in Figures 5 and 6 for the θ1 fixed values of 0° and 90°.

FIG. 5.

The original C6H–H2 PES (dashed) and the regularized PES (solid) for θ1 = θ2 = 0.

FIG. 5.

The original C6H–H2 PES (dashed) and the regularized PES (solid) for θ1 = θ2 = 0.

Close modal
FIG. 6.

The original C6H–H2 PES (dashed) and the regularized PES (solid) for θ1 = 90°, θ2 = 0°.

FIG. 6.

The original C6H–H2 PES (dashed) and the regularized PES (solid) for θ1 = 90°, θ2 = 0°.

Close modal

One can see in these figures that the regularized PES smoothly saturates therefore eliminating the cusp that is problematic for the lambda terms. The regularization procedure simply amounts to applying a slightly lower ceiling to the PES and with a smoother onset. The ceiling of the regularized potential is sufficiently high enough in energy for our scattering study here, but it should be noted that this modified PES should not be used for high collision energies.

Collisional cross sections were calculated using the quantum non-reactive molecular scattering code MOLSCAT.35 The interaction potential was expanded over the complete set of functions of the angular coordinate. In the case of C6H–H2, the PES was expanded as follows:

V(R,θ1,θ2,ϕ)=l1,l2,lvl1,l2,l(R)Al1,l2,l(θ1,θ2,ϕ),
(5)

where Al1,l2,l are contracted normalized spherical harmonics.36 The angular dependence of C6H was expanded up to order l1,max = 48, while H2 was expanded up to order l2,max = 4.

The C6H–He PES was expanded according to

V(R,θ)=l1=0l1maxvl1(R)Pl1(cosθ),
(6)

where Pl1 are the Legendre polynomials. Radial coefficients up to order l1 = 50 were considered in the expansion.

For collisions between two linear rigid rotators, the fully quantal close-coupling (CC) approach of Green37 was used to determine the integral cross sections. The log-derivative propagator of Manolopoulos38 was used to solve the coupled-channel equations from 4 to 40 a0. The rotational basis included all open channels and several closed channels to secure convergence of the cross sections to within 10% compared to fully converged calculations for a given PES. For H2, only the j2 = 0 and j2 = 1 levels were included in the basis set; the additional inclusion of j2 = 2 and j2 = 3 for the para- and ortho- calculations, respectively, led to computed cross sections within 5%, of the previous more affordable result.

For the C6H–He system, the CC approach for the scattering of a rigid rotor and an atom of Arthurs and Dalgarno39 was used to determine the cross sections. The same log-derivative propagator of Manolopoulos38 was used to solve the coupled-channel equations from 4.8 to 50 a0. The resulting cross sections are converged to within 10%.

De-excitation cross sections for Δj1 = − 1 are presented in Figures 7 and 8 for collisions of C6H with para-H2 and ortho-H2, respectively, while de-excitation cross sections for collisions with He are shown in Figure 9.

FIG. 7.

The largest state-to-state de-excitation cross sections for C6H in collisions with para-H2 (j2 = 0).

FIG. 7.

The largest state-to-state de-excitation cross sections for C6H in collisions with para-H2 (j2 = 0).

Close modal
FIG. 8.

The largest state-to-state de-excitation cross sections for C6H in collisions with ortho-H2 (j2 = 1).

FIG. 8.

The largest state-to-state de-excitation cross sections for C6H in collisions with ortho-H2 (j2 = 1).

Close modal
FIG. 9.

The largest state-to-state de-excitation cross sections for C6H in collisions with He.

FIG. 9.

The largest state-to-state de-excitation cross sections for C6H in collisions with He.

Close modal

The ortho- and para-H2 cross sections are similar in magnitude and in fact there is no significant difference between the two moieties. The energy-dependent de-excitation cross sections for C6H in collision with ortho-H2 (j2 = 1) appear to have a smoother energy dependence than the cross section for collision with para-H2 (j2 = 0). This is a result of the fact that there are many more, and hence overlapping, resonances for ortho-H2 than for para-H2. Hence, the resonance features are mostly washed out for ortho-H2. The 1/v dependence of the state-to-state cross sections can be clearly seen as the values decrease with increasing kinetic energy. Therefore, one can expect a slow temperature variation of the rate coefficients for these low temperature collisions in agreement with Langevin theory for ion-neutral interactions.

The C6H–He cross sections have a similar magnitude to those of H2. The He resonances disappear after 20–30 cm−1, however, while the H2 resonances sustain until ∼80 cm−1.

The propensities for the de-excitation of C6H out of initial level j1 = 5 in collisions with H2 and He are shown for 10 and 300 cm−1 in Figures 10 and 11, respectively. These plots confirm that both the ortho and para-H2 data exhibit similar values for the state-to-state cross sections. Similar trends are generally found for ionic systems (HCO+, N2H+, CN,…). This is because the inelastic process is mostly governed by long range interactions that are weakly anisotropic with respect to H2 rotation. The plots also show that the magnitude of the cross sections decrease with increasing Δj1. Note that for the C6H, the decrease is quite slow because of the small energy spacing between the rotational states and because of the large well depth of in the PESs that will significantly couple states with large Δj1. At low collision energies, the C6H–He cross sections slowly decrease with increasing Δj1 as well. However, at high collision energies, a propensity rules in favor of transitions with even Δj1 due to the near-homonuclear symmetry of the potential energy surface.

FIG. 10.

Propensities for the de-excitation of C6H from initial level j1 = 5 in collisions with H2 and He for EK = 10 cm−1.

FIG. 10.

Propensities for the de-excitation of C6H from initial level j1 = 5 in collisions with H2 and He for EK = 10 cm−1.

Close modal
FIG. 11.

Propensities for the de-excitation of C6H from initial level j1 = 5 in collisions with H2 and He for EK = 300 cm−1.

FIG. 11.

Propensities for the de-excitation of C6H from initial level j1 = 5 in collisions with H2 and He for EK = 300 cm−1.

Close modal

Finally, it is interesting to compare the new C6H–H2 results with the CN–H2 one of Kłos and Lique.10 At low collision energies, both C6H and CN exhibit similar behavior with the largest cross section being the Δj = − 1 transition.

We present the first PESs for the C6H–H2 and C6H–He systems. Fortran codes of the PESs are available upon request. Both ab initio data were computed at the CCSD(T*)-F12b level of theory. 4D and 2D PESs, respectively, were constructed using the IMLS method. Both surfaces similarly exhibit large anisotropies for linear versus T-shape approaches and have deep van der Waals wells. We furthermore performed inelastic scattering calculations on these surfaces in order to describe the collisions of C6H with para-H2, ortho-H2, and He for low internal excitations of C6H. The H2 and He state-to-state de-excitation cross sections are similar in magnitude, although the He cross sections are slightly smaller. At high collisional energies, the C6H–He system exhibits propensity rules in favor of even Δj1 transitions contrarily to the C6H–H2 system. For odd Δj1 transitions, the difference between He and H2 results can be up to an order of magnitude. Such behavior prevents from using He rate coefficients as a template for the H2 ones. Collisional excitation data with H2 are especially important for modeling astrophysical environments containing anions, such as TMC-1 or IRC +10216, and determining the physical and chemical conditions in the regions. In the future, this PES can be used to obtain a complete set of collisional rate coefficients, including highly excited rotational states, for the C6H–H2 and C6H–He systems.

This research was supported by the French National Research Agency (ANR) through a grant to the Anion Cos Chem Project (No. ANR-14-CE33-0013), the CNRS-INSU Programme National de Physique et Chimie du Milieu Interstellaire, and the National Science Foundation (Grant No. CHE-1300945 to R.D.). This work was granted access to the HPC resources of IDRIS under the Allocation No. 2016-047659 made by GENCI (Grand Equipement National de Calcul Intensif).

1.
M. C.
McCarthy
,
C. A.
Gottlieb
,
H.
Gupta
, and
P.
Thaddeus
,
Astrophys. J., Lett.
652
,
L141
(
2006
).
2.
J.
Cernicharo
,
M.
Guélin
,
M.
Agúndez
,
K.
Kawaguchi
,
M.
McCarthy
, and
P.
Thaddeus
,
Astron. Astrophys.
467
,
L37
(
2007
).
3.
S.
Brünken
,
H.
Gupta
,
C. A.
Gottlieb
,
M. C.
McCarthy
, and
P.
Thaddeus
,
Astrophys. J., Lett.
664
,
L43
(
2007
).
4.
M.
Agúndez
,
J.
Cernicharo
,
M.
Guélin
,
C.
Kahane
,
E.
Roueff
,
J.
Kłos
,
F. J.
Aoiz
,
F.
Lique
,
N.
Marcelino
,
J. R.
Goicoechea
 et al,
Astron. Astrophys.
517
,
L2
(
2010
).
5.
P.
Thaddeus
,
C. A.
Gottlieb
,
H.
Gupta
,
S.
Brünken
,
M. C.
McCarthy
,
M.
Agúndez
,
M.
Guélin
, and
J.
Cernicharo
,
Astrophys. J.
677
,
1132
1139
(
2008
).
6.
J.
Cernicharo
,
M.
Guélin
,
M.
Agúndez
,
M. C.
McCarthy
, and
P.
Thaddeus
,
Astrophys. J.
688
,
L83
(
2008
).
7.
V.
Wakelam
,
E.
Herbst
,
J.-C.
Loison
,
I. W. M.
Smith
,
V.
Chandrasekaran
,
B.
Pavone
,
N. G.
Adams
,
M.-C.
Bacchus-Montabonel
,
A.
Bergeat
,
K.
Béroff
 et al,
Astrophys. J., Suppl. Ser.
199
,
21
(
2012
).
8.
F.
Dumouchel
,
A.
Spielfiedel
,
M. L.
Senent
, and
N.
Feautrier
,
Chem. Phys. Lett.
533
,
6
(
2012
).
9.
D.
Hauser
,
S.
Lee
,
F.
Carelli
,
S.
Spieler
,
O.
Lakhmanskaya
,
E. S.
Endres
,
S. S.
Kumar
,
F.
Gianturco
, and
R.
Wester
,
Nat. Phys.
11
,
467
(
2015
).
10.
J.
Kłos
and
F.
Lique
,
Mon. Not. R. Astron. Soc.
418
,
271
(
2011
).
11.
N.
Sakai
,
T.
Sakai
,
Y.
Osamura
, and
S.
Yamamoto
,
Astrophys. J., Lett.
667
,
L65
(
2007
).
12.
H.
Gupta
,
C. A.
Gottlieb
,
M. C.
McCarthy
, and
P.
Thaddeus
,
Astrophys. J.
691
,
1494
(
2009
).
13.
M. A.
Cordiner
,
S. B.
Charnley
,
J. V.
Buckle
,
C.
Walsh
, and
T. J.
Millar
,
Astrophys. J., Lett.
730
,
L18
(
2011
).
14.
M.
Wernli
,
L.
Wiesenfeld
,
A.
Faure
, and
P.
Valiron
,
Astron. Astrophys.
464
,
1147
(
2007
); e-print arXiv:physics/0611258.
15.
T. B.
Adler
,
G.
Knizia
, and
H.-J.
Werner
,
J. Chem. Phys.
127
,
221106
(
2007
).
16.
J. G.
Hill
,
S.
Mazumder
, and
K. A.
Peterson
,
J. Chem. Phys.
132
,
054108
(
2010
).
17.
H.-J.
Werner
,
P. J.
Knowles
,
G.
Knizia
,
F. R.
Manby
, and
M.
Schütz
,
Wiley Interdiscip. Rev.: Comput. Mol. Sci.
2
,
242
(
2012
).
18.
F.
Lique
,
J.
Kłos
, and
M.
Hochlaf
,
Phys. Chem. Chem. Phys.
12
,
15672
(
2010
).
19.
R.
Dawes
,
X.-G.
Wang
,
A. W.
Jasper
, and
T.
Carrington
, Jr.
,
J. Chem. Phys.
133
,
134304
(
2010
).
20.
R.
Dawes
,
X.-G.
Wang
, and
T.
Carrington
, Jr.
,
J. Phys. Chem. A
117
,
7612
(
2013
).
21.
K. B.
Gubbels
,
S. Y.
van de Meerakker
,
G. C.
Groenenboom
,
G.
Meijer
, and
A.
van der Avoird
,
J. Chem. Phys.
136
,
074301
(
2012
).
22.
R. A.
Kendall
,
T. H.
Dunning
, Jr.
, and
R. J.
Harrison
,
J. Chem. Phys.
96
,
6796
(
1992
).
23.
F.
Weigend
and
R.
Ahlrichs
,
Phys. Chem. Chem. Phys.
7
,
3297
(
2005
).
24.
R.
Dawes
,
A. F.
Wagner
, and
D. L.
Thompson
,
J. Phys. Chem. A
113
,
4709
(
2009
).
25.
M.
Majumder
,
S. A.
Ndengue
, and
R.
Dawes
,
Mol. Phys.
114
,
1
(
2015
).
26.
J.
Brown
,
X.-G.
Wang
,
T.
Carrington
, Jr.
,
G.
Grubbs
II
, and
R.
Dawes
,
J. Chem. Phys.
140
,
114303
(
2014
).
27.
X.-G.
Wang
,
T.
Carrington
,
R.
Dawes
, and
A. W.
Jasper
,
J. Mol. Spectrosc.
268
,
53
(
2011
).
28.
J.
Brown
,
X.-G.
Wang
,
R.
Dawes
, and
T.
Carrington
, Jr.
,
J. Chem. Phys.
136
,
134306
(
2012
).
29.
S. A.
Ndengué
,
R.
Dawes
, and
F.
Gatti
,
J. Phys. Chem. A
119
,
7712
(
2015
).
30.
M.
Majumder
,
S. E.
Hegger
,
R.
Dawes
,
S.
Manzhos
,
X.-G.
Wang
,
C.
Tucker
, Jr.
,
J.
Li
, and
H.
Guo
,
Mol. Phys.
113
,
1823
(
2015
).
31.
I. M.
Sobol
,
USSR Comput. Math. Math. Phys.
16
,
236
(
1976
).
32.
U.
Fink
,
T. A.
Wiggins
, and
D. H.
Rank
,
J. Mol. Spectrosc.
18
,
384
(
1965
).
33.
J. V.
Foltz
,
D. H.
Rank
, and
T. A.
Wiggins
,
J. Mol. Spectrosc.
21
,
203
(
1966
).
34.
S.
Chapman
and
S.
Green
,
J. Chem. Phys.
67
,
2317
(
1977
).
35.
J. M.
Hutson
and
S.
Green
, Distributed by Collaborative Computational Project No. 6, Swindon, UK: Engineering and Physics Science Resoruce Council, 1994.
36.
H.
Rabitz
,
J. Chem. Phys.
57
,
1718
(
1972
).
37.
S.
Green
,
J. Chem. Phys.
62
,
2271
(
1975
).
38.
D. E.
Manolopoulos
,
J. Chem. Phys.
85
,
6425
(
1986
).
39.
A. M.
Arthurs
and
A.
Dalgarno
,
Proc. R. Soc. London A
256
,
540
(
1960
).