Polycyclic Aromatic Hydrocarbons (PAHs) in Interstellar Ices: A Computational Study into How the Ice Matrix Influences the Ionic State of PAH Photoproducts

It has been experimentally observed that water–ice-embedded polycyclic aromatic hydrocarbons (PAHs) form radical cations when exposed to vacuum UV irradiation, whereas ammonia-embedded PAHs lead to the formation of radical anions. In this study, we explain this phenomenon by investigating the fundamental electronic differences between water and ammonia, the implications of these differences on the PAH–water and PAH–ammonia interaction, and the possible ionization pathways in these complexes using density functional theory (DFT) computations. In the framework of the Kohn–Sham molecular orbital (MO) theory, we show that the ionic state of the PAH photoproducts results from the degree of occupied–occupied MO mixing between the PAHs and the matrix molecules. When interacting with the PAH, the lone pair-type highest occupied molecular orbital (HOMO) of water has poor orbital overlap and is too low in energy to mix with the filled π-orbitals of the PAH. As the lone-pair HOMO of ammonia is significantly higher in energy and has better overlap with filled π-orbitals of the PAH, the subsequent Pauli repulsion leads to mixed MOs with both PAH and ammonia character. By time-dependent DFT calculations, we demonstrate that the formation of mixed PAH–ammonia MOs opens alternative charge-transfer excitation pathways as now electronic density from ammonia can be transferred to unoccupied PAH levels, yielding anionic PAHs. As this pathway is much less available for water-embedded PAHs, charge transfer mainly occurs from localized PAH MOs to mixed PAH–water virtual levels, leading to cationic PAHs.


INTRODUCTION
The presence of polycyclic aromatic hydrocarbons (PAHs) in the interstellar medium has been proven through the presence of strong infrared (IR) emission bands between 3.3 and 11.3 μm. 1 For quite some time, these IR emission features could not be explained and therefore were named the unidentified IR features. In the 80s, these features were hypothesized to originate from various stretching and bending modes of vibrationally excited hydrocarbon systems, such as PAHs, after excitation by photons from the interstellar radiation field. 1b −d,g,2 Nowadays, the existence of interstellar PAHs is generally accepted, specifically after the recent unambiguous astronomical identification of a number of aromatic species. 3 Organic molecules such as PAHs play an important role in interstellar chemistry. They are formed in the outflows of dying stars, contribute upon fragmentation to the molecular inventory in space, and are expected to freeze out onto cold dust grains. Here, they get embedded in an ice matrix comprising mainly water but that also comprises other species that are expected to form upon the solid-state hydrogenation of atomic precursors, such as ammonia. 4,5 Although water−ice is the main component of interstellar ice analogues, observations with the Spitzer Space Telescope showed that ammonia can be present in water-rich ices, with an abundance typically 10−20 times lower with respect to water. 6 Upon exposure to high-energy irradiation, reactions of PAHs or PAH derivatives with the surrounding matrix are expected to occur. In fact, the vacuum UV (VUV) irradiation and particle bombardment of these interstellar ices produces complex organic molecules that possibly played a vital role in the origin of life. 7 Previous in situ spectroscopic studies have shown that PAHs embedded in a water−ice matrix are readily ionized. 8,9 For example, Gudipati and Allamandola showed that 4-methylpyrene, naphthalene, and quaterrylene in water−ice can be ionized by Lyman-α (10.2 eV) irradiation to their cationic forms. 8 Similarly, Kofman et al. and Bouwman et al. observed the formation of cationic triphenylene and pyrene, respectively, in water−ice at 20 K when exposed to VUV broadband irradiation (120−160 nm) by recording the time-resolved electronic spectra of these PAHs as a function of fluence. 9,10 In the study of Bouwman et al., the authors state that the photolysis of pyrene is the result of direct single-photon ionization of the neutral species. In 2012, however, Cuylle et al. showed that both pyrene and benzo [ghi]perylene embedded in ammonia−ice become anionic when exposed to Lyman-α irradiation. 11 The formation of anionic PAHs becomes more prevalent as the ratio of ammonia to water in mixed ices increases. Cuylle et al. hypothesized that an electron can be transferred to the PAH through ammonia-related photoproducts, and therefore, anionic PAHs are formed. They did, however, not provide experimental evidence for how this process occurs.
The purpose of our study is to achieve an understanding of the underlying physical factors responsible for the observation that upon VUV irradiation, the PAHs in water−ice turn into their radical cations, whereas the PAHs in ammonia−ice rather transform into the corresponding radical anions. Under interstellar and laboratory conditions, PAHs are fully surrounded by matrix molecules. Herein, we aim at disentangling the origin of the observed phenomena from other bulk effects and investigate if, already at the level of the interaction between a single representative PAH and a single matrix molecule, 12 the two opposite tendencies can be recovered and related to a causal physical mechanism.
To understand the effect of the matrix on the photochemical behavior of the PAH in its essence, we need to go back to the fundamental difference in the interaction between one PAH molecule and one water versus one ammonia molecule. We first quantum-chemically analyze the fundamental differences between the electronic structures of water and ammonia in relation to those of three representative PAHs, benzene (Ben, C 6 H 6 ), pyrene (Py, C 16 H 10 ), and benzo[ghi]perylene (BgP, C 22 H 12 ), using Kohn−Sham molecular orbital (MO) theory, as contained in Kohn−Sham density functional theory (DFT). Second, we analyze how these differences affect the nature of the PAH−water versus PAH−ammonia interaction with respect to the MO electronic spectra of the resulting PAH− matrix complexes and the associated potential ionization routes. Third, we compute the actual excitations and their oscillator strength using time-dependent DFT (TDDFT).
Here, we anticipate that we find fundamental differences between water and ammonia molecules that translate into opposite charge-transfer excitations in the corresponding complexes with PAHs, leading preferentially to PAH-to-water excitation (i.e., PAH cation formation in the case of water) and ammonia-to-PAH excitation (i.e., PAH anion formation in the case of ammonia).

COMPUTATIONAL DETAILS
2.1. General Procedure. Calculations were carried out with the Amsterdam Density Functional (ADF) program, version 2013.01. 13 The generalized gradient approximation functional BLYP was employed to find the equilibrium structures of all species using the basis set TZ2P, a small frozen core, and Grimme's DFT-D3 correction with Becke− Johnson (BJ) damping to correct for dispersion effects. 14 The TZ2P basis set is of valence triple-ζ quality and has been augmented by two polarization functions. 15 This basis set is known from previous work to yield relatively small basis-set superposition errors that can be safely neglected; see ref 15b. Scalar relativistic corrections were included using the zerothorder regular approximation (ZORA) to allow for the possible expansion of this work in the future. 16,17 Numerical integration was performed with the Becke grid set to "verygood". 18 All geometry optimizations were performed in the gas phase and without imposing geometric constraints. The equilibrium geometries were confirmed to be in a (local) minimum-energy state using vibrational frequency analyses (i.e., no imaginary frequencies were found).
2.2. Fock Matrix-Symmetrized Fragment Orbitals (FMATSFOs). In the ADF program, the electronic structure of molecular systems, such as our PAH−matrix complexes, can be analyzed in terms of the interaction between two meaningfully chosen fragments and the bonding mechanism between the corresponding fragment MOs. The effective orbital energies of fragments can and, in general, do change once these fragments are at their positions in the final complex and experience each other's electrostatic potential. 19,20 Often, but not always, the effective fragment MO energies are stabilized when they are exposed to the potential of the other fragment. Note that this is the actual starting situation from which orbital interactions take place. These shifted effective fragment orbital energies can be computed in ADF using the keyword FMATSFO. 19 2.3. TDDFT Computations. Possible charge-transfer excitations between the PAHs and matrix molecules were investigated using TDDFT. 21 All optically allowed (i.e., singlet−singlet) vertical excitations were calculated for the PAH−water and the PAH−ammonia complexes. The rangeseparated functional CAMY-B3LYP was used, along with a TZ2P basis set. 22 Relativistic effects were accounted for using ZORA. 16,17 The CAMY-B3LYP functional has shown to provide accurate results for charge-transfer excitations in naphthalene derivatives and is therefore a suitable functional to study the PAHs in this work. 23 For more detailed information on the CAMY-B3LYP functional, we refer to ref 22b.

RESULTS AND DISCUSSION
3.1. MO Analysis of Ammonia, Water, and PAHs. Understanding of the ionization pathways of PAHs in water and ammonia−ice starts with the quantum chemical analysis of the electronic structure of the two different-matrix molecules (H 2 O and NH 3 ). In this analysis, we focus on the highest occupied MOs (HOMOs) and lowest unoccupied MOs (LUMOs) in particular since these levels are most likely involved in electronic excitations. Figure 1 shows the optimized geometries and the HOMO and LUMO orbitals and corresponding energies of water and ammonia.
The HOMO of ammonia (HOMO A ) is at a relatively high energy (−6.1 eV) compared to the HOMO of water (HOMO W , −7.2 eV), resulting in an energy difference between the two HOMOs of 1.1 eV. These HOMOs represent the lone pair (LP) of ammonia and the highest-energy LP of water (the HOMO − 1 W of water, i.e., the other LP, lies at −9.2 eV). In contrast to the HOMOs, the energy difference between the LUMOs of the two matrix molecules is relatively small, only 0.4 eV, where the LUMO of ammonia lies at −0.5 eV and the LUMO of water at −0.9 eV. Qualitatively, similar results are obtained when the water and ammonia molecule are in a medium of other like molecules (the HOMO−LUMO gaps increase by ca. 0.5 eV; see Figure S1).
Next, we compare the relative positioning of the HOMO and LUMO of ammonia and water to the frontier MOs of pyrene (Py) and benzo[ghi]perylene (BgP), which are representative for the PAHs used in the experiments by Cuylle et al. [see Figure 2 for the energy levels of Py, W, and A (left) and BgP, W, and A (right)]. 11 The HOMO and LUMO have πand π*-like characters, respectively, for both Py and BgP. The smallest intermolecular HOMO−LUMO gap with regard to water and ammonia is encircled in red (see Figure 2). The energy gap between the HOMO A and the LUMO PAH is smaller than the energy gap between the HOMO W and the LUMO PAH due to the higher-lying HOMO A . For ammonia, it amounts to 3.7 and 3.5 eV between HOMO A and, respectively, the LUMO Py and LUMO BgP . For water, the values are 4.8 and 4.6 eV between the HOMO W and, respectively, the LUMOs of Py and BgP. The energy differences between the HOMO of the PAHs and the LUMOs of ammonia and water are less pronounced as the LUMOs of the matrix molecules differ less. The smallest energy gap here is found between the HOMO PAH and the LUMO W , being 4.1 and 4.0 eV for Py and BgP, respectively (see also Figure 2).
Let us consider an electron transfer in these model systems assuming that the electron transfer occurs via the smallest energy gap (see encircled gaps in Figure 2). In a PAH−water complex, the electron is excited over the smaller energy gap from the PAH to water (4.1 and 4.0 eV for Py and BgP) than from water to the PAH (4.8 and 4.6 eV for Py and BgP). This suggests that it is more likely to form cationic PAHs in a PAH−water complex than the anionic species. In the case of PAH−ammonia complexes, it is the other way around. The electron is more easily excited over the smaller energy gap from ammonia to the PAH (3.7 and 3.5 eV for Py and BgP, respectively) than from the PAH to ammonia (4.5 and 4.4 eV for Py and BgP, respectively). Here, the suggestion is that it is more likely to form anionic PAHs in a PAH−ammonia complex than cationic PAHs. These findings are consistent with the observations of Cuylle et al. in their in situ experiments. 11 The ground-state energetics of the individual molecules give a good indication of the electronic differences between ammonia and water and how the relative HOMO and LUMO energies of ammonia, water, and the PAHs are responsible for the size of the intermolecular energy gaps. Note, however, that this simplified picture is an incomplete representation as charge-transfer excitations that occur in interstellar ices involve PAHs embedded in them, thus interacting with ammonia and water. The individual groundstate HOMO−LUMO levels do not take into account any change in the MOs as a result of the interaction between the PAHs and the matrix molecules. This coupling or mixing of the PAH levels with the matrix electronic structure is crucial for  understanding the actual factors determining the transition probability for PAHs in water−ice versus that of PAHs in ammonia−ice (vide inf ra).
3.2. MO Interaction Diagrams of PAH−Matrix Complexes. In the previous section, we analyzed the ground-state frontier MOs of the individual PAH, water, and ammonia molecules. For a basic analysis of the MOs of small molecules, we refer the reader to the book "Orbital Interactions in Chemistry" by T. A. Albright et al., see ref 24. Next, we quantify the interactions in the PAH−ammonia and PAH−water complexes to see how the HOMO−LUMO levels change as a result of the interaction. Interacting PAH− matrix complexes were created by placing a single matrix molecule above the π-electronic system of the aromatic rings and subsequently optimizing the complex (see Figure 3). For each of the final six complexes (labeled BenW, BenA, PyW, PyA, BgPW, and BgPA), energy minima were found. In the optimized geometries, the ammonia and water molecules direct the N−H or O−H bond(s) toward the π-system of the PAH, see Figure 3. The Cartesian coordinates of the optimized PAH−matrix complexes can be found in the Supporting Information. Note that we find several local minima for each PAH−matrix combination. These complexes are very similar, both energetically and in terms of their interaction mechanism but only differ in the exact position at which the matrix molecule binds to the π-system of the PAH. Herein, we highlight one of the obtained complexes per PAH−matrix combination that clearly demonstrates the fundamental differences between the PAH−ammonia and PAH−water complexes. For the geometries of the other computed PAH− matrix complex minima, the reader is referred to Figure S3.
The interaction between the PAH and matrix molecule when forming the PAH−matrix complexes is small, with the bond energies ranging between −2.27 and −3.27 kcal mol −1 , which can be considered as a weak N(H)···π or O(H)···π hydrogen bond interaction. In general, the bond energy for the PAH−water complexes (−3.02 to −3.27 kcal mol −1 ) is more stabilizing than that for the PAH−ammonia complexes (−2.27 to −2.66 kcal mol −1 ).
Next, the interaction between the PAH and the matrix molecules in the PAH−matrix complexes is further investigated. The Py−matrix complexes are highlighted here to provide general insights on the PAH−matrix interaction as similar results were found for the Ben−matrix and BgP− matrix complexes. Additional orbital interaction figures for the Ben−matrix, Py−matrix, and BgP−matrix complexes can be found in Figures S4−S9 of the Supporting Information.
The interaction between the occupied MOs of pyrene and the matrix molecules water and ammonia is given in Figure 4. The extent of mixing between the PAH and matrix-molecule orbitals is determined not only by the differences in energy between the orbitals of the PAHs and the matrix molecules but also for an important extent by the orbital overlap. Now, we first inspect how much the extent of mixing is in the case of a water and an ammonia matrix. For the PyW complex ( Figure  4a), there is very little-to-no mixing that occurs in the occupied MOs of the PyW complex, with roughly 98% of the HOMO − 2 PyW originating from the HOMO − 2 Py , 100% of the HOMO − 3 PyW originating from the HOMO − 3 Py , and 98% of the HOMO − 4 PyW originating from the HOMO W . Thus, the occupied π-orbitals of pyrene (HOMO − 2 Py and HOMO − 3 Py ) and the LP of water (HOMO W ) mix minimally and remain mostly localized on the original fragments.
In contrast, the LP of ammonia (HOMO A ) and a π-orbital of pyrene (HOMO − 2 Py ) mix strongly in the PyA complex ( Figure 4b). Two delocalized MOs (HOMO − 2 PyA and HOMO − 3 PyA ) are formed, with roughly 14% of the total MO originating from one monomer while the remaining 85−86% originating from the other monomer (see Figure 4b). It is a manifestation of the Pauli exclusion principle due to the overlap between occupied orbitals of both monomers (see Tables S1 and S2). The HOMO − 3 PyA is a bonding combination, and the HOMO − 2 PyA is an antibonding combination between the HOMO − 2 Py and the HOMO A . This difference in interaction between the PAH−water and PAH−ammonia complexes can be attributed to the HOMO energies of the two matrix molecules and orbital overlap between HOMO of the matrix molecules (ammonia or water) and π-orbitals of the PAH. Now, we return to the rationale behind the stronger PAH− matrix orbital mixing in the case of ammonia compared to water. The amount of mixing between an orbital on one monomer and an orbital on the other monomer is approximately proportional to the overlap between the two orbitals on both monomers S (more precisely, it is proportional to the matrix element between the two orbitals) and inversely proportional to the energy difference between those orbitals, Δε. The absolute energy gap between the HOMO A of ammonia and the HOMO − 2 Py (0.24 eV) is much smaller than the gap between the HOMO W of water and the HOMO − 2 Py (0.84 eV), as can be seen in Figure 4. Furthermore, the spatial orientation of the LP of ammonia relative to the PAH allows for better orbital overlap with the PAH MOs in comparison to water (see also Figures 1 and 3). The reason for this is that the ammonia LP is essentially the nitrogen 2p z (with a slight admixture of hydrogen 1s character), and due to ammonia's trigonal pyramidal geometry, as it points with its N−H bond(s) to the PAH, the LP also points directly to the πsystem of the PAH and, thus, overlaps well. In contrast, the oxygen 2p-type LP is arranged parallel to the PAH surface as the water molecule points with its O−H bond(s) to the PAH (see also Figures 1 and 3). This orbital does not point with a lobe toward the PAH and therefore has a smaller overlap with the π-system. For example, in the PyA complex, the orbital overlap between the LP of ammonia (HOMO A ) and the πorbital of pyrene (HOMO − 2 Py ) is larger (S = 0.0266) than the overlap between the LP of water (HOMO W ) and the πorbital of pyrene (HOMO − 2 Py ) in the PyW complex, which is only 0.0175 (see Table S2). The orbital overlap between the HOMO W and HOMO − 3 Py is even 1 order of magnitude smaller (S = 0.0017). These smaller orbital overlaps explain why the MOs of pyrene do not mix significantly with the HOMO W .
The strong mixing in the PyA complex and the poor mixing in the PyW complex is a clear fundamental difference between these two complexes. As will be shown in the following section, the formation of these mixed PAH−ammonia MOs plays an important role in the available charge-transfer pathways.
3.3. Charge-Transfer Excitations in Interacting PAH− Matrix Complexes. In this section, we analyze the possible charge-transfer excitation pathways in the PAH−matrix complexes (Figure 3) using TDDFT calculations. All chargetransfer excitations below 10.2 eV (Lyman-α radiation) with an oscillator strength greater than 0.2 are considered in the exploration of the possibility of charge transfer upon VUV irradiation. The calculated charge-transfer excitations are reported in Table 1 for the complexes of Figure 3. See Table  S4 and Figure S11 for the charge-transfer excitations of other local minima of these PAH−matrix complexes. In this section, we focus mostly on the Py−matrix complexes as these complexes clearly explain the difference between PAH−water and PAH−ammonia complexes. All relevant MOs involved in the Py−matrix excitations are shown in Figure 5. The MOs of the Ben−matrix and BgP−matrix complexes can be found in Figure S10.
In the PyW complex, charge-transfer excitations can occur from the localized π-orbitals of pyrene (HOMO − 1 PyW or HOMO − 2 PyW ) to a delocalized unoccupied MO consisting of the mixed PAH−water σ*-orbital with a large coefficient on the water molecule (LUMO + 4 PyW ), see Figure 5. This direction of charge-transfer excitations suggests that pyrene donates electronic density to water and therefore can become positively charged upon irradiation.
Interestingly, the direction of charge-transfer excitations in the PyA complex is different from the water-embedded variant PyW. In the PyA complex, charge-transfer excitations are observed from a mixed MO with ammonia and pyrene πorbital character (HOMO − 3 PyA ) to two different unoccupied MOs with a π*-like character localized on the pyrene molecule (LUMO + 2 PyA or LUMO + 3 PyA ), see Figure 5. This direction of charge-transfer excitation leads to a shift of charge density from the occupied mixed MO of the complex toward the unoccupied orbitals of pyrene and indicates a net electronic density flow from ammonia to pyrene, leading to a negatively charged PAH upon irradiation. These results on the dimers of PAH−water and PAH−ammonia show the same trends as the experimental findings by Cuylle et al., who reported the formation of anionic PAHs in an ammonia matrix, while cationic PAHs were formed in a water matrix. 11 The results for the BenW and BgPW complexes are similar to those for the PyW complex: charge-transfer excitations are observed from localized π-orbitals (HOMO − 1 BenW for BenW and HOMO BgPW or HOMO − 1 BgPW for BgPW) toward mixed PAH−water σ*-MOs (LUMO + 2 BenW for BenW and LUMO + 5 BgPW or LUMO + 8 BgPW for BgPW). This suggests the formation of cationic PAH radicals upon the irradiation of PAH−water complexes. For the BenA and BgPA complexes, all observed charge-transfer excitations originate from mixed MOs of the LP of ammonia and the π-orbitals of the aromatic molecule (HOMO BenA or HOMO − 2 BenA for BenA and HOMO − 2 BgPA or HOMO − 3 BgPA for BgPA) toward localized π*-orbitals of the aromatic molecule (LUMO BenA or LUMO + 1 BenA for BenA and LUMO + 2 BgPA , LUMO + 3 BgPA , or LUMO + 4 BgPA for BgPA). Accordingly, this direction of charge transfer suggests the formation of anionic PAH radicals. These results are recovered in all other PAH−ammonia complexes, which all go with predominant ammonia to PAH excitations. In the case of other PAH−water complexes, the situation is not always that consistent, and in some cases, besides PAH-to-water excitations, there also exist water-to-PAH excitations.
In the study by Cuylle et al., it was shown that changing the ratio of water/ammonia in the ice matrix gives rise to a turning point in the ionization direction, that is, raising the ratio water/ ammonia led to the formation of cationic PAHs over anionic PAHs. Our TDDFT findings nicely explain this phenomenon: raising the ratio water/ammonia implies strengthening the mechanism that leads to PAH cation formation while lowering the ratio water/ammonia corresponds to reinforcing the mechanism of PAH anion formation. Note however that our archetypal model systems do not allow for a prediction of the water/ammonia ratio at which the turning point occurs.

CONCLUSIONS
We have quantum-chemically developed a rationale for the experimental observation that water−ice-embedded PAHs become cationic when exposed to high-energy irradiation, Lone pairs are abbreviated as LP. The MOs involved in the charge-transfer excitations (columns 5 and 6) can be found in Figure 5.
whereas ammonia−ice-embedded PAHs lead to the formation of anionic PAH photoproducts. Interestingly, this turns out to be not so much an average medium effect but rather the consequence of discrete quantum chemical interactions that are already effective upon bonding one matrix molecule in PAH−water and PAH−ammonia complexes. Our analyses, based on Kohn−Sham MO theory and TDDFT, reveal that the nitrogen LP-type HOMO of ammonia is optimally oriented toward the aromatic π-system as this matrix molecule binds with its N−H bond(s) to the PAH. This orientation, and the relatively high energy of the ammonia LP HOMO, guarantees good overlap and strong mixing between occupied ammonia and PAH orbitals. Indeed, our TDDFT analyses confirm that excitations from occupied orbitals of the PAH−ammonia complex have a strong component of charge transfer from the ammonia matrix into virtual orbitals that are largely on the PAH. In other words, this generates PAH radical anions.
At variance, the oxygen LP-type HOMO of water is suboptimally oriented, namely, in parallel to the aromatic πsystem as this matrix molecule binds with its O−H bond(s) to the PAH. This orientation, and the somewhat lower energy of the water LP HOMO, impair the buildup of overlap and the mixing between occupied water and PAH orbitals. This yields a situation in which excitations take place from occupied orbitals that remain mainly localized on the PAH into virtuals with water and PAH character, as perceived from TDDFT analyses. This comes down to generating PAH radical cations.

* sı Supporting Information
The Supporting Information is available free of charge at h t t p s : / / p ub s . a c s. o r g / d o i / 1 0 .1 0 2 1 / a c se a r th s p a c echem.1c00433.
Bond energies, orbital interaction diagrams, and MOs of the PAH−matrix complexes and energy decomposition analysis, TDDFT results, and Cartesian coordinates of the optimized PAH−matrix complexes (PDF) Notes