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 H_{2} molecules. In this paper, we focus on the collisional excitation of the first observed molecular anion, C_{6}H^{−}, by He and H_{2}. 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 C_{6}H^{−}–H_{2} and C_{6}H^{−}–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 H_{2} cross sections. Conversely, no significant differences exist between the collisions of C_{6}H^{−} with the two species of H_{2} (*para*- and *ortho*-H_{2}). We expect that these new data will help in accurately determining the abundance of the C_{6}H^{−} anions in space.

## I. INTRODUCTION

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 C_{6}H^{−}, 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, C_{4}H^{−}, C_{8}H^{−}, CN^{−}, C_{3}N^{−}, and C_{5}N^{−}, 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 C_{2}H^{−} 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^{−}–H_{2} collisional system computed by Kłos and Lique.^{10} The CN^{−}–H_{2} rate coefficients were obtained from a new reliable *ab initio* potential energy surface of the CN^{−}–H_{2} complex and the quantum scattering calculations were performed to investigate rotational energy transfer in collisions of the CN^{−} molecule with both *para*-H_{2} (*j* = 0) and *ortho*-H_{2} (*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, C_{6}H^{−} is of particular interest. C_{6}H^{−} has been found in several star-forming regions^{1,11–13} and is one of the most abundant anionic species. He and H_{2} are generally the dominant colliders in interstellar and circumstellar media, the predominant one being H_{2}.

Our aim is then to compute collisional rate coefficients for the C_{6}H^{−}–He and C_{6}H^{−}–H_{2} 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 H_{2} 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 C_{6}H^{−}–H_{2} and C_{6}H^{−}–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.

## II. COMPUTATIONAL DETAILS

### A. Electronic structure calculations

For the C_{6}H^{−}–H_{2} system, both monomers were held rigid in collinear arrangements with the internuclear bond distance for H_{2} fixed at *r*_{HH} = 0.766 65 Å (the vibrationally averaged bond distance for *para*- hydrogen, *j* = 0). The geometry of C_{6}H^{−} 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.

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 others^{21} 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 T_{1}-diagnostic was monitored and found to be roughly 0.017 for all geometries in the high-level data set. For the C_{6}H^{−}–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 C_{6}H^{−}–H_{2} system out to 22 Å separation, some small numerical discrepancies were noted beyond 17 Å for C_{6}H^{−}–He. Tests with various other auxiliary bases did not improve the behavior. Thus, to determine an analytic long-range for the 2D C_{6}H^{−}–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).

Method . | Well 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 |

Method . | Well 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 |

### B. Analytic representation

The four dimensional (4D) intermolecular C_{6}H^{−}–H_{2} potential includes 4691 symmetry-unique high-level *ab initio* energy data and is represented analytically by the interpolating moving least squares (IMLS) method^{19,24,25} using a weight function to interpolate between local fitting basis expansions. The C_{6}H^{−}–H_{2} system is similar to the previously studied CO_{2}–CS_{2},^{26} (NNO)_{2},^{19,27} (OCS)_{2},^{28} and (CO)_{2}^{20,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 H_{2} 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 H_{2} 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 CO_{2}–CS_{2}, (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 C_{6}H^{−}–H_{2} 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 sequence^{31} 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.

The long length of the C_{6}H^{−} fragment introduces anisotropy with respect to approach of the H_{2} fragment either towards the ends or side of C_{6}H^{−}. In addition, there is a remarkably strong anisotropy with respect to the orientation of the H_{2} fragment (see Figures 2 and 3). End-on approach (*θ*_{2} fixed at 0 or *π*) of H_{2} includes geometries such as the completely collinear global minimum (−712.1 cm^{−1}) with H_{2} interacting with the C-end of the C_{6}H^{−} fragment. The side-on approach of H_{2} (*θ*_{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 C_{6}H^{−}–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 C_{6}H^{−}–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°.

## III. INELASTIC CROSS SECTIONS

### A. Methods

We consider collisions of C_{6}H^{−} with *para*-H_{2} and *ortho*-H_{2} such as

where *j*_{1} and *j*_{2} designate the rotational levels of C_{6}H^{−} and H_{2}, respectively.

The calculations are restricted to low/moderate collisional energies (*E _{c}* < 500 cm

^{−1}). The rotational energy levels of the C

_{6}H

^{−}and H

_{2}molecules were computed using the experimental spectroscopic constants of $BeC6H\u2212=0.045\u2009927$,

*α*

^{C6H−}= 3.3356 × 10

^{−5}, $DeC6H\u2212=1.079\xd710\u22129$, and $BeH2=60.853$,

*α*

^{H2}= 3.062, $DeH2=0.0471$.

^{1,32,33}

We also consider collisions of the C_{6}H^{−} anion with helium,

where *j* is the orbital angular momentum quantum number.

Since the C_{6}H^{−}–H_{2} PES is extremely anisotropic, the surface was further adapted from the analytical fits. Chapman and Green^{34} found that extremely anisotropic interactions, specifically HC_{3}N–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 C_{6}H^{−}–H_{2} PES. This modification was not necessary for the less expensive calculations on the 2D C_{6}H^{−}–He PES.

For low energies the new regularized potential *v*_{reg}(*V*(*R*, *θ*_{1}, *θ*_{2}, *ϕ*)) remained identical to the original potential, but at a threshold energy *V*_{a} it begins to smoothly saturate along the repulsive curve until reaching a limiting value based on the parameter *V*_{b}. A function *S _{f}* 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

where the switching function is defined as

The threshold energy chosen for *V*_{a} was 300 cm^{−1}, while *V*_{b} = 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 *V*_{b} 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°.

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 C_{6}H^{−}–H_{2}, the PES was expanded as follows:

where *A*_{l1,l2,l} are contracted normalized spherical harmonics.^{36} The angular dependence of C_{6}H^{−} was expanded up to order *l*_{1,max} = 48, while H_{2} was expanded up to order *l*_{2,max} = 4.

The C_{6}H^{−}–He PES was expanded according to

where *P*_{l1} are the Legendre polynomials. Radial coefficients up to order *l*_{1} = 50 were considered in the expansion.

For collisions between two linear rigid rotators, the fully quantal close-coupling (CC) approach of Green^{37} was used to determine the integral cross sections. The log-derivative propagator of Manolopoulos^{38} was used to solve the coupled-channel equations from 4 to 40 a_{0}. 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 H_{2}, only the *j*_{2} = 0 and *j*_{2} = 1 levels were included in the basis set; the additional inclusion of *j*_{2} = 2 and *j*_{2} = 3 for the *para*- and *ortho*- calculations, respectively, led to computed cross sections within 5%, of the previous more affordable result.

For the C_{6}H^{−}–He system, the CC approach for the scattering of a rigid rotor and an atom of Arthurs and Dalgarno^{39} was used to determine the cross sections. The same log-derivative propagator of Manolopoulos^{38} was used to solve the coupled-channel equations from 4.8 to 50 a_{0}. The resulting cross sections are converged to within 10%.

### B. Results

De-excitation cross sections for Δ*j*_{1} = − 1 are presented in Figures 7 and 8 for collisions of C_{6}H^{−} with *para*-H_{2} and *ortho*-H_{2}, respectively, while de-excitation cross sections for collisions with He are shown in Figure 9.

The *ortho*- and *para*-H_{2} 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 C_{6}H^{−} in collision with *ortho*-H_{2} (*j*_{2} = 1) appear to have a smoother energy dependence than the cross section for collision with *para*-H_{2} (*j*_{2} = 0). This is a result of the fact that there are many more, and hence overlapping, resonances for *ortho*-H_{2} than for *para*-H_{2}. Hence, the resonance features are mostly washed out for *ortho*-H_{2}. 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 C_{6}H^{−}–He cross sections have a similar magnitude to those of H_{2}. The He resonances disappear after 20–30 cm^{−1}, however, while the H_{2} resonances sustain until ∼80 cm^{−1}.

The propensities for the de-excitation of C_{6}H^{−} out of initial level *j*_{1} = 5 in collisions with H_{2} 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*-H_{2} data exhibit similar values for the state-to-state cross sections. Similar trends are generally found for ionic systems (HCO^{+}, N_{2}H^{+}, CN^{−},…). This is because the inelastic process is mostly governed by long range interactions that are weakly anisotropic with respect to H_{2} rotation. The plots also show that the magnitude of the cross sections decrease with increasing Δ*j*_{1}. Note that for the C_{6}H^{−}, 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 Δ*j*_{1}. At low collision energies, the C_{6}H^{−}–He cross sections slowly decrease with increasing Δ*j*_{1} as well. However, at high collision energies, a propensity rules in favor of transitions with even Δ*j*_{1} due to the near-homonuclear symmetry of the potential energy surface.

Finally, it is interesting to compare the new C_{6}H^{−}–H_{2} results with the CN^{−}–H_{2} one of Kłos and Lique.^{10} At low collision energies, both C_{6}H^{−} and CN^{−} exhibit similar behavior with the largest cross section being the Δ*j* = − 1 transition.

## IV. CONCLUSIONS

We present the first PESs for the C_{6}H^{−}–H_{2} and C_{6}H^{−}–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 C_{6}H^{−} with *para*-H_{2}, *ortho*-H_{2}, and He for low internal excitations of C_{6}H^{−}. The H_{2} 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 C_{6}H^{−}–He system exhibits propensity rules in favor of even Δ*j*_{1} transitions contrarily to the C_{6}H^{−}–H_{2} system. For odd Δ*j*_{1} transitions, the difference between He and H_{2} results can be up to an order of magnitude. Such behavior prevents from using He rate coefficients as a template for the H_{2} ones. Collisional excitation data with H_{2} 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 C_{6}H^{−}–H_{2} and C_{6}H^{−}–He systems.

## Acknowledgments

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).