Fingerprinting fragments of fragile interstellar molecules: dissociation chemistry of pyridine and benzonitrile revealed by infrared spectroscopy and theory

The cationic fragmentation products in the dissociative ionization of pyridine and benzonitrile have been studied by infrared action spectroscopy in a cryogenic ion trap instrument at the Free-Electron Lasers for Infrared eXperiments (FELIX) Laboratory. A comparison of the experimental vibrational fingerprints of the dominant cationic fragments with those from quantum chemical calculations revealed a diversity of molecular fragment structures. The loss of HCN/HNC is shown to be the major fragmentation channel for both pyridine and benzonitrile. Using the determined structures of the cationic fragments, potential energy surfaces have been calculated to elucidate the nature of the neutral fragment partner. In the fragmentation chemistry of pyridine, multiple non-cyclic structures are formed, whereas the fragmentation of benzonitrile dominantly leads to the formation of cyclic structures. Among the fragments are linear cyano-(di)acetylene˙+, methylene-cyclopropene˙+ and o- and m-benzyne˙+ structures, the latter possible building blocks in interstellar polycyclic aromatic hydrocarbon (PAH) formation chemistry. Molecular dynamics simulations using density functional based tight binding (MD/DFTB) were performed and used to benchmark and elucidate the different fragmentation pathways based on the experimentally determined structures. The implications of the difference in fragments observed for pyridine and benzonitrile are discussed in an astrochemical context.


Spectroscopic measurements
Experiments were performed in a cryogenic 22-pole ion trap instrument coupled to the FELIX freeelectron lasers described in more detail earlier 1 .Vapors of pyridine and benzonitrile were led into the ion source where the molecules were fragmented using dissociative electron impact ionization with 15 or 50 eV electrons.Under the experimental conditions, we expect a mixture of 80 % canonical and a maximum of 20 % of the alpha-distonic isomer for pyridine 2,3 The benzonitrile cation is determined here to be purely the canonical structure (Supplementary Fig. 2).For most of the fragment masses, a direct electron impact ionization source was used.For the m/z 51, 53 and 78 fragments from the dissociative ionization of pyridine, a Gerlich-type electron impact storage ion source was used 4 .The comparability of the two sources was assessed by measuring the infrared spectrum of m/z 52 from pyridine formed in both sources (Supplementary Fig. 4).From these measurements, similar abundances for the main isomer were measured.More details are shown in Supplementary Fig. 5 and 6.The fragments of interest were mass-selected by a quadrupole mass spectrometer and led into the cryogenic 22-pole ion trap 53 .The ions were trapped and cooled down to around 6-7 K by a pulse of He:Ne (3:1 ratio).The formed ion-Ne complexes were studied using IRPD spectroscopy using the freeelectron laser at the FELIX Laboratory 5 .The ions were irradiated with intense (up to 30 mJ) and tunable infrared radiation provided by FEL-2 operating at 10 Hz in the 550-2200 cm -1 range with a typical 0.5% FWHM bandwidth.When the laser is resonant with a vibrational mode of the ion-Ne complex, the amount of depletion is measured as a function of the wavelength to yield an infrared spectrum.The wavelength is calibrated using a grating spectrum analyser and the ion signal is normalized to the laser pulse energy (E), and number of FELIX pulses (N) to determine the relative intensity (I) according to: I * , with S the observed ion counts and B the baseline ion count number.To increase the match with theoretical calculations, some of the spectra have been also calibrated for the photon flux by multiplying the normalized intensity (I) with the photon energy (hν), were ν is the frequency of the light in wavenumbers.
A saturation depletion measurement 6 has been performed to measure the abundancies of the two isomers in the m/z 52 channel.Multiple laser pulses, resonant with a vibrational mode of one isomer, were used to fully deplete this active isomer.The analysis of this depletion as a function of the deposited energy yielded the abundance of the active isomer (Supplementary Fig. 5 and 6).

Quantum chemical calculations
Structures of relevant cations, inspired from literature, where optimized to their lowest energy structure using density functional theory (DFT) using Gaussian 16 7 .
For most of the molecules the functional/basis set combination B3LYP-GD3/N07D was used to perform the geometry optimization and the vibrational frequency calculations as it showed to be a robust method for various cationic molecules [8][9][10][11] .A typical scaling factor of 0.976 was used to account for anharmonicity.For some of the ions anharmonic treatment of the vibrational modes, using the VPT2 functionality of Gaussian 16, showed a significant improvement and is shown as the calculated spectrum.For the non-linear closed-shell species CH 2 CCCH + with m/z 51, the B2PLYPD3/aug-cc-pVTZ functional and basis set combination was used 12 .
The Renner-Teller and spin-orbit splitting patterns of the bending modes of HC 3 N •+ ( 2 П) and HC 2 NC •+ ( 2 П) with m/z=51 and HC 5 N •+ ( 2 П) with m/z 75 were calculated with an effective Hamiltonian approach similar to that of Steenbakkers et al. 13 For this work no cross-mode Renner-Teller term was assumed so that the splitting of each of the bending modes could be calculated separately.The Hamiltonian for mode  then reduces to: where  , is the vibrational Hamiltonian,  the spin-orbit Hamiltonian and  , the Renner-Teller Hamiltonian.The expressions are given below: Here the nomenclature of Steenbakkers et al. 13 was followed.The spectroscopic parameters  ,  and  , represent the Renner Teller constant, the harmonic vibrational frequency and the spinorbit constant, respectively.For HC 5 N •+ all spectroscopic constants were taken from Gans et al. 14 For HC 3 N •+ and HC 2 NC •+ the harmonic vibrational frequencies and their corresponding Renner-Teller constants were determined based on frequency calculations on the RCCSD(T)-F12a/cc-pVTZ-F12 level of theory.The spin-orbit constant used for HC 3 N •+ was taken from Steenbakkers et al. 13 and was assumed to be equal for HC 2 NC •+ .The intensity ratio of the bending modes of all three species was calculated at the B2PLYPD3/aug-cc-pVTZ level of theory.
The minima structures and transition states for the PESs were calculated at the B3LYP-GD3/N07D level of theory and the corresponding minima were connected using intrinsic reaction coordinate (IRC) calculations.All energies have been corrected for the zero-point vibrational energy.

Molecular dynamic simulations
In order to explore the observed dissociation pathways and kinetics, extensive on-the-fly Born Oppenheimer molecular dynamics (MD) simulations with the electronic structure computed at the SCC-DFTB level of theory -hereafter quoted MD/DFTB -were run using the deMonNano code. 15We showed that using the original mio set of parameters 16 led to discrepancies with experiments regarding the H loss vs. the C 2 H 2 loss for the dissociation of PAHs 17 , leading us to scale the C-H atomic integrals by a scaling factor (0.95) in order to decrease the strength of the C-H bond with respect to that of the C-C bond.Using such a Hamiltonian led to better agreement with experimental results 18,19 Similarly, in this work, we consistently adjusted the parameters for the N-H bond as explained in the next paragraph.
Two sets of parameters, Set1 and Set2, were used.The default atomic integral values  ℎ  and   , where  are the atomic orbitals of X=C,N,H and ℎ the monoelectronic Hamiltonian, were scaled by 0.95 for X=C,N,H in the Set2 parameters set, consistently with our previous work 19 .However, the energy difference between cationic benzonitrile and the corresponding alpha-distonic ion was found too high compared to DFT results (84 kJ/mol vs. 49 kJ/mol 20 ).In order to recover the DFT energy difference,  ℎ  and   were scaled by 0.95 for X=C,H and 0.98 for X=N.Besides, charge model3 (CM3) charges were used instead of the original Mulliken charges in the SCC-DFTB Hamiltonian to improve the description of the polarity of the bonds as in previous works (q cm3 =q Mull + d XY with d NH =0.120 21 , d CH =0.09 22 , and dN C =0.08).With this final set of parameters named hereafter Set1, the energy difference between cationic benzonitrile and the corresponding alpha-distonic ion was computed to be 41 kJ/mol.Similarly, Set1 parameters provide a better description of the energy difference between cationic pyridine and the corresponding alpha-distonic ion (+4.4 kJ/mol) whereas with Set2, an energy difference of 44 kJ/mol is obtained, our reference being the DFT value (with ZPE corrections) of -5.8 kJ/mol 2 .However, with both sets of parameters HNC becomes too low in energy compared to HCN (with Set2: E (HNC-HCN) = +9.9kJ/mol while with Set1 E (HNC-HCN) = -56 kJ/mol), the reference value being +62 kJ/mol 23 .This illustrates the limits of transferability of the SCC-DFTB parameters when different types of bonds are involved.However, for these systems accurate energetics of the minima and the corresponding barriers of the entire ion are more important in order to simulate the fragmentation correctly.
Several hundreds to 1800 simulations of 300 ps to 1 ns with initial random atomic velocities were ran in the microcanonical (NVE) ensemble for several internal energies.The SCC-DFTB energy and its gradient were computed every 0.1 fs.The lowest internal energy for which fragmentation was observed for cationic benzonitrile in these timescales was determined to be 8.0 eV.Several energy values were investigated up to 8.8 eV.The dissociation of cationic pyridine was also investigated using Supplementary Fig. 18: Kinetics for the dissociation of canonical pyridine Supplementary Tab.1: Stoichiometry of the fragments resulting from the dissociation of benzonitrile •+ (pool of 1800 simulations) at the end of the MD/DFTB simulations for the two sets of parameters at different internal energies of the respective cation.
•+ using Set1 parameters (1800 simulations at 6.74 eV internal energy of the respective cation).