Complete, Theoretical Rovibronic Spectral Characterization of the Carbon Monoxide, Water, and Formaldehyde Cations

New high-level ab initio quartic force field (QFF) methods are explored which provide spectroscopic data for the electronically excited states of the carbon monoxide, water, and formaldehyde cations, sentinel species for expanded, recent cometary spectral analysis. QFFs based on equation-of-motion ionization potential (EOM-IP) with a complete basis set extrapolation and core correlation corrections provide assignment for the fundamental vibrational frequencies of the A˜ 2B1 and B˜ 2A1 states of the formaldehyde cation; only three of these frequencies have experimental assignment available. Rotational constants corresponding to these vibrational excitations are also provided for the first time for all electronically excited states of both of these molecules. EOM-IP-CCSDT/CcC computations support tentative re-assignment of the ν1 and ν3 frequencies of the B˜ 2B2 state of the water cation to approximately 2409.3 cm−1 and 1785.7 cm−1, respectively, due to significant disagreement between experimental assignment and all levels of theory computed herein, as well as work by previous authors. The EOM-IP-CCSDT/CcC QFF achieves agreement to within 12 cm−1 for the fundamental vibrational frequencies of the electronic ground state of the water cation compared to experimental values and to the high-level theoretical benchmarks for variationally-accessible states. Less costly EOM-IP based approaches are also explored using approximate triples coupled cluster methods, as well as electronically excited state QFFs based on EOM-CC3 and the previous (T)+EOM approach. The novel data, including vibrationally corrected rotational constants for all states studied herein, provided by these computations should be useful in clarifying comet evolution or other remote sensing applications in addition to fundamental spectroscopy.


Introduction
Small molecular cations would be potent markers for comet evolution in the Solar System, as these species are believed to be abundant in cometary comae [1]. Understanding the rovibronic transitions of these small molecules would aid in their potential inclusion as markers of photochemical processes in these environments [2]. Previously unidentified lines in the spectra of the Hyakutake and Ikeya-Zhang comets [3] have been attributed to transitions of higher-excited levels of H 2 O + [4], exemplifying the importance of these rovibronic transitions. The "ROSINA ion zoo" observed by the Rosetta spacecraft at comet 67P/Churyumov-Gerasimenko underscores the importance of these species as sentinels for important astrophysical phenomena [5] through in situ mass spectroscopic detection of many small cation species such as CH 3 + , H 2 O + , HCO + and CH 3 OH 2 + . Additionally, characterization of the rovibronic and photoionization spectra of such species would help connect unassigned spectral lines with these transitions, thereby "pulling the weeds" of spectra observed near various Solar System bodies. This would indirectly aid in the detection of new molecules by clarifying unknown spectra of known molecules. Many cometary lines from the near-UV and visible regions of the electromagnetic spectrum are unattributed, representing significant knowledge gaps in the characterization of these bodies and their rovibronic spectra [3,[6][7][8][9][10][11][12][13]. The detection of H 2 O + and CO + in situ by the Rosetta spacecraft as well as the presence of peaks in low resolution double focusing mass spectrometer results corresponding to the mass of H 2 CO + , also during this mission to comet 67P, make an understanding of the full rovibronic spectra of these three molecules desirable for subsequent remote sensing within planetary science beyond fundamental spectroscopic clarification of these three molecular cations.
While some experimental rovibronic spectral data are known for these molecules, many constants are still missing. Modern quantum chemical analysis can provide a complete set of such data. However, no good means for quantum chemically computing anharmonic vibrational parameters of electronically excited states exist. This is necessary to produce a full set of vibrationally-corrected rotational constants for each electronic state. Although this is not all that is necessary for a complete model of electronic spectra, these data are a portion of what is required for a complete UV/Vis model of these electronic excitations.
Looking at the above set of small molecular cations, CO + has already been extensively explored with theoretical approaches [14][15][16] and many of its spectroscopic constants are well characterized experimentally [17,18]. However, because it is a small system similar to H 2 O + and H 2 CO + , it should act as a reasonable benchmark for electronically excited state methods that may then be applied to the latter molecules. Thus, rovibrationally probing the electronically excited states of CO + with the methods outlined in the present work will validate the utility of said methods for small molecular radical cations as well as producing useful and novel spectroscopic data for CO + , such as predictions for the H e distortion constant. TheX 2 Σ + , Ã 2 Π, andB 2 Σ + states of CO + , i.e. the lowest-lying states, are the primary electronic states of astrochemical interest and will be explored in the present study.
H 2 O + and H 2 CO + have several stable, low-lying excited electronic states, showcasing these molecules as prime candidates for rovibronic characterization. These excited electronic states of H 2 O + have been explored quite extensively experimentally [19][20][21][22][23][24][25] and theoretically [26,27]. Similarly, H 2 CO + has also been explored with both experiment and theory [28][29][30][31][32]. However, significant gaps still remain for their spectral classification [2]. H 2 O + has a stable Ã 2 A 1 state and a meta-stableB 2 B 2 state [33]. The Ã 2 A 1 state of H 2 O + is linear and exhibits classic Renner-Teller behavior [34]. Its equilibrium geometry is linear, at which point the species is properly labeled as a Σ state but is degenerate with theX 2 B 1 state at linearity. The equilibrium geometry of theB 2 B 2 state is, however, highly bent to even less than 60 • . TheB 2 B 2 state has complex observed spectra due to coupling with theX and Ã surfaces [20,27] most notably in the form of a conical intersection between the Ã 2 A 1 andB 2 B 2 states [26]. Consequently, these species have been treated in the past with multireference methods [35] as well as diabatic treatments [27]. However, there is considerable practical difficulty with such approaches, where the use of more simplified methods would be desired, if possible, in order to provide additional data. Another complication is that the low intensity transitions of the lower vibrational excitations of the Ã 2 A 1 state [23,35] result in poor experimental characterization of this state. As a result, further analysis of these transitions with other methods is warranted.
H 2 CO + has transitions into the Ã 2 B 1 andB 2 A 1 states in the visible range, which share a similar planar geometry with the ground electronic state [31]. Characterization of these states is relatively more straightforward than the H 2 O + system, although many vibrational modes have not been experimentally characterized. This represents an opportunity for high-level theoretical spectroscopy to provide novel insights.
Even though experimental and theoretical spectroscopic data are available for several of these species, there is significant room for improvement in their characterization. Most experimental data come from photoelectron spectra, where there is difficulty in accurate assignment, and some data are of fairly low resolution [20,28,29,36,37]. For instance, Feller and Davidson [38] argue for a reassignment of the fundamental vibrational frequencies of theB 2 B 2 state of H 2 O + attributed from photoelectron spectra. The highly accurate ab initio adiabatic ionization potentials computed in [38] imply that the fundamental frequencies are likely mislabeled, highlighting the role that theory can play in clarifying the spectra of such species beyond providing corroboration with purely experimental inferences. Additionally, providing theoretical data for low intensity vibrational frequencies would assist in thermodynamic characterization of these species and complete spectral modeling as would be required for the classification of cometary observations. Herein, a simple, but highly accurate [39] adiabatic approach is utilized in an attempt to provide fundamental vibrational frequencies for the electronically excited states for CO + , H 2 O + and H 2 CO + .
Quartic force fields (QFFs) can be used readily and accurately to characterize vibrational frequencies, rotational constants, and other spectroscopic parameters for molecules in their ground electronic states [40]. These have previously been extended to electronically excited states [41,42] using EOM-CC3 [43,44] as well as the recently proposed (T)+EOM approach [39]. The last method has achieved mean absolute differences as low as 1.6 cm −1 for anharmonic frequencies relative to the established benchmark CcCR [45] approach (defined below). Although the (T)+EOM method seems promising, it has not yet been applied beyond the initial set of test cases. Thus, formulating alternative electronically excited state QFF approaches is necessary in order to produce accurate high-level spectroscopic data for these small molecular cations for assisting in their role as markers of cometary phenomenon.
In this work, the application of quartic force fields based on the ionization potential variant of equation-of-motion (EOM-IP) [46] is undertaken. This allows for usage of a closed-shell reference state to describe the open-shell target electronic states. These closedshell references are in many cases more well-behaved than their open-shell counterparts, making such a treatment less prone to computational pitfalls. Additionally, EOM-IP based approaches allow for a broader choice of quantum chemical codes, which account for higher order correlation since many excited state methods are not implemented for openshell references. The application of EOM-IP to these systems is not wholly novel [14]; however, the QFFs explored in the present study provide anharmonic data at a higher level than has been previously available. These QFF approaches may also be useful for astrochemical applications beyond the present study. Beyond EOM-IP, the straightforward equation-of-motion excitation energy (EOM-EE) based QFFs [47] are also employed to provide potentially corroborating data for the molecules of interest if the reference states are sufficiently well-behaved. EOM-EE-CC3 is one of the only methods with higher order correlation widely available for open-shell references [44], and it forms the basis for one set of the QFFs employed here.
Thus, EOM-IP based QFFs together with (T)+EOM are used to explore low-lying excited electronic states of CO + , H 2 CO + and H 2 O + , with the intention of providing what is necessary for a complete theoretical model of these species' spectra in the UV/Vis region. Together with other forms of characterization beyond the scope of the present study (e.g., oscillator strengths), the data presented herein should be fruitful for photophysical applications such as the study of cometary phenomena, especially upon perhelion.

Computational Methods
Spectroscopic data for ground state species are computed with two primary QFF approaches. The first QFF uses explicitly correlated coupled cluster singles, doubles, and perturbative triples energies within the F12 formalism [CCSD(T)-F12b] [48] conjoined to the cc-pVTZ-F12 basis set [49][50][51]. This is dubbed the F12-TZ approach from here on [52,53]. The second uses energies consisting of a three-point aug-cc-pVXZ (X = T, Q, 5) complete basis set extrapolation (CBS: C) [54] with additive core correlation (cC) corrections using the Martin-Taylor basis set [55] and an additional scalar relativistic correction using the Douglas-Kroll (R) formalism [56]: this is known as the CcCR approach [45]. The CBS extrapolation used is given by the following formula, where A, B, and C are the 5Z, QZ, and TZ basis set computations, respectively, and l is the highest angular momentum included in the given basis set: (1) For other QFFs which use a two-point extrapolation, the following formula is used: The Douglas-Kroll correction is computed at the CCSD(T) level using a triple-zeta Douglas-Kroll basis set [57] as the difference between the energy computed with relativity enabled and disabled. For F12-TZ, optimized reference geometries are obtained using CCSD(T)-F12b/cc-pVTZ-F12 energies. CcCR uses reference geometries obtained at the CCSD(T)/augcc-pV5Z level corrected by a Martin-Taylor core correlation correction. All computations for the CcCR and F12-TZ QFFs are performed using MOLPRO 2020 [58]. The QFF is formed by displacing the reference geometry at step sizes of 0.005 Å or radians along a set of symmetry internal coordinates (SICs). The SIC schemes employed for each species are given in Appendix A. The single point energies described above are computed at each of the displaced geometries. Force constants are then generated by performing a least squares fitting of these energies followed by a refitting to the numerically exact minimum geometry. The resulting force constants of the refit are then converted to Cartesian coordinates by the INTDER program [59]. The SPECTRO software package [60] then uses these force constants in second-order vibrational and rotational perturbation theory (VPT2) to produce spectroscopic data [61][62][63].
(Ro)vibrational variational configuration interaction (VCI) calculations are also performed for some QFFs for validation of the VPT2 results. For VCI calculations, the symmetry internal coordinate QFFs are converted to simple internal coordinate QFFs so that α values can be determined [64] for the Morse function of C/O-H and C-O bond stretches. Then, the simple internal coordinate QFFs are analytically converted [65] to simple Morsecosine(-sin) QFFs, which are adopted in VCI calculations using VTET [66] for H 2 O + and MULTIMODE (MM) [67,68] for H 2 CO + . Later for verification, symmetry Morse-cosine(sin) QFFs are directly fitted from the same geometry and energy sets of theB 2 B 2 state of H 2 O + and theB 2 A 1 state of H 2 CO + , and their MM calculations confirm the reliability and consistency of VCI results acquired on simple Morse-cosine(-sin) QFFs. Tests show the VCI values are converged to better than 0.05 cm −1 (H 2 O + ) or 0.2 cm −1 (H 2 CO + ).
For electronically excited states, (T)+EOM/CcCR QFFs [39] are constructed by approximating the energy of a target higher electronic state as a combination of a ground-state CCSD(T) energy and an EOM-CCSD excitation energy to the target state. These (T)+EOM energies are used with the same core correlation and scalar relativistic corrections as in the ground state CcCR QFF to form the (T)+EOM/CcCR approach. The optimized geometry for the (T)+EOM/CcCR QFF is obtained following the scheme in Equation (3), where R (T)+EOM/CcCR is a given geometric parameter for the final geometry. The scheme for the single point energies (E (T)+EOM/CcCR ) is shown in Equation (4). The geometry optimizations for the (T)+EOM/CcCR QFFs are constructed in NWChem [69] due to its ability to optimize molecular geometries for the full (T)+EOM energy for open-shell molecules in a straightforward manner without the need for additional, hand-written wrapper programs. The CCSD(T) energies for these QFF points are calculated using MOLPRO 2020.1 [58]. The EOM-CCSD energies are calculated in PSI4 [70] since this program has restricted open-shell Hartree-Fock EOM methods available [44].
A second electronically excited state approach employed here is based on EOM-EE-CC3. Two approaches are investigated: one uses a two-point CBS extrapolation at the triple-zeta and quadruple-zeta levels ("TQ"); the second uses the same three-point CBS extrapolation ("C") defined for (T)+EOM/CcCR and CcCR. Both of these levels of theory include an additional core correlation correction term ("cC") using the Martin-Taylor basis set and are, thus, named EOM-EE-CC3/TQcC and EOM-EE-CC3/CcC, respectively. These computations are performed in PSI4.
Finally, EOM-IP based QFFs are constructed utilizing two iterative perturbative triples methods: EOM-IP-CCSDT-3 and EOM-IP-CC3 [71]. In addition, QFFs at the EOM-IP-CCSDT level are computed. The latter is employed at the EOM-IP-CCSDT/TQcC and EOM-IP-CCSDT/CcC levels, similar to the EOM-EE-CC3 QFFs. EOM-IP-CCSDT-3 and EOM-IP-CC3 employ the same approaches but with an additional corrective term for higher order correlation at the triple zeta level (T): These approaches are labelled EOM-IP-X/TQcCT and EOM-IP-X/CcCT, where X is CC3 or CCSDT-3. The EOM-IP family of approaches are computed using CFOUR [72]. These QFFs use the aug-cc-pCVTZ basis set [73] rather than the Martin-Taylor basis set used for the other QFFs as a matter of convenience for working with CFOUR. The effect of additive corrections to the EOM-IP based QFFs are examined by providing spectroscopic data from QFFs which use only some components that make up the total EOM-IP-CCSDT-3/CcCT and EOM-IP-CCSDT-3/TQcCT QFFs. These include QFFs using a quadruple-zeta and quintuple-zeta Dunning basis set with no further corrections (EOM-IP-CCSDT-3/QZ and EOM-IP-CCSDT-3/5Z), QFFs using a three-point extrapolation (EOM-IP-CCSDT-3/TQ5) and a three-point extrapolation with a core correlation correction (EOM-IP-CCSDT-3/TQ5+cC). The additive corrections here are chosen partly for parity with other, established QFFs. Further corrections such as quantum electrodynamics effects may ultimately be necessary for spectroscopic accuracy [74,75], however, such treatment is beyond the scope of the present work.
Adiabatic excitation energies (AEEs) for the electronically excited state QFFs are provided by taking the difference of the excited state energy and the corresponding ground state energy (e.g., (T)+EOM/CcCR for the excited state and CcCR for the ground state) while also including the obtained anharmonic zero-point corrections and refitting energies from the QFFs.
EOM-IP based QFFs are computed for ground electronic states of the cationic species as well as the electronically excited states. High-resolution experimental data are available for these ground electronic states and can be used to validate the EOM-IP based QFFs. This ensures that the method's predictions for hitherto unassigned frequencies of the higher states will be reliable. Additionally, ground state type CcCR and F12-TZ QFFs are undertaken for the Ã 2 A 1 H 2 O + state with the goal of providing highly accurate spectroscopic predictions from these reliable QFFs. This is feasible here because the drop in symmetry to C s still results in a different electronic symmetry label for theX and Ã states. The caveat is that such results in uneven treatment of the Renner-Teller pair of this species. While this is not the case for the other species, computation of the available states in this way will allow for internal benchmarking of the pure electronically excited state QFFs when available.

CO+
The results for all three states of CO + studied herein-theX 2 Σ + , Ã 2 Π andB 2 Σ + states-are given in Tables 1-3, respectively. In addition to its role as a marker of cometary phenomena, CO + is a useful test case for this group of open-shell cations and helps to understand the behavior of the electronically excited state QFFs investigated in the present work. Thus, this data will clarify which set of QFF data is reliable for spectroscopic prediction.  Anharmonic vibrational fundamental frequencies from all QFFs match well with experimental values for theX state. No QFF here deviates more than 14 cm −1 from the literature values [17,18] for the C-O stretch fundamental frequency. Notably, the F12-TZ QFF matches with experiment to less than 1.0 cm −1 . Rotational constants also compare reasonably with experimental values, with the closest QFF being CcCR, which places the B 0 constant at 59,148 MHz compared to the experimental value of 59,270.5 MHz. EOM-IP-CCSDT/CcC, the highest-level EOM-IP based QFF employed here, differs from experiment by 16.8 cm −1 for ν 1 and 149.5 MHz for B 0 . The other EOM-IP based approaches are somewhat closer, although this is likely fortuitous and overall they behave similarly to EOM-IP-CCSDT/CcC. All QFFs agree with the experimental value for the D e distortion constant to within 3.0 kHz. Thus, both reference ground-state type QFFs and the EOM-IP family of QFFs, which may be conveniently applied to this ground electronic state for benchmarking, behave reasonably for this system.
The predictions for the H e constant here may also be useful for fully fleshing out line list models for this system. CcCR places this constant at 122.642 mHz, whereas EOM-IP-CCSDT/CcC places it at 133.314 mHz. All other QFFs employed herein fall between these two values. The predictions for the B 1 constant here may also be valuable owing to the high level of theory employed. QFF performance for the electronically excited states of CO + is not so straightforward. The (T)+EOM/CcCR method seems to have significant trouble in treating the Ã 2 Π state of CO + ( Table 2). The fundamental vibrational frequency is nearly 100 cm −1 lower than in experimental data. Investigation of the T 1 diagnostic [76] shows an abnormally high value of 0.035 near the equilibrium geometry of the Ã 2 Π state (approximately 1.25 Å). Figure 1 depicts the T 1 diagnostic of variationally accessible states of CO+ and of neutral CO. This high T 1 diagnostic indicates that theX 2 Σ + state has significant multi-reference character at higher C-O bond lengths, and, thus, lays a poor foundation for (T)+EOM/CcCR's treatment of the higher electronic states. By contrast, the ground electronic state of neutral CO is much better described by a single reference determinant at the equilibrium geometries of the Ã 2 Π state. EOM-EE-CC3 exhibits similar issues as (T)+EOM/CcCR here and has a similar magnitude of error as (T)+EOM/CcCR compared to experiment for the ν 1 frequency. For both of these QFFs, the H e constant, with large negative values on the order of −400 mHz, disagrees significantly with the other QFFs, which place this constant around 30 mHz.
On the other hand, the EOM-IP approaches, which access the Ã state by way of the closed-shell neutral's ground electronic state, behave quite well in comparison to known values, with the high level EOM-IP-CCSDT/CcC QFF placing the result about 22.8 cm −1 higher than experimental data, as seen in Table 2. Most of the QFFs overestimate the experimental frequency to a similar magnitude, but this may be due to a weakness in treating diatomic species with the VPT2 approach used. The problem may also arise from the higher order C-O bond, which is somewhat problematic for QFFs [52,77,78]. This effect is likely exacerbated by the diatomic system, which has no counterbalancing contributions from other, well-behaved types of bonds. The B 0 constant demonstrates similar behavior across QFFs as the fundamental vibrational frequency: significant divergence from experiment for (T)+EOM/CcCR and EOM-EE-CC3 based QFFs are in stark contrast to the reasonable agreement seen for EOM-IP based QFFs. The D e and H e constants for the former two are also significantly different than for the rest of the QFFs. However, EOM-IP does not perform as well for theB 2 Σ + state. EOM-IP QFFs all overestimate this frequency by more than 100 cm −1 . This may indicate potential pitfalls for EOM-IP based QFFs, possibly due to lower quality dynamic correlation in this instance. EOM-EE-CC3 based QFFs, however, perform quite well for this state.
The shorter equilibrium C-O bond length for theB 2 Σ + state results in a lower T 1 diagnostic for the corresponding reference state and, hence, lower multi-reference character. This results in better behavior for the EOM-EE based treatments, as the EOM-EE energies are dependent on accurate ground-state coupled cluster amplitudes for the reference state. Despite this, (T)+EOM still exhibits pathological failure here. The excited electronic states of CO + , thus, offer an important lesson: an intelligent choice of method is likely necessary for studying properties of electronically excited states, and the application of QFFs to electronically excited states will likely depend on having a varied toolbox of methods. The right tool should be applied judiciously for a given system.
Comparison between TQcC QFFs and CcC QFFs is also necessary as the latter group of approaches is likely infeasible for larger systems. The two point T-Q extrapolation does not appear to harm performance much for these species. However, T-Q extrapolation without the quintuple-zeta contribution is less computationally demanding by a fair margin: walltimes for the Ã 2 Π state for EOM-IP-CCSDT are 105 h at the CcC level compared to 22 h at the TQcC level. This reduction of more than 75% produces a difference of less than 1.0 cm −1 in the fundamental frequency between the two methods. In addition, there is minimal difference between EOM-IP-CC3 and EOM-IP-CCSDT-3 QFFs for CO + : 1555.3 cm −1 for the former vs 1554.9 cm −1 for the latter for the fundamental C-O stretch frequency. EOM-IP-CC3 is somewhat less expensive due to the inclusion of less triples terms. Both approximate the highly expensive EOM-IP-CCSDT QFF reasonably well with a difference of less than 2.0 cm −1 for the CcC variants.
Tables 1-3 also provide spectroscopic data from QFFs which use only some pieces of the overall composite EOM-IP-CCSDT-3/CcCT and EOM-IP-CCSDT-3/TQcCT QFFs. First examining basis set size, the ν 1 fundamental frequency forX 2 Σ + CO + changes from 2151.0 cm −1 for a QZ basis set to 2155.6 cm −1 for a 5Z basis set. A three-point T-Q-5 extrapolation continues this trend, with a value of 2160.3 cm −1 . This hierarchical convergence is also seen for the Ã 2 Π state. TheB 2 Σ + state shows a similar trend going from QZ to 5Z, however the T-Q-5 extrapolation lowers the anharmonicity to 1798.9 cm −1 from the 5Z value of 1793.3 cm −1 . The core-correlation correction contributes no more than 10 cm −1 difference to any of these frequencies. The higher-order correlation correction ("T"), results in a large increase of more than 30 cm −1 for ν 1 for theX 2 Σ + state and results in closer agreement with experiment and ground state QFFs compared to not including this correction.
Based on experimental comparisons, CcCR and F12-TZ data are trustworthy for thẽ X 2 Σ + and Ã 2 Π states. For theB 2 Σ + state, EOM-EE-CC3 QFFs seem well-behaved although they fail catastrophically for the Ã 2 Π state. The present CO + data do not clearly favor EOM-IP or EOM-EE based QFFs for future application to electronically excited states which are not variationally accessible and suggest that intelligent method application is important for complex open shell cationic species. Finally, novel high-level data is provided for some spectroscopic constants of these electronic states of CO + in addition to the benchmarking outlined above: H e is predicted for each for each electronic state, which was not available from the existing high resolution electronic studies. Across all non-pathological levels of theory, the H e constant differs substantially across the three electronic states present herein: EOM-IP-CCSDT/CcC places it at 133.314 mHz, 8.936 mHz, and −0.002 mHz for theX 2 Σ + , Ã 2 Π, andB 2 Σ + states, respectively. Thus, models utilizing this constant may provide additional clarity for identification of CO + spectral signatures. B 1 is also provided herein. Referring to EOM-IP-CCSDT/CcC values as qualitatively representative, B 1 is 58558 MHz for theX 2 Σ + state. TheB 2 Σ + state value is 53604 MHz, while there is a more significant shift for the Ã 2 Π state at 46983.5 MHz. Although data from previous theory are available [15] since CO+ is a simple diatomic, the high level anharmonic treatment from established QFF approaches provided here should be useful in fully modeling CO + spectra. Table 4 for F12-TZ, CcCR, and the family of EOM-IP QFFs, agree excellently with the given experimental data for theX 2 B 1 ground state of H 2 O + . Of these, the CcCR QFF is expected to perform the best due to the high-quality basis set extrapolation used as well as the additional corrective terms included. In addition, CcCR being a variational ground state method results in better treatment of dynamic correlation. CcCR places the ν 1 frequency at 3260.2 cm −1 , which agrees to within 1.2 cm −1 of 3259.04 cm −1 from gas phase experiments of Dinelli et al. and Huet et al. [79,80]. CcCR also agrees well the experimental value of 3212.86 cm −1 for ν 2 given in several studies [20,[79][80][81]. Literature values for the ν 3 frequency is reasonably unanimous in assignment at 1408.42 cm −1 [20,36,37]. CcCR predicts this frequency to be 1408.3 cm −1 . Agreement is generally quite close between CcCR and F12-TZ, as is often the case, with a mean absolute difference (MAD) of 2.6 cm −1 for F12-TZ vs. CcCR for the fundamental vibrational frequencies.

Anharmonic vibrational frequencies, provided in
Notably, EOM-IP-CCSDT/CcC achieves quite reasonable agreement with experiment as well, with a mean absolute error (MAE) of 8.8 cm −1 . Although this approach uses costly full triples, the treatment of dynamic correlation for EOM based energies tends to be worse. Hence, EOM-IP-CCSDT/CcC underperforming compared to CcCR is not surprising. This QFF could also potentially be improved by the inclusion of scalar relativistic corrections in future work. Similar to CO + , the EOM-IP-CCSDT-3 and EOM-IP-CC3 approaches approximate the EOM-IP-CCSDT method quite well with EOM-IP-CCSDT-3/CcC almost exactly matching the EOM-IP-CCSDT values. Vibrational fundamentals generated from both VCI and VPT2 calculations for the EOM-IP-CCSDT/CcC QFF agree within a few cm −1 . Thus, the VPT2 calculations here are likely behaving appropriately. Breakdowns of the individual contributions to the EOM-IP-CCSDT-3 QFFs are also given in Table 4, which do not show any anomalous behavior and the overall effect of each of the corrections appears to be small. Increasing basis set quality from QZ, 5Z, to a three-point extrapolation appears to result in systematic convergence of frequencies, as is to be expected. The CcCR harmonic frequencies also agree to within 6.0 cm −1 of the previous theoretical work by Feller and Davidson [38], which uses a similar CCSD(T) approach with additive corrections. EOM-IP-CCSDT/CcC, the most high-level EOM-IP based QFF employed here, corroborates with these harmonic frequencies reasonably well with no more than an 11.0 cm −1 deviation from the CcCR values. The less expensive EOM-IP approaches differ by no more than a few additional cm −1 . F12-TZ reproduces CcCR fairly well. It agrees to within 5.0 cm −1 for all three fundamental frequencies. The anharmonic assignments given here for the ground state of H 2 O + are therefore reliable.
The rotational constants for theX 2 B 1 state of H 2 O + are given in Table 5. CcCR has an MAE of 3415 MHz, or 0.4% error with the assignments from Muller et al. [83]. F12-TZ presents an MAE of 2812 MHz, surprisingly closer to experiment. The core correlation and relativistic effects included in CcCR are often necessary for accurate description of rotational constants [77]. Consequently, F12-TZ being more successful for the ground state of H 2 O + is unusual, which may indicate some minor issues with CcCR here. EOM-IP-CCSDT/CcC has an MAE of 3807 MHz for rotational constants compared to values from Muller et al., which is a respectable performance compared to CcCR at only a 392 MHz difference. Furthermore and once more, there is marginal difference between the very costly EOM-IP-CCSDT/CcC and cheaper approximations. EOM-IP-CCSDT-3/CcCT is actually marginally closer to experiment with an MAE 3755 MHz. There is also a marginal difference here between EOM-IP-CCSDT-3/TQcCT and EOM-IP-CCSDT-3/CcT: an MAD of 643 MHz between the two. Likewise, EOM-IP-CC3/CcCT and EOM-IP-CCSDT-3/CcCT have an MAD of only 725 MHz. Lastly, variations between rotational constants at different levels of theory are most severe in the A rotational constants [77]; EOM-IP-CCSDT-3/CcCT has an MAE of only 643 MHz for the B and C rotational constants compared to experiment. Distortion constants and geometrical parameters for this state are given in Table A1 in Appendix B.
Overall, the ground state QFFs describe theX 2 B 1 state of H 2 O + quite well, as does the family of EOM-IP QFFs. The latter's reasonable performance here, where high quality experimental data is available, is encouraging for its further application to systems and states herein with missing data. Table 6 gives harmonic and anharmonic vibrational frequencies for the Ã 2 A 1 state of H 2 O + . The variational F12-TZ and CcCR approaches agree moderately well with results of Truong et al. [25] in consideration of the wide range of error reported for experiment. CcCR places the ν 1 symmetric stretch at 3232.1 cm −1 compared to the experimental value of 3153 ± 169 cm −1 . The CcCR value for the ν 2 bending frequency, at 995 cm −1 , is slightly outside of the reported margin of error at 903 ± 80 cm −1 . A similar case is seen for the ν 3 anti-symmetric stretch, placed at 3420.7 cm −1 for CcCR compared to 3331 ± 24 cm −1 from Truong et al. The value from Reutt et al. [20] for the ν 1 frequency of 3547 ± 16 cm −1 disagrees both with the values of Truong et al. and theory presented herein. CcCR and F12-TZ agree to within 5.0 cm −1 for ν 1 and ν 2 , with a difference of 41.0 cm −1 for the ν 3 bending motion, with F12-TZ predicting the higher frequency. These methods also compare well with the harmonic frequencies provided by previous authors [38], validating the approach in the present study. Mode Special mention should be given to these methods' uneven treatment of the Renner-Teller pair: uneven occupation of two π orbitals, which are degenerate at linearity, in each of these methods leads to a breaking of said degeneracy. EOM-EE-CC3 QFFs are not investigated for the Ã 2 A 1 state for this reason. However, based on moderately close experimental comparisons given above, the uneven treatment of these orbitals does not seem to have catastrophic consequences. This is probably because the deviation from degeneracy is small.

Ã 2 A 1 H 2 O +
Moving on to the AEEs (i.e., the barrier to linearity from theX 2 B 1 state to the Ã 2 A 1 state), the MRCI work of Brommer et al. [35] calculates a value of 7886 cm −1 , which agrees reasonably with the values given in Table 6 for most levels of theory, except for EOM-IP-CC3. The closest value is 7950.1 cm −1 from CcCR, highlighting the utility of these variationally-accessible QFFs in providing reliable estimates of missing spectroscopic data.
The work herein overall supports the qualitative positions of the band origins of this state given by Truong et al. [25], with potential guidance as to the more exact positions of the band origins. The full set of rotational constants, distortion constants, and geometrical parameters predicted by these and other levels of theory herein are given in Table A2 for the Ã 2 A 1 state. As these are not widely available at high levels of theory [2], these data could prove useful in modeling complete UV/Vis spectra of H 2 O + .
Turning to analysis of the pure electronically excited state QFFs, EOM-IP-CCSDT/CcC matches within reported margins of error for experimental frequencies from Truong et al. [25] for the symmetric stretch as well as the bending motion of Ã 2 A 1 H 2 O + : ν 1 (the symmetric stretch) is 3222.0 cm −1 at the EOM-IP-CCSDT/CcC level compared to the experimental value of 3153 ± 169 cm −1 , while ν 2 (the bend) is placed at 970.2 cm −1 compared to 903 ± 80 cm −1 from experiment. Placement of the ν 3 anti-symmetric stretch is somewhat higher at 3405.3 cm −1 compared to Truong et al.'s value of 3331 ± 24 cm −1 . However, EOM-IP-CCSDT/CcC still has the best performance of all the QFFs employed for comparison to known data, even better than the variationally accessible F12-TZ and CcCR approaches. Application of EOM-IP-CCSDT/CcC to theB 2 B 1 state, which is not variationally accessible, should therefore be trustworthy. Replacing the costly EOM-IP-CCSDT method with CCSDT-3 and CC3 variants seems to harm the description of the ν 2 bending frequency significantly, unlike in other systems where the two are well behaved. EOM-IP-CCSDT-3/CcCT places the ν 2 frequency higher at 1067.1 cm −1 , while EOM-IP-CC3/CcCT places it significantly lower at 690.4 cm −1 . Upon closer examination, the higher placement from EOM-IP-CCSDT-3/CcCT appears to result from the "T" higher-order correlation correction, suggesting that this correction may have issues with the highly anharmonic bending motion. Basis set convergence is also not straightforward for this motion, with EOM-IP-CCSDT-3/QZ placing the motion at 906.9 cm −1 , EOM-IP-CCSDT-3/5Z placing it at 900.5 cm −1 and EOM-IP-CCSDT-3/TQ5 placing it at 912.9 cm −1 .
(T)+EOM/CcCR and variationally-accessible ground state-style CcCR QFFs show reasonably close agreement in the vibrational frequencies for the Ã 2 A 1 state of H 2 O + , given in Table 6. There is only a 0.4 cm −1 difference between the two for the ν 2 frequency, and a difference of 4.5 cm −1 for the ν 1 frequency. The discrepancy in the ν 3 frequency is more significant, with (T)+EOM/CcCR placing this nearly 60 cm −1 higher in frequency at 1056.6 cm −1 . (T)+EOM/CcCR thus seems well behaved for this state, in contrast to its performance for the higher electronic states of CO + .
EOM-IP-CCSDT/CcC has an MAD of 911.18 MHz compared to CcCR for rotational constants (Table A2). Again, this deviation is possibly due to the latter's uneven Renner-Teller treatment. Reasonable agreement is shown between (T)+EOM/CcCR and CcCR with an MAD of 564.2 MHz or 0.2% error. Overall and again, the data for this state shows that EOM-IP-CCSDT/CcC seems to be a valid approach for treating this state, as EOM-IP provides a balanced treatment of the Renner-Teller pair. However, this approach is very costly and not easy to extrapolate to larger systems. Hence, the large difference between CCSDT-3 and CC3 variants in their description of the bending motion is disappointing. This issue may be due to the highly anharmonic nature of this mode, however, and the perturbative triples variants may perform better with more well behaved states.
The ν 2 bending motion of this state is worth some further consideration. This motion is highly anharmonic with a positive anharmonicity of more than 500 cm −1 at the EOM-IP-CCSDT/CcC level with VPT2 (424.1 cm −1 for ω 2 compared to 970.2 cm −1 for ν 2 ). VCI calculations give an anharmonic frequency of 557.535 cm −1 at the EOM-IP-CCSDT/CcC, markedly different from the VPT2 result of 970.2 cm −1 for the ν 2 bend. Thus, although VPT2 and VCI results agree closely for the other two modes, a VPT2 based treatment of the bend for the Ã 2 A 1 state of H 2 O + is likely inadequate. The highly anharmonic nature of this mode, combined with difficulties in properly handling the Renner-Teller effect's breakdown of the Born-Oppenheimer approximation, means that semi-global potential surfaces with multi-reference methods may be necessary for confident characterization of this state. There is considerable difficulty in determining experimental fundamental frequencies for the bending motion due to the vibrational progression and low intensity of the fundamental and lower overtones of this mode [24,35]. Truong et al.
give fairly large margins of error at as much as 169 cm −1 , and values from the previous work by Reutt et al. [20] differ by as much as 400 cm −1 from Truong et al. for the symmetric H-O stretch frequency. Thus, ambiguity in proper assignment of this mode remains, which could possibly be clarified with a more thorough theoretical treatment with semi-global potential surfaces.

3.2.3.B 2 B 2 H 2 O +
Large discrepancies exist between previous experimental assignments and the theoretical results given here for theB 2 B 2 state of H 2 O + in Table 7. However, Feller and Davidson [38] point out that the assignments given by Truong et al. [25] for this state likely result from misassignment of the origin band, which is supported by the discrepancy between all of the QFFs presented here as well as other experimental data both for the fundamental frequencies and the AEE. The AEE given by Reutt et al. [20] of 36,757 cm −1 is about 3000 cm −1 higher than that of all QFFs, a difference of approximately 8%. This lends credence to the assertions of Feller and Davidson of possible misassignment of the origin band. Harmonic frequencies at the CCSD(T)/aug-cc-pV5Z level with additional corrections calculated by Feller and Davidson [38] agree quite closely with those of (T)+EOM/CcCR QFF as well as those for the EOM-EE-CC3/TQcC QFF. The EOM-IP family places these harmonic frequencies relatively higher, e.g. ω 3 is 1985.1 cm −1 for EOM-IP-CCSDT/CcC compared to 1945.1 cm −1 from Feller and Davidson. Frequencies from the other EOM-IP based QFFs closely approximate the EOM-IP-CCSDT/CcC values. Examining the effect of the individual components of the composite EOM-IP QFFs shows that the "T" correction has the largest effect here for the EOM-IP-CCSDT-3 QFFs, lowering the ω 2 harmonic frequency for the bend by more than 10 cm −1 and bringing the value more closely in line with the previous work of Feller and Davidson.
As experimental comparisons are problematic, estimating which set of data and, hence, which method likely produces the most meaningful novel constants is more accurate is difficult. Also, issues with the present adiabatic approach may exist due to the complex potential surface near the well of theB 2 B 2 state [27]. VCI calculations at the EOM-IP-CCSDT/CcC level show some noteworthy discrepancies. The ν 1 symmetric stretch differs 23.8 cm −1 between VPT2 and VCI for the EOM-IP-CCSDT/CcC QFF, which appears to be due to differences in Fermi resonances with the 2ν 2 and 2ν 3 overtone bands. The ν 2 bending motion differs more severely, at 1478.9 cm −1 for VPT2 compared to 1561.7 cm −1 for VCI, a difference of 83 cm −1 . This suggest unusual anharmonicity for this frequency. These VCI results thus suggest some caution in interpreting the VPT2 QFFs for this state. Nevertheless, the present work provides an alternative to a costly and complex diabatic treatment and may be useful in guiding interpretation of future experimental work. Based on adequate performance for CO + and the other electronic states of H 2 O + , the EOM-IP-CCSDT/CcC QFF computed for this state likely offers reasonable estimations of spectroscopic data from the present group of adiabatic approaches.
All QFF approaches, as well as previous theory [38], report a significant lowering of the symmetric and antisymmetric stretching frequencies for theB 2 B 2 state of H 2 O + relative to the values of theX 2 B 1 state. For example, EOM-IP-CCSDT/CcC places the ν 1 and ν 3 frequencies at 2409.3 cm −1 and 1785.7 cm −1 , respectively, for theB 2 B 2 state, compared to 3203.0 cm −1 and 3248.7 cm −1 for theX 2 B 1 state. The geometrical parameters for this state, given in Table A4, showcase an increase in the r 0 (H-O) bond length of 0.125 Å. This is a reasonable explanation for the dramatic decrease in frequency of these motions when combined with the dramatically smaller H-O-H bond angle of 57.662 • and the change in electronic configuration-singly occupying the out-of-plane b 1 π orbital rather than the b 2 π orbital as in the ground state. The shift in frequency for the bending motion, ν 2 , is not as severe but is still marked: 1597.1 cm −1 for theB 2 B 2 state compared to 1402.3 cm −1 for the ground state. These large shifts should provide a clear disambiguation for the presence of this electronic state of H 2 O + in cometary spectra. Rotational constants for all QFFs for thẽ B 2 B 2 state of H 2 O + are given in Table A3, while geometrical parameters and distortion constants are given in Table A4 and should prove useful in modeling spectral features. Lastly, EOM-IP-CCSDT-3/CcCT seems to act as a reasonable, cost-effective approximation to EOM-IP-CCSDT/CcC here. Anharmonic fundamental vibrational frequencies differ no more than 4 cm −1 . EOM-IP-CCSDT/CcC ran with a walltime of 1253 h compared to 108.5 h for EOM-IP-CCSDT-3/CcCT. As a result, the latter may be a prudent choice for application to larger systems where prohibitive scaling will make the former infeasible. EOM-IP-CCSDT-3/TQcCT may also be a viable choice for a further shaving-off of computational cost, as its frequencies differ by no more than 7 cm −1 from EOM-IP-CCSDT/CcC.

3.3.1.X 2 B 2 H 2 CO +
Reasonable agreement is exhibited between the CcCR and F12-TZ QFFs and available experimental frequencies (Table 8) for theX 2 B 2 state of H 2 CO + . EOM-IP-CCSDT-3/TQcCT actually outperforms both ground-state QFFs in comparison to several experimental frequencies: ν 4 is computed to have a value of 1039.2 cm −1 at the EOM-IP-CCSDT-3/TQcCT level, compared to the experimental value of 1036 ± 4 cm −1 . F12-TZ places this frequency at 1013.6 cm −1 while CcCR places it at 1007.9 cm −1 . EOM-IP-CCSDT-3/TQcCT also outperforms the ground state methods for the ν 5 anti-symmetric C-H stretch. Additionally, VCI calculations agree within several cm −1 of the VPT2 results for EOM-IP-CCSDT-3/TQcCT, suggesting the VPT2 QFFs are well behaved for this state. The core correlation correction here makes a large difference in the ν 2 O-C stretching frequency, increasing it by about 6.6 cm −1 . The "T" correction has, again, the largest effect overall and brings many of the fundamental frequencies in closer agreement with experiment, with the only exception being the ν 1 symmetric C-H stretch. EOM-IP-CCSDT-3/TQcCT is the only EOM-IP based QFF employed here, owing to its adequate performance for smaller test cases as well as time and cost constraints. CcCR places the ν 6 frequency nearly 60 cm −1 lower than the experimental value. F12-TZ places this much closer to the experimental value (823.7 cm −1 ) at 825.0 cm −1 . This may be due to known issues with the CcCR approach at handling low-frequency, large amplitude motions [86]. Overall, comparison with experimental data for this state suggests that the EOM-IP-CCSDT-3/TQcCT QFF likely provides a valid prediction of spectroscopic constants for the electronically excited states. Rotational constants and distortion constants are given in Tables A5 and A6, respectively, where there is an MAD of 136 MHz between EOM-IP-CCSDT-3/TQcCT and CcCR in the rotational constants. The vibrationally-corrected rotational constants provided herein are not available experimentally, and the high-level CcCR computations given here are, thus, likely the most reliable estimates produced for these constants.

Ã 2 B 1 H 2 CO +
Frequencies for the Ã 2 B 1 state of H 2 CO + are given in Table 9. Experimental frequencies are available for the O-C stretch and the H-C-O symmetric bending motions, given at 1488 cm −1 [29,30] and 1250 cm −1 [28][29][30], respectively. EOM-IP-CCSDT-3/TQcCT matches these frequencies well: predicting the O-C stretch to be 1489.3 cm −1 , only 1.3 cm −1 higher than the literature value and well within the reported margin of error. The H-C-O symmetric bend is computed to be 1262 cm −1 by this method, only a bit higher than the literature value of 1250 cm −1 . The AEE computed by EOM-IP-CCSDT-3/TQcCT is 26,469.6 cm −1 , 540.6 cm −1 higher than the literature value of 25929 ± 5 cm −1 [28][29][30], a relatively small error compared to CO + and H 2 O + .
VCI and VPT2 results for EOM-IP-CCSDT-3/TQcCT are generally in close agreement, with the largest difference being the ν 2 frequency. VCI places this frequency at 1505.7 cm −1 , 16.4 cm −1 higher than the VPT2 value of 1489.3 cm −1 . The "T" correction here brings both the ν 2 and ν 3 frequencies much closer to the experimental values. The core-correlation correction also has a large effect on several modes, notably the ν 3 H-C-O symmetric bend, which drops by nearly 34 cm −1 when comparing EOM-IP-CCSDT-3/TQ to EOM-IP-CCSDT-3/TQ+cC. Several other motions exhibit large changes due to the core correlation correction, and it is likely important here due to the O-C double bond.
(T)+EOM/CcCR performs reasonably for the ν 2 O-C stretch at 1470.1 cm −1 . However, its H-C-O symmetric bend is 1226.4 cm −1 or 23.6 cm −1 lower than in experiment. The T 1 diagnostic at the CCSD/aug-cc-pV5Z level for the equilibrium geometry of the Ã 2 B 1 state is relatively high at 0.037. This suggests that (T)+EOM/CcCR may have issues here due to the reference state's high multi-reference character, as was seen with CO + .
Experimental frequencies are unavailable for the other vibrational motions. The novel data predicted by EOM-IP-CCSDT-3/TQcCT, which seems well behaved for H 2 CO + , should, therefore, make valuable predictions for these band origins. Several vibrational motions differ significantly from theX 2 B 1 state: the ν 1 , ν 2 , ν 4 , and ν 5 frequencies all differ by more than 150 cm −1 between electronic states. The Ã 2 B 1 state is formed by a π → n transition from the doubly occupied b 1 orbital into the singly occupied b 2 orbital resulting in a lengthening of the C-O bond to 1.338 Å from 1.193 Å from the ground electronic state. This explains the drop in the O-C stretching frequency from 1677.3 cm −1 to 1489.3 cm −1 . As the other geometrical parameters do not change much, the shift in electronic density to the non-bonding orbital, which is primarily centered on the oxygen, is likely responsible for the changes in the other frequencies.
Rotational constants, distortion constants and geometrical parameters for this state predicted herein are given in Tables A7 and A8. The EOM-IP-CCSDT-3/TQcCT values for these parameters will be useful in producing full vibronic spectra of H 2 CO + and enhance understanding of its potential role in cometary phenomena. The A 0 constant increases by 1760 MHz from theX 2 B 1 state value of 267989 MHz while the B 0 changes more significantly, decreasing by 7444 MHz from the ground state value. Several distortion constants also change to an appreciable degree and should help in matching features based on qualitative assessment of spectral signatures.  Table 10 contains the fundamental vibrational frequencies of theB 2 A 1 state of H 2 CO + . Experimental assignment is available for the ν 2 frequency, which is placed at 1304 cm −1 [28][29][30]. There has been controversy with which totally symmetric mode this fundamental corresponds [87]. Although Niu et al. [29] assign this band to the ν 2 C-O stretching motion, the present results support assignment of the experimental 1304 cm −1 band origin to the ν 3 H-C-O symmetric bend, instead. EOM-IP-CCSDT-3/TQcCT places ν 3 at 1326.4 cm −1 and ν 2 at 1360.0 cm −1 . Looking at the AEEs, EOM-IP-CCSDT-3/TQcCT is a mere 46.2 cm −1 (0.11%) off from the literature value of 39,928 cm −1 . By contrast, (T)+EOM/CcCR is off by 955 cm −1 , a much larger percent error of 2.4 % but still better than many of the other molecular states reported above in this work. Therefore, EOM-IP-CCSDT-3/TQcCT is likely the more trustworthy of the two methods employed here, in line with the results for the Ã 2 B 1 state of H 2 CO + .
VCI and VPT2 comparisons for this state at the EOM-IP-CCSDT-3/TQcCT level show some notable discrepancies. The frequencies for ν 2 differ by 14 cm −1 while the ν 3 frequencies differ by 20 cm −1 . These differences may be explained by the presence of a Darling-Dennison type resonance between these two modes which are not present in the VPT2 calculations. ν 1 shows a larger difference between the two approaches, with VPT2 placing the frequency at 2734.8 cm −1 compared to 2653.8 cm −1 for VCI. This difference appears to arise from differences in the Fermi resonance polyad of ν 2 = ν 3 + ν 4 = 2ν 3 = 2ν 4 , with the contributions to the VCI polyad being higher in energy. Looking at the breakdown of individual contributions to the EOM-IP-CCSDT-3/TQcCT composite QFF, the core correlation ("cC") correction does not appear to have as significant a contribution here as with the Ã 2 B 1 state, although the contributions of the "T" correction appear to be of similar magnitude. Since the EOM-IP-CCSDT-3/TQcCT results here seem trustworthy based on experimental comparison for the ν 2 frequency and also taking into consideration its performance for the other states of H 2 CO + , the novel data for the missing fundamental frequencies and rotational constants should prove valuable for modeling the UV/Vis transitions of H 2 CO + in cometary environments. The ν 1 frequency for theB 2 A 1 state is closer to theX 2 B 1 state value at 2733.8 cm −1 compared to the shift seen for the Ã 2 B 1 state. The ν 2 frequency at 1360.0 cm −1 , however, is even lower than the Ã 2 B 1 state value. The r 0 (O-C) bond length at 1.285 Å is shorter than the value for the Ã 2 B 1 state (1.338 Å), marking this behavior as unexpected. TheB 2 A 1 state is formed by a σ → n excitation, where the σ a 1 orbital is centered on the C-O bond, and the non-bonding b 2 orbital is centered on the oxygen. This weakening of the σ bond explains the larger decrease in vibrational frequency. Other frequencies, such as the ν 5 anti-symmetric C-H stretch, are, like ν 1 , generally closer to theX 2 B 1 values.
The feature at 1304 cm −1 , discussed above and attributed here to the H-C-O bend with a value of 1326.4 cm −1 at the EOM-IP-CCSDT-3/TQcCT level, is at least 50 cm −1 away from any band origins from other electronic states, and is likely a useful indicator of thẽ B 2 A 1 in analyzing H 2 CO + spectra. Rotational constants for theB 2 A 1 state are given in Table A9. Distortion constants and geometrical parameters are given in Table A10.

Conclusions
Several high level approaches were explored for evaluating the spectroscopic data of CO + , H 2 O + and H 2 CO + . Newly employed EOM-IP based quartic force fields with higherorder correlation correction and core correlation perform well in comparison to available experimental benchmarks in many cases. EOM-IP-CCSDT-3/TQcCT outperforms groundstate QFFs in its treatment of several modes of theX 2 B 1 state of the formaldehyde cation. New data are provided using this QFF for missing fundamental vibrational frequencies of several states of H 2 CO + , as well as vibrationally-excited rotational constants and distortion constants for the electronically excited states of all molecules studied herein. Additionally, work done in the present study as well as by previous authors supports reassignment of the fundamental vibrational frequencies for theB 2 B 1 state of H 2 O + to approximately 2409.3 cm −1 , 1478.9 cm −1 , and 1785.7 cm −1 .
Many of the EOM-IP based approaches herein, especially the high-level EOM-IP-CCSDT/CcC approach, are infeasible for larger systems due to prohibitive computational cost. The EOM-IP-CCSDT-3/TQcCT QFF seems to be a promising avenue for treating small systems of astrochemical interest such as the small molecular cations in the present study. Even so, this method fails for theB Σ + state of CO + . Additionally, VCI calculations highlight some issues with the VPT2 based treatment for several highly anharmonic motions of the electronically excited states of H 2 O + . These concerns highlight the need for a judicious choice of method when dealing with troublesome electronic states. The previously explored (T)+EOM/CcCR seems to behave poorly for the molecules in the present study, and its usage may be limited to more well-behaved systems. Pathological systems like the Ã 2 A 1 state of H 2 O + or the higher electronic states of CO + may necessitate semi-global potential surfaces for truly accurate treatment but the QFFs presented herein represent a reasonable attempt at making a simple, useful approximation. The new data provided by this paper will be of astrochemical interest for modeling molecular transitions with application to cometary phenomena. Examples of such data include the fundamental vibrational frequencies of the Ã 2 B 1 andB 2 A 1 states of H 2 CO + (e.g., 1489.3 cm −1 and 1360.0 cm −1 for the O-C stretch of both states, respectively) and the tentative reassignments forB 2 B 1 H 2 O + . These new data are necessary for a complete spectral model of electronic excitations, and should be useful for high resolution vibrational studies using ground-based observatories and even the Hubble Space Telescope.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: For the Ã 2 A 1 state, which is linear in its equilibrium configuration, the following scheme is used for 57 symmetry-unique points: For all three states of H 2 CO + , the following scheme consisting of 413 points is used: All three states of CO + use a trivial one-dimensional scheme displacing along the C-O bond length coordinate.