11institutetext: Departamento de Física and IUdEA, Universidad de La Laguna (ULL), E-38200 Tenerife, Spain
11email: [email protected]
22institutetext: De Vinci Higher Education, De Vinci Research Center,92 916 Paris, France 33institutetext: Institut Supérieur des Arts Multimédia de la Manouba, Université de la Manouba, 2010 la Manouba, Tunisia 44institutetext: Faculté des Sciences de Tunis, Département de Physique, (LPMC), Université de Tunis El Manar, 2092 Tunis, Tunisia 55institutetext: Sorbonne Université, Observatoire de Paris, Université PSL, CNRS, LERMA, F-92195 Meudon, France 66institutetext: Instituto de Astrofísica de Canarias, C/ Via Láctea s/n, E-38205 La Laguna, Spain 77institutetext: Departamento de Astrofísica, Universidad de La Laguna (ULL), E-38206 La Laguna, Spain 88institutetext: Institut de Quimica Computacional i Catalisi (IQCC), Departament de Quimica, Universitat de Girona, Girona E-17003, Catalonia, Spain 99institutetext: CIGUS CITIC - Department of Nautical Sciences and Marine Engineering, University of A Coruña, Paseo de Ronda 51, E-15011 A Coruña, Spain 1010institutetext: Consejo Superior de Investigaciones Científicas (CSIC), Spain 1111institutetext: Departament de Física Quàntica i Astrofísica (FQA), Universitat de Barcelona (UB), c. Martí i Franqués, 1, 08028 Barcelona, Spain 1212institutetext: Institut de Ciències del Cosmos (ICCUB), Universitat de Barcelona (UB), c. Martí i Franqués, 1, 08028 Barcelona, Spain 1313institutetext: Department of Physics, College of Khurma University, Taif University, P.O. Box 11099, Taif 21944, Saudi Arabia

Peculiarities in the infrared emission of PAH-C60 adducts

R. Barzaga    B. Kerkeni,,,    D. A. García-Hernández,    X. Ribas    T. Pelachs    M. Manteiga    A. Manchado,,    M. A. Gómez-Muñoz,    T. Huertas-Roldán,    G. Ouerfelli
(Received ??; Accepted ??)

The coexistence of polycyclic aromatic hydrocarbons (PAHs) and the C60 fullerene in different astrophysical environments can give rise to the formation of new complex species denoted as PAH-C60 adducts, which may contribute to the infrared (IR) emission observed. These PAH-C60 adducts have been previously reported experimentally due to the high reactivity between PAHs and C60. From the astrophysical point of view, however, they have not been considered in detail yet. Here we have performed a combined experimental and theoretical study in order to characterize the IR spectra of PAH-C60 adducts, including multiple adducts. By using new advanced experimental techniques, we have been able to synthesize some specific PAH-C60 adduct isomers, and measured their IR spectra. These experimental data are used to correct their harmonic scaled spectra, as obtained from quantum-chemistry calculations performed at the density functional theory (DFT) level under the B3LYP-GD3/6-31+G(d) approach. This way, we simulate the IR (\sim3-25 μ\mum) spectra of multiple PAH-C60 adducts, composed by a different number of PAH units: mostly one or two units. In addition, the chemical kinetics data available in the literature are used to tentatively estimate the possible order of magnitude of the abundances of these PAH-C60 adducts using the available observational data. Essentially, our results reveal a possible strong modification of the IR spectra when astronomically estimated abundances are considered. Several spectral peculiarities are observed, such as a broad \sim3.4-3.6 μ\mum feature, and important modifications in the 6-10 and 12-16 μ\mum spectral regions together with contributions to the C60 features at 7.0 and 18.9 μ\mum. Interestingly, these PAH-C60 adducts lack aliphatic CH bonds, but they display IR features around 3.4 μ\mum, challenging previous interpretations of this astronomical feature.

Key Words.:
Stars: carbon – Infrared: stars – Abundances – Astrochemistry

Introduction

Fullerenes, like C60, C+60{}_{60}^{+} and C70, are the biggest molecules detected in space so far; either by their infrared (IR) emission bands or their visible absorption bands (Cami et al., 2010; Campbell et al., 2015). Specially neutral C60, which has been detected via its IR transitions mainly in circumstellar envelopes of evolved stars such as young planetary nebulae (PNe) (Cami et al., 2010; García-Hernández et al., 2011b), but also in diverse astrophysical environments like R Coronae Borealis (RCB) stars (García-Hernández et al., 2011a), reflection nebulae (Sellgren et al., 2010), young stellar objects (Roberts et al., 2012), the diffuse interstellar medium (ISM, Berné et al., 2017) and more recently, in regions related to proto-planetary disks (Iglesias-Groth, 2019; Arun et al., 2023). The spectral signature of C60 is recognizable by its four strongest mid-IR emission bands at \sim7.0, 8.5, 17.4 and 18.9 μ\mum. According to these four bands, thermal models can be used to predict the C60 temperature in different astrophysical objects; specially in the case of circumstellar envelopes of evolved stars where it has been estimated that an average temperature of \sim300 K is representative of the C60 IR emission (Cami et al., 2011; García-Hernández et al., 2012b). In many of these astrophysical environments C60 is accompanied by other organic species, which potentially could react with C60, forming more complex molecular systems or C60 derivatives. In fact, Barzaga et al. (2023a) show that the variation of the C60 17.4μ\mum/18.9μ\mum band ratio observed in different C60-rich circumstellar envelopes (PNe and RCB stars) may imply that there are other species like metallofullerenes, among possibly others, contributing at these wavelengths. This reinforces the idea of the possible presence of C60 derivatives in space environments. One of the most obvious candidate molecules to form C60 derivatives are polycyclic aromatic hydrocarbons (PAHs), which are detected in conjunction with C60 in very different fullerene-rich sources (e.g. García-Hernández et al., 2010; Sellgren et al., 2010; Otsuka et al., 2013; Arun et al., 2023).

The PAHs have been largely assumed to be responsible for the discrete unidentified IR (UIR) emission bands at \sim3.3, 6.2, 7.7, 8.6 and 11.2 μ\mum widely observed in the Universe; from our Solar System to old stars and very distant galaxies, among others (see e.g. Peeters et al., 2021, 2024; Chown et al., 2024, for a recent review and the comprehensive study of the Orion Bar). According to astrophysical models, the UIR emission bands can be attributed to large dimension PAHs composed by \sim25-100 C atoms (e.g. Allamandola et al., 1989; Puget & Léger, 1989; Chown et al., 2024). Indeed, a family of large PAHs (neutrals and ions) have to be randomly combined in order to reproduce the UIR bands (e.g. in terms of the broadness of the features, Rosenberg et al., 2014)111Note that the PAH model has been proved to reproduce the UIR emission bands, but it is not the only explanation; e.g. a family of PAH model spectra can also, in some sources, fit the spectra of several other species, both organic and inorganic ones (see e.g. Zhang & Kwok, 2015).. However, such large PAHs have not been unambiguously detected yet in space; neither in the ISM nor any other astrophysical environment. In contrast, small cyano-PAHs like cyanonaphtalene (McGuire et al., 2021) as well as the pure PAH indene (Burkhardt et al., 2021; Cernicharo et al., 2021) have been recently detected for the first time in space, towards the cold dark cloud TMC-1. These detections have been only possible by radioastronomy, combined with microwave spectroscopy measurements, being the first convincing evidences of the existence of small PAHs in astronomical environments222Very recently, cyano derivatives of the PAHs acenaphthylene (C12H8), pyrene (C16H10) and coronene (C24H12) have been radioastronomically detected in TMC-1, being the largest PAH derivatives presently detected in space (Cernicharo et al., 2024; Wenzel et al., 2025b, a)..

The C60 and PAH molecules can easily react to form new hybrid species denoted as C60 adducts333The definition of C60 adducts implies the binding of another organic molecule and/or functional group, e.g. alkyl chains, alcohols, among others., which are obtained by well-known organic chemistry experimental procedures (e.g. Sarova & Berberan-Santos, 2004; García-Hernández et al., 2013; Dunk et al., 2013; Barham et al., 2018). In particular, these methods have been employed to synthesize C60 adducts with acenes (anthracene, tetracene, pentacene), as well as with indene. The formation reaction of such C60 adducts occurs smoothly, yielding single to multiple PAHs attached to the C60 cage (see e.g. He et al., 2010; Cataldo et al., 2014). This provokes that a mixture of regioisomers444Regioisomers are molecules with the same chemical composition and functionality, but with different spatial arrangement. is produced, which is a main drawback, especially for their IR spectra measurements. However, a recent study reports the use of a supramolecular mask strategy in order to perform regioselective synthesis of anthracene and pentacene C60 bis-adducts (Pujals et al., 2022). This novel strategy allows to pinpoint the selection of the desired PAH isomers, even for multiple addend PAH units to the carbon cage (see Pujals et al., 2022, for more details). In short, all previous experimental studies clearly show the trend of C60 to easily form adducts under the presence of these small PAHs, something that it is very likely to also occur in astrophysical environments under UV-shielding conditions such as evolved stars, where both types of organic ingredients may coexist (see Sect.Astrophysical relevance).

Herein, we present a quantum-chemical study, supported by novel experimental data, of multiple C60 adducts with the small PAHs indene, indenyl, anthracene, tetracene and pentacene. Mono-adducts models are reported for all these PAHs, while the bis- and tris-adducts models are only included for the small PAHs indene, indenyl and anthracene. Their IR spectra have been simulated in order to capture the structural effects of the multiple addend PAH units and their relation with the IR spectral features observed. For this purpose, reliable scaling factors have been previously obtained and validated using the advanced experimental IR spectra measurements available for these C60 adducts. Several spectral peculiarities are seen in the simulated IR spectra of multiple PAH-C60 adducts, which are made publicly available to the astrophysical community, especially to users of the James Webb Space Telescope (JWST), for potential comparisons with their astronomical observations.

Methods

Quantum-chemistry models of PAH-C60 adducts require an accurate description of the numerous electrons system composing these molecules. The approach selected in order to carry out these quantum-chemical calculations is mainly determined by the number of electrons in the system. The density functional theory (DFT) is the natural choice to model large molecules containing more than 100 electrons. The PAH-C60 adducts presented here contain up to \sim650 electrons, which corresponds to the maximum PAH units of anthracene added to C60. The large system size makes thus unfeasible the treatment of our models with more accurate quantum-chemical approaches other than DFT; e.g. ab-initio methods.

Computational details

The Gaussian16 code (Frisch et al., 2016) has been used to perform all the DFT calculations, in conjunction with the double zeta 6-31+G(d) basis set (Petersson et al., 1988; Petersson & Al‐Laham, 1991). Similarly to previous works on other PAH-C60 adducts, the B3LYP functional (Stephens et al., 1994) has been chosen to describe the exchange and electronic correlation (Beheshtian et al., 2012; Khodam Hazrati & L. Hadipour, 2016; El Bakouri et al., 2018). Large aromatic molecules like PAHs and C60 can exhibit strong long-range forces (London forces) when they interact with each other (Ehrlich et al., 2013). In order to account for such forces, the third-level empirical dispersion correction of Grimme (GD3) has been also included in the calculations (Grimme et al., 2010). Each PAH-C60 adduct, after full optimization, was verified to be a stationary point with non-imaginary vibrations to be considered as a minimum. A thermodynamic analysis for all the PAH-C60 adducts was conducted at the temperature of 300 K, along with the corresponding zero-point energy correction.

The IR vibrational spectra of all species were obtained within the framework of the harmonic oscillator and, subsequently, the harmonic frequencies were adjusted by applying a triple-scaling factor scheme to account for anharmonicity, vibro-rotational couplings, etc. (see next Section). The IR intensity has been modeled by a Lorentzian function as peak profile and with a full-width at half-maximum (FWHM) of 0.02 μ\mum. This approximately reproduces the average spectral resolution (R\sim1700) of the Mid-IR Instrument (MIRI) onboard the JWST, working in the \sim5-30 μ\mum spectral range. One of our goals is to provide the theoretically simulated IR spectra of PAH-C60 adducts to the astronomical community as an useful tool to interpret JWST observational data.

Scaling factors

The predicted harmonic frequencies from quantum-chemical calculations must be corrected with scaling factors to avoid the overestimation with respect to the experimental fundamentals (Zapata Trujillo & McKemmish, 2022a). Generally, scaling factors are determined by a fitting procedure of the theoretical IR spectral data against the experimental ones (see e.g. Xu et al., 2023). In our case, we lack experimental data for all PAH-C60 adducts under study and the standard procedure for the scaling factors fitting is not possible. Thus, we have applied an extrapolation method on existing scaling factors, depending of the PAH type bonded to the C60 cage and based on the equation:

νiSFaωi\nu_{i}\simeq\mathrm{SF}_{a}\cdot\omega_{i} (1)

where νi\nu_{i} and ωi\omega_{i} stand for experimental and theoretical frequencies, respectively, with the ithi^{th} suffix indicating that they correspond to the same normal mode and SFa\mathrm{SF}_{a} is the scaling factor. The Eq. 1 above can be reformulated as:

νiSFaωi=SFaωi\nu_{i}\simeq\mathrm{SF}_{a}\cdot\omega_{i}=\mathrm{SF}_{a}^{{}^{\prime}}\cdot\omega_{i}^{{}^{\prime}} (2)

where Eq. 2 describes the relation between harmonic frequencies obtained from different quantum-chemical calculations and is fulfilled if the two sets of harmonic frequencies (ωi\omega_{i}, ωi\omega_{i}^{{}^{\prime}}) are obtained at the same level of theory, producing thus scaling factors (SFa\mathrm{SF}_{a},SFa\mathrm{SF}_{a}^{{}^{\prime}}) that achieve a good agreement with the experimental values (νi\nu_{i}). Formally, the accuracy of harmonic frequencies not only depends on the level of theory used, but also on the basis set applied in the description of the quantum-chemistry model (Zapata Trujillo & McKemmish, 2022a, b). However, recent theoretical studies demonstrate that the median error in the obtained harmonic frequencies with respect to their experimental fundamental vibrations depends more on the level of theory selected than on the basis set applied (Zapata Trujillo & McKemmish, 2023). On the other hand, the scaling factors can correct the harmonic frequencies by applying a global (single number) scaling factor or a frequency-range-specific (multiple numbers) scaling factor. Frequency-range-specific scaling factors exhibit a higher accuracy than their global counterparts, which is due to the frequency range division performed. The frequency range is usually divided into three different frequency regions, each one of them with its corresponding scaling factor (Zapata Trujillo & McKemmish, 2023).

Refer to caption
Figure 1: Structure of the several PAHs considered to build the PAH-C60 adducts. The chemical formula, name and notations are also displayed.

Following the concepts mentioned above to obtain new scaling factors for the PAH-C60 adducts, we have applied Eq. 2 in combination with the frequency-range-specific scaling factors (at << 5 μm\mu\rm m, 5-10 μm\mu\rm m and >> 10 μm\mu\rm m) from the NASA Ames PAH IR Spectroscopic Database (e.g. Bauschlicher et al., 2018; Boersma et al., 2014; Mattioda et al., 2020, AMES PAH Database hereafter). The frequency-range-specific scaling factors of the AMES PAH Database have been computed at the 4-31G/B3LYP level for a large variety of PAH molecules, including the set of small PAHs of our interest (see Figure 1)555For the particular case of indene, we only obtained the adduct formed by 2H-indene and C60. According to our calculations, the 1H-indene, which is the most stable isomer, transforms into 2H-indene upon bonding to the C60 cage.. In addition to our standard 6-31+G(d)/B3LYP+GD3 approach, we have also computed the resulting harmonic frequencies at the 6-31+G(d)/B3LYP level only. The main purpose is to account for the effect of the dispersion correction (GD3) over the scaling factors.

Table 1: Frequency-range-specific scaling factors for the PAHs in Figure 1.666The scaling factors from the AMES PAH Database (SFAMES3{}_{3}^{\mathrm{AMES}}) are compared to those from this work: 6-31+G(d)/B3LYP (SFB3LYP3{}_{3}^{\mathrm{B3LYP}}) and 6-31+G(d)/B3LYP+GD3 (SFD33{}_{3}^{\mathrm{D3}}). The maximum relative errors (Max.r.e) of the scaled frequencies are also included for both SFB3LYP3{}_{3}^{\mathrm{B3LYP}} and SFD33{}_{3}^{\mathrm{D3}}. The specific frequency (or wavelength) range of harmonic frequencies for each scaling factor in this work (SF3B3LYP,SF3D3\mathrm{SF_{3}^{B3LYP}},\mathrm{SF_{3}^{D3}}) is according to the standard range set proposed by Zapata Trujillo & McKemmish (2022a, b). The wavelength ranges are highlighted for indene (In) as example: , , . Max.r.e has been computed using wavelength units.
PAH SF3AMES\mathrm{SF_{3}^{AMES}} SF3B3LYP\mathrm{SF_{3}^{B3LYP}} Max.r.e.B3LYP\mathrm{\overset{B3LYP}{Max.r.e.}} SF3D3\mathrm{SF_{3}^{D3}} Max.r.e.D3\mathrm{\overset{D3}{Max.r.e.}}
0.9563aaaaaa>> 10 μm\mu\rm m 0.9834 0.0132 0.9852 0.0269
In 0.9523bbbbbb5-10 μm\mu\rm m 0.9677 0.0124 0.9675 0.0151
0.9595cccccc<< 5 μm\mu\rm m 0.9612 0.0025 0.9614 0.0023
\arrayrulecolorgray 0.9563 0.9830 0.0099 0.9838 0.0093
Iyl 0.9523 0.9665 0.0074 0.9664 0.0114
0.9595 0.9621 0.0026 0.9627 0.0028
0.9794 0.9815 0.0109 0.9807 0.0190
An 0.9691 0.9734 0.0103 0.9711 0.0114
0.9597 0.9603 0.0001 0.9610 0.0003
0.9563 0.9830 0.0103 0.9818 0.0118
Tn 0.9523 0.9677 0.0120 0.9679 0.0037
0.9595 0.9608 0.0012 0.9614 0.0009
0.9563 0.9854 0.0402 0.9848 0.0500
Pn 0.9523 0.9676 0.0102 0.9648 0.0113
0.9595 0.9612 0.0013 0.9618 0.0010
\arrayrulecolorblack

Table Scaling factors summarizes the comparison of frequency-range-specific scaling factors obtained through Eq. 2 with respect to those of the AMES PAH Database (Bauschlicher et al., 2018; Boersma et al., 2014; Mattioda et al., 2020). The new scaling factors SF3B3LYP\mathrm{SF_{3}^{B3LYP}} and SF3D3\mathrm{SF_{3}^{D3}} provide the same level of accuracy, which is clearly reflected in the maximum relative error (Max.r.e.)777The Max.r.e is defined as the maximum value from the relative error of the new scaling factors with respect to the AMES PAH Database ones. The Max.r.e. is calculated using the scaled harmonic frequencies for each method according to the equation: Max.r.e<SF3newSF3AMESSF3AMES\rm Max.r.e<\frac{SF_{3}^{new}-SF_{3}^{AMES}}{SF_{3}^{AMES}}., being very similar for both methods. The Max.r.e indicates the deviation of the scaled frequencies obtained from SF3B3LYP\mathrm{SF_{3}^{B3LYP}} and SF3D3\mathrm{SF_{3}^{D3}} with respect to those of SF3AMES\mathrm{SF_{3}^{AMES}}. Our new scaling factors mostly have a maximum relative error of the order of \sim10-2, with a slightly larger value (\sim5×\times10-2) at wavelengths longer than 10 μ\mum (see Table Scaling factors). However, the spectral resolution of the experimental IR spectra used to obtain the AMES PAH Database scaling factors SF3AMES\mathrm{SF_{3}^{AMES}} is significantly higher than our Max.r.e. values (see e.g. Mattioda et al., 2020). So, the error introduced by the new scaling factors SF3B3LYP\mathrm{SF_{3}^{B3LYP}} and SF3D3\mathrm{SF_{3}^{D3}} can be considered as negligible. In short, we conclude that both scaling factors SF3B3LYP\mathrm{SF_{3}^{B3LYP}} and SF3D3\mathrm{SF_{3}^{D3}} can reproduce fairly well the vibrational frequencies (or IR features) of the small PAHs displayed in Figure 1. In the following, we only use the scaled frequencies as obtained from the new scaling factors SF3D3\mathrm{SF_{3}^{D3}}.

Scaling factors validation

Despite the scaling factors SF3D3\mathrm{SF_{3}^{D3}} reliably reproduce the experimental IR features of small PAHs (Fig.1), here we are interested in their adducts with the C60 molecule. Therefore, it is not appropriate to directly extrapolate our SF3D3\mathrm{SF_{3}^{D3}} scaling factors to PAH-C60 adducts before a proper validation. In order to validate our scaling factors, we have used the laboratory IR spectra measurements of three different types of PAH-C60 adducts (see Figure 2). These three archetypal systems have been synthesized via a novel supramolecular mask strategy (Pujals et al., 2022), which allows us to selectively obtain the bis-adducts (two PAH units) or mono-adducts (one PAH unit). It is worth noting that the synthesis route via supramolecular mask strategy can pinpoint to the desired PAH-C60 adduct isomer. Thus, we are completely confident that the laboratory IR spectra are not contaminated by a mixture of regioisomers.

The experimental IR spectra of the three PAH-C60 adducts in Figure 2 were obtained by a Bruker ALPHA II FT-IR spectrometer with a spectral resolution of ±\pm2 cm-1 equipped with a DTGS (deuterated triglycine sulfate) detector. All measurements were carried out in the solid state of the sample by using a Bruker Platinum ATR Adapter (attenuating total reflectance module with a diamond crystal support), allowing a direct measurement without prior treatment of the sample. All spectra were recorded at room temperature and under air atmosphere. The FT-IR measurements were recorded between 2.5 and 25 μ\mum (4000-400 cm-1). Each final spectrum was obtained by averaging 24 recorded scans.

Refer to caption
Figure 2: The three archetypal PAH-C60 adducts synthesized by supramolecular mask strategy and used for the scaling factors validation: 1AnC60C_{60} (mono-anthracene), 2AnC60C_{60} (bis-anthracene) and 1PnC60C_{60} (mono-pentacene) (see Pujals et al., 2022, for more details). The notation follows the legend in Figure 1, adding the number of PAH units and the C60 identification.

In order to quantitatively assess the agreement between the experimental and theoretical IR spectra, we have followed the procedure of the cosine similarity score recently introduced in the literature (Fu & Hopkins, 2018; Kempkes et al., 2019; Müller et al., 2020). This similarity score method calculates the cosine of the angle θ\theta between two n-dimensional vectors using their normalized Euclidean dot product according to:

Similarity=cos(θ)=ABAB=i=1nAiBii=1nAi2i=1nBi2\rm{Similarity=cos(\theta)=\frac{A\cdot B}{||A||\,||B||}=\it{\frac{\sum_{i=1}^{n}A_{i}B_{i}}{\sqrt{\sum_{i=1}^{n}A_{i}^{2}}\cdot\sqrt{\sum_{i=1}^{n}B_{i}^{2}}}}} (3)

Where the elements AiA_{i} and BiB_{i} are the experimental (A) and theoretical (B) IR intensity at the same ithi^{th} frequency. Concretely, the procedure is to compute the similarity score at the same reference frequency; both experimental and theoretical spectra have thus the same xx-axis. In this case, the experimental spectrum is taken as reference and the similarity score was calculated at the exact experimental wavenumber. A similarity score close to 1 indicates a higher resemblance between the experimental and theoretical spectra.

On the other hand, Kempkes et al. (2019) introduce an alternative definition of the IR intensity to make the cosine similarity score more sensitive to the overlap of frequencies in the spectra AA and BB of Eq. 3, but diminishing the importance in the deviations of their IR intensities (see Kempkes et al., 2019, for more details). The Kempkes et al. (2019) IR intensity definition stand as:

Airev=log(AiAm+c)A_{i}^{rev}=log\left(\frac{A_{i}}{A_{m}}+\it{c}\right) (4)
Refer to caption
Figure 3: Experimental (red) versus scaled harmonic B3LYP-D3/6-31+G* (black) spectra for the three archetypal PAH-C60 adducts used for the scaling factors validation: 1AnC60C_{60}, 2AnC60C_{60} and 1PnC60C_{60} (Fig.2). The experimental spectra of the pristine PAHs (light blue) are included to ease the discussion. (A): Spectra built from scaling factors in Table Scaling factors. (B): Spectra built from empirical-corrected scaling factors for wavelengths ¡ 5 μ\mum ( ¿ 2000 cm-1) (see the text for more details). The similarity score value for each case is also indicated. The theoretical spectra have been constructed in the wavenumber scale by the convolution of a Lorentzian function with FWHM = 4 cm-1, in order to ease the comparison with the experimental spectra.

The variable AirevA_{i}^{rev} denotes a new recalculated intensity, which is computed from the normalized intensity up to a maximum of 1 \equiv (AiAm\frac{A_{i}}{A_{m}}) and weighted by a factor cc. We apply this procedure (Eq. 4) to both the experimental (A) and theoretical (B) IR spectra, compromising the similarity score with the sensitivity to low-intensity bands in the spectrum, but avoiding the experimental noise. The factor cc is a constant for both A and B, which has to be derived from the correlation between the experimental and theoretical spectra. However, due to the lack of a significant number of resolved IR features (i.e. poor statistics) in our PAH-C60 adducts experimental spectra (Fig.3), we have not correlated them with the theoretical ones. Instead, we have used a factor of cc = 0.71 recently reported by Xu et al. (2023), as obtained from a significant number of experimentally resolved IR features on smaller C60 adducts like CO+60\rm{}_{60}O^{+} and CH+60\rm{}_{60}H^{+}. This way, the recalculated experimental and theoretical intensities, AirevA_{i}^{rev} and BirevB_{i}^{rev}, are used in Eq. 3 to compute the similarity score values for the three archetypal PAH-C60 adducts and indicated in the legend of Fig.3.

In Figure 3a we show that the scaled harmonic theoretical IR spectra of our three archetypal PAH-C60 adducts are in good agreement with the experimental spectra regarding the vibrational frequencies. This is evidenced by the similarity scores obtained, which resemblance to those reported in previous experimental-theoretical studies (Müller et al., 2020; Kempkes et al., 2019). However, according to our experimental spectra, the archetypal PAH-C60 adducts do not display IR features below \sim 3.3 μ\mum (above 3000 cm-1). Apparently, our scaling factors produced from free PAHs cannot reproduce this spectral region for PAH-C60 adducts. Thus, we have opted to apply an empirical correction to the scaling factors in Table Scaling factors, but only in the ¡ 5 μ\mum spectral range. The procedure imply a shift of the theoretical spectra until a maximum in the similarity score was obtained with respect to the experimental data. Figure 3b illustrates the empirical-corrected theoretical spectra this way, showing \sim 3.3 μ\mum region in better agreement with the experimental data. Using these three archetypal PAH-C60 adducts we have found an empirical wavelength correction of 0.1984 μ\mum, which was applied to the scaling factors for the rest of PAH-C60 adducts. The corresponding scaling factors finally adopted can be found in Appendix A. It is worth noting that the experimental spectra of the pristine PAHs and PAH-C60 adducts are strongly different in the \sim3.3-3.6 μ\mum region (see Figure 3b); something that is discussed in the following sections.

However, deviations in the intensities are noticeable; specially in the case of 2AnC60, which shows the largest intensity deviations of the three PAH-C60 adducts (e.g. at wavelengths longer than 7 μ\mum but also in the \sim3.3-3.5 μ\mum spectral region). The difference in intensity between the theoretical and experimental IR spectra is likely due to several reasons. The collection of the experimental spectra is performed in solid state, which can modify the IR intensity with respect to gas-phase measurements like IR multiple-photon dissociation (IRPMD, Polfer, 2011; Palotás et al., 2020) or even more sophisticated techniques (Gerlich et al., 2018; Roithová et al., 2016). In addition, we can not completely discard whether a better fitting of the experimental spectra baseline and/or a reduction of the noise contribution could improve the intensity match with the theoretical spectra. Nevertheless, we note that the experimental measurements of the three archetypal PAH-C60 adducts are only used to validate our new derived scaling factors (see Table Scaling factors and Appendix A), which, according to the similarity score values, perform quite well regarding the vibrational IR frequencies. The intensity mismatch is a long-standing and well-known problem in the comparison of experimental and theoretical IR spectra of complex organic compounds (see e.g. Katsyuba et al., 2013), and its resolution is, of course, beyond the scope of the present work.

Multiple PAH-C60 adducts models

In the following sections we present the theoretical models and simulated IR spectra of several PAH-C60 adducts containing multiple PAH units; in many cases they were described by more than one spatial configuration (isomers). The mono-adducts notation refers to those C60 adducts formed by the addition of one PAH unit only, while bis-adducts and tris-adducts stand for two and three PAH units, respectively. For the largest PAHs considered here, tetracene (Tn) and pentacene (Pn) (see Figure 1), only a maximum of two units and one unit, respectively, were attached to the C60 cage; mainly due to our computational limitations.

PAH-C60 adducts comparison

Figure 4 displays the theoretical IR spectra of the different mono-adducts models together with their molecular structure. At first glance, the IR spectra of the mono-adducts are richer than those corresponding to the isolated PAH and C60 molecules (Figure 4). Especially, in those spectral regions (e.g. the \sim12-13 and 14-17.5 μ\mum regions, among others, see Figure 4 and Appendix F) where no contributions from the corresponding PAH and C60 are observed. In addition, the typical \sim3.3 μ\mum feature of the pristine PAH appears red-shifted in the case of PAH-C60 adducts. The PAH binding to the C60 cage is the cause for all new emission features observed in the IR spectra, denoting several unique spectral regions to distinguish PAH-C60 mono-adducts. The black arrows in Figure 4 also indicate those specific features that are free of contribution from C60, the pristine PAH and even from the other adducts. In some adducts, like those with In and Iyl, it is very difficult to identify unique (specific of the adduct) features since they contribute throughout all the IR spectra. The peculiarities of In or Iyl adducts are due to the multiple possible ways of binding to the C60 cage. Single and double-bond bindings have also been found in previous studies about C60-coronene+ adducts, but in this case inducing the loss of one or two H atoms, which can be another types of binding between PAHs and C60 reacting under energetic conditions (see Figure 2 in Dunk et al., 2013). In contrast, the rest of mono-adducts could only be obtained under a double-bonded geometry (see Figure 4). Nevertheless, all mono-adducts have common emitting spectral regions with features at \sim3.3-3.6, 6-10, 12-16, and 17-19 μ\mum.

Refer to caption
Figure 4: Scaled theoretical IR spectra of the mono PAH-C60 adducts (one PAH unit: C60\rm C_{60}-CxHy\rm C_{x}H_{y}) modeled at the B3LYP-D3/6-31+G* level. All spectra were convolved with a Lorentzian function of FWHM = 0.02 μ\mum. In the case of 2H-indene (In) and indenyl (Iyl), the blue and red arrows highlight the binding of the PAH with a pentagon and hexagon ring to the C60 cage, respectively. For indenyl, with a single bond connecting its pentagon to C60, a green arrow has been used instead. In all panels, the red dashed spectra correspond to the isolated pristine PAH, while the blue dashed lines mark the four strongest C60 features (\sim7.0, 8.5, 17.4 and 18.9 μ\mum). The black arrows on the spectra indicate those characteristic features (specific of each adduct), which are free of contribution from C60, the pristine PAH and even from the other adducts.

The 3.3-3.6 μ\mum feature can become broader (or splitted into various peaks) and less intense as the size of the PAH attached to the C60 cage decreases (see Appendix F). The former trend is only an indication of the reduction in the number of CH bonds and a breakdown of the symmetry as the PAH is smaller. Less CH bonds reduce the intensity of the 3.3-3.6 μ\mum feature, while asymmetric PAHs like In and Iyl make the feature broader (or even with additional resolved features at slightly longer wavelengths of \sim3.5-3.6 μ\mum) because of the different chemical environments affecting the CH bonds.888In the case of pure hexagon PAHs like anthracene (An), tetracene (Tn) and pentacene (Pn), the maximum different CH bonds in the mono-adducts is two, which are those CH close and far from the binding. In contrast, 2H-indene (In) and indenyl (Iyl) can have more than four different CH bonds due to the presence of pentagons and hexagons, broadening the 3.3-3.6 μ\mum feature. This is quite noticeable for In, which contains CH bonds connected to a Csp3\rm{}_{sp_{3}}, creating multiple chemically different CH stretching. In fact, In and Iyl show a clear increase in the IR intensity of the CH stretching region when the PAHs are bonded through the pentagon ring; as a consequence of a more marked dipole moment change for the CH vibrations inside the pentagon than for the hexagon.

At 6-10 μ\mum the trend is similar in all spectra, a number of new low-intensity IR features; with 1InC60{}^{{}^{\prime}} being the spectrally richer model (see Figure 4). The C-C vibrations lay in this region, which, as expected, shows the major coincidence (in terms of features) with its progenitors the free PAHs and the C60 cage. In particular, all mono-adducts strongly contribute to the 7.0 μ\mum C60 feature, with almost a negligible emission contribution at 8.5 μ\mum. The enrichment of the IR spectra in the 6-10 μ\mum spectral region is smaller than for other fullerene based species like metallofullerenes (Barzaga et al., 2023a, b). This indicates that the dipole moment change related to C-C vibrations is weakly modified by the binding between the PAHs and C60. Consequently, charge reordering and charge transfer processes in the mono-adducts are less important than for metallofullerenes (Barzaga et al., 2023b).

The 12-16 μ\mum spectral region is very interesting because it is free from C60 emission and the isolated PAHs usually only display one strong feature around \sim13-14 μ\mum, but their mono-adducts counterparts display richer spectra. This region is mainly dominated by the out-of-plane CH vibrations, which is less IR active when the PAH is alone. In the particular case of 2H-indene (In) and indenyl (Iyl), the PAH can bind to the C60 cage through its hexagon or pentagon ring999In the case of tetracene (Tn), we have also considered another possible way of binding, but it is highly thermodynamically unstable and also lacks of chemical principles (see Appendix E for more details). (see Figure 4). Both In and Iyl bind more favorably to the C60 cage through the pentagon ring with energies of -1.62 eV (1InC60{}^{{}^{\prime}}) and -1.41 eV (1IylC60), respectively. This specific binding is also reflected in the change of the IR spectra in the 12-16 μ\mum region, with the appearance of new IR active features. On the other hand, in the case of Iyl single-bonded to the C60 (1IylC60, marked with a green arrow in Figure 4), multiple out-of-plane CH vibrations become IR active. This is due to the loss of symmetry compared to the doubled-bonded 1IlylC60{}^{{}^{\prime}} model; the vibrational modes appear splitted, producing a richer spectra in the 12-16 μ\mum range.

Finally, another interesting spectral region is the one from 17 to 19 μ\mum. This region looks quite similar in the mono-adducts, with several new weak features emerging around 17.4 μ\mum and a broader 18.7-19.2 μ\mum peak whose strongest intensity is at 18.9-19.0 μ\mum. Clearly, this implies an important contribution to the 18.9 μ\mum C60 feature, with only a marginal emission contribution at 17.4 μ\mum where there is a marked decay in intensity (see Appendix F). The 18.9 μ\mum feature of pristine C60 corresponds to a combined carbon cage vibration, but the symmetry is destroyed by PAH binding.101010This occurs for all the PAH-C60 adducts under study, with the only exception of Iyl, which shows a broadness reaching up to 19.2 μ\mum due to the combination with the intrinsic PAH vibrations. Curiously, the largest mono-adducts from Tn and Pn are the only ones exhibiting noticeable IR features at wavelengths >>20 μ\mum; in particular, 1PnC60 with a relatively strong feature at 20.7 μ\mum. These features at long wavelengths are intrinsic to the free PAHs, although they are perturbed by binding to C60, causing a shift in wavelength up to almost \sim0.8 μ\mum in the case of Tn (see Figure 4).

Refer to caption
Figure 5: Theoretical scaled IR spectra of PAH-C60 bis-adducts from 2H-indene (2InC60, a-c) and indenyl (2IylC60, d-e). The red dashed spectra correspond to the mono-adduct analogues, while the blue dashed lines mark the four strongest C60 features (\sim7.0, 8.5, 17.4 and 18.9 μ\mum). Note that the convolution parameters are the same as in Figure 4. For a detailed discussion on the models and spectral features we refer to Appendix B.

Figure 5 shows a comparison between the theoretical scaled IR spectra of mono-adducts and bis-adducts, indicating that no significant differences are observed. The changes from mono- to bis-adducts are only in terms of IR intensity caused by the increment of PAHs units, but regarding the features position, the spectra remain unchangeable. The only exception is for the model 2IylC60 h and is provoked by a strong structural deformation of the C60 cage (see Appendix B for more details and also Appendix F). Furthermore, the same spectral behavior is observed for the tris-adducts when they are compared with the bis-adducts (see Appendix C). The spectral differences among all these PAH-C60 adducts are more appreciated at low IR intensity, which, in most cases, are expected to be well below the detection threshold of astronomical IR observations (see Appendices B and C).

PAH-C60 adducts stability and kinetics

Thermodynamic stability

The common mechanism for the formation of C60 adducts is by exohedral addition of new molecular species, specially alkyl or benzyl groups, which is due to the high reactivity of C60 (Hirsch, 1995; Krusic et al., 1991; Li, 2022). We can not confirm or refute that exohedral addition is the appropriate chemical process taking place in C60-rich astrophysical environments like young PNe, characterized by both low temperatures (\sim300 K; see e.g. Cami et al., 2010; García-Hernández et al., 2012b) and densities; but it is the most likely one taking into consideration the well-known C60 reactivity. In our case, multiple PAHs addends may undergo this exohedral cyclo-addition through sequential reaction pathways:

  1. I.

    PAH+C60{}_{60}\xrightarrow[]{} PAH-C60

  2. II.

    PAH+PAH-C+PAH60{}_{60}\xrightarrow[]{\rm+PAH} [PAH]2\rm[PAH]_{2}-C60

  3. III.

    PAH+[PAH]2\rm[PAH]_{2}-C+PAH60{}_{60}\xrightarrow[]{\rm+PAH} [PAH]3\rm[PAH]_{3}-C60

The three pathways (I-III) imply the sequential (step-by-step) addition of PAHs to C60 from mono- to tris-adducts. These reactions depend on different factors that determine the mechanism of formation, but also lack the description of UV radiative processes, which are commonly observed in many astrophysical sources, including C60-rich environments. Such analysis requires a detailed description of the several reaction pathways with the inclusion of photo-dissociation and photo-reaction, which involves the development of precise reactivity models; something that is extensively time-consuming and beyond of the scope of the present work. Thus, we have opted here for the derivation of the thermodynamic stability in order to identify the most spontaneous reactions following the sequential pathways (I-III). Assuming this approach implies that reactions take place under thermodynamical equilibrium. The former could be in contradiction with the observational data that indicate uncertainties in the estimation of physical conditions of astrophysical environments (see e.g. Brieva et al., 2016), but our models can still be used as a rough prediction. Thus, we have chosen the difference in Gibb’s free energy of formation per PAH unit (ΔGf\rm\Delta G_{f}), defined by:

ΔGf=GproductsGreactantsn\rm\Delta G_{f}=\frac{\sum G_{products}-\sum G_{reactants}}{n} (5)

where Gproducts\rm G_{products} and Greactants\rm G_{reactants} are the Gibb’s free energy of products and reactants, respectively, which have to be summed according to the reactions pathways I-III mentioned above. The variable n indicates the number of PAH molecules attached to the C60 cage and is used to normalize the energy with respect to the PAHs units. A negative and positive value of ΔGf\rm\Delta G_{f} indicates a spontaneous and non-spontaneous process, respectively. It is important to clarify that ΔGf\rm\Delta G_{f} is a thermodynamic quantity that describes the probability that reactants become products; however, there is no information on the activation energy to carry out these reactions. In this sense, even though ΔGf\rm\Delta G_{f} provides insight into the binding energy of the PAH on C60, it can not be used as a threshold for any reaction mechanism (e.g. dissociation, photodissociation, etc.). Theoretical studies on the C60-anthracene+ adduct suggest that the energy barriers in PAH binding to C60 can be as small as 0.08 eV, which implies that the activation energy is almost equal to the binding energy (see Figure 6 in Zhen et al., 2019b). In principle, this could be used to infer the range in activation energies of our PAH-C60 adducts because they also follow a Diels-Alder binding. However, it is likely that differences in the thermodynamics could appear due to the fact that our models are neutral systems. Table Thermodynamic stability summarizes the thermodynamic data for all PAH-C60 adducts under study; all values have been computed assuming a temperature of 300 K and adequately corrected by the zero-point energy.

Table 2: Thermodynamic values of all the PAH-C60 adducts111111Gibb’s free energy of formation per PAH unit (ΔGf\rm\Delta G_{f}) for each of the PAH-C60 adducts studied here (Figures 4, 7 and 8). The binding energies (Eb\rm E_{b}) for the mono-adducts are also indicated. Each ΔGf\rm\Delta G_{f} value corresponds to the sequential reaction pathways (I-III). All energies are described in eV units.
Mono-adduct Bis-adduct Tris-adduct
PAH Model Eb\rm E_{b} ΔGf(I)\rm\Delta G_{f}(I) Model ΔGf(II)\rm\Delta G_{f}(II) Model ΔGf(III)\rm\Delta G_{f}(III)
In 1InC60{}^{{}^{\prime}} -0.11 -1.015 2InC60 a -0.976 3InC60 a -0.953
2InC60 b –1.010
2InC60 c -0.978 3InC60 b -0.974
\arrayrulecolorgray Iyl 1IylC60{}^{{}^{\prime}} 0.35 1.000 2IylC60 g 1.820 3IylC60 c 1.052
2IylC60 f 1.022 3IylC60 d 1.075
1IylC60 -0.46 0.102 2IylC60 e -0.056
2IylC60 d 0.075 3IylC60 e -0.064
An 1AnC60 -0.49 0.121 2AnC60 h 0.124 3AnC60 f 0.183
3AnC60 g 0.154
Tn 1TnC60 -0.33 -3.261 2TnC60 i -2.389
Pn 1PnC60 -1.05 -0.426
\arrayrulecolorblack

According to the values in Table Thermodynamic stability, tetracene (Tn) exhibits the most spontaneous (exergonic) reactions for any of the reaction pathways assumed. Increasing the number of Tn addends makes the reactions less exergonic. Such trend in ΔGf\rm\Delta G_{f} is seen for almost all PAHs, with the exception of some C60 adducts with In and Iyl. This generally indicates that the inclusion of more PAH units to C60 is energetically hindered by the geometrical arrangement of multiple PAHs to the carbon cage, which is known as steric effect. In contrast, the thermodynamics of the indenyl (Iyl) adducts is the most complex one since there is no clear relationship between the ΔGf\rm\Delta G_{f} values and the reaction pathways and/or regioisomer (Table Thermodynamic stability). However, it is clear that the single-bonded Iyl-C60 adducts models are more exergonic than the double-bonded ones (see all adducts models from 1IylC60 in Table Thermodynamic stability). Indenyl is a planar and highly unsaturated radical that tends to change from Csp2{}_{\rm sp_{2}} to Csp3{}_{\rm sp_{3}} when it binds to C60. A similar behavior is known to occur for other aromatic addends (e.g. Mas-Torrent et al., 2002).

Our predictions also show that the formation of full hexagon PAH-C60 adducts with anthracene (An) is an exothermic process (see 𝐄𝐛\rm\mathbf{E_{b}} in Table Thermodynamic stability); similar theoretical results have been reported for C60-anthracene+ derivatives in recent studies (Zhen et al., 2019a; Wu et al., 2024). Steric and electrostatic effects should be the main factors for the thermodynamic changes when the amount of PAH addends increases, but our models contain both symmetric and asymmetric PAHs, and the thermodynamics information is thus not straightforward to interpret. Seemingly, there is no correlation between the increment in the size of fused-hexagon PAHs (An << Tn << Pn) and the observed thermodynamics (see Table Thermodynamic stability). Future intensive efforts would be necessary to include several factors in the mechanistic description of these chemical reaction pathways.

Kinetics of PAH-C60 adducts

Although in the present study no precise reactivity models have been developed to predict the kinetics of our PAH-C60 adducts, previous experimental works have determined these parameters for some of the adducts presented here. This Section covers the state-of-the-art experimental kinetics, concerning the rate constants of the mono- and bis-adducts of C60 with anthracene, tetracene and indene.

Kinetics analyses of PAH-C60 adducts formation have been mainly made during Diels-Alder synthesis, which is the common experimental method to obtain them. Sarova & Berberan-Santos (2004) study the chemical kinetics of C60 mono-adducts from anthracene and tetracene combining experimental measurements and theoretical predictions. By using a fitting procedure, they determine rate constants following a second-order kinetics with values of 1.6×\times10-4M-1s-1 and 3.1×\times10-2M-1s-1 for anthracene and tetracene, respectively, at T\sim298 K (Sarova & Berberan-Santos, 2004). Such values indicate that the tetracene mono-adduct kinetics is much faster than for anthracene at least by a factor of 200. As we will see in the next Section, these kinetics data can be used, at least when the small PAHs like anthracene and tetracene are equally abundant, for the construction of more reliable theoretical IR spectra representative of PAH-C60 adducts.

Refer to caption
Figure 6: The DFT simulated IR (\sim3-25 μ\mum) spectra of the mixture of PAH-C60 adducts. (A): the total (summed) mixture spectrum; (B): the total abundance-weighted spectrum using the indene abundance from Burkhardt et al. (2021); and (C) the total abundance-weighted spectrum using the indene abundance from Cernicharo et al. (2021). In all panels, the red dashed lines mark the most common astronomical PAH features (at \sim3.3, 6.2, 7.7, 8.6, 11.2 and 12.7 μ\mum), instead blue dashed lines mark the four strongest C60 features (\sim7.0, 8.5, 17.4 and 18.9 μ\mum). In the case of total the (summed) mixture spectrum, the intensity has been normalized to the maximum intensity peak, while the abundance-weighted spectra have been normalized with respect to the higher intensity peak from the weighted spectrum (C); i.e. the one using the indene abundance from Cernicharo et al. (2021).

More recently, the kinetics of solid indene-C60 mono- and bis-adducts decomposition have been reported by Rodrigues et al. (2023). The most remarkable result of the previous study is the relative rate constants between the mono- and bis-adducts, being one of the few scarce data of this type reported in the literature. Accordingly, a first-order kinetics defines the decomposition of mono-adducts carrying out a slower process with respect to bis-adducts (at least by a factor of 2). Arguably, solid-state kinetics can differ from the gas-phase since different forces may rule out the reaction. We are thus aware that the kinetics data by Rodrigues et al. (2023) could be far from the chemical processes occurring in astrophysical environments. However, they can be still used to build more reliable and chemical intuitive preliminary models, describing the IR emission of these PAH-C60 adducts and elaborated below.

Astrophysical relevance

A comparison between our theoretical spectra of PAH-C60 adducts and the observational ones is not straightforward because it is well known that the relative abundances of the several species may play a key role in the landscape of the IR spectra (see e.g. Barzaga et al., 2023a, b). Also, all previous IR spectroscopic observations of fullerene-rich astrophysical environments were carried out by the Spitzer Space Telescope, which did not cover the \sim3-5 μ\mum region. As we will see through the current Section, this is a key spectral region where the PAH-C60 adducts strongly contribute to the IR emission.

The main possible drawback of the kinetics information in the previous Section could be that the rate constants are obtained from experimental reactions where anthracene and tetracene always exceeded C60. The abundances of such small PAHs and C60 in astrophysical environments are not precisely known; mainly because their more specific UV/Visible electronic transitions have not been detected yet and only upper limit estimates to their abundances121212Note that the C60 abundance estimates from the IR emission usually assume thermal emission as excitation mechanism (see e.g. Cami et al., 2010; García-Hernández et al., 2012b, and references therein) and the IR bands attributed to C60 may actually comprise contributions by analogous species like metallofullerenes (see e.g. Barzaga et al., 2023a, b). In fact, Brieva et al. (2016) concluded that thermal excitation by itself cannot explain the IR spectra attributed to C60 in some planetary nebulae, and also that a possible explanation is that other molecules than C60 may contribute to the observed spectra. are available (see Rouillé et al., 2021, for a recent review about C60). Another difficulty is that present abundance estimates are not coming from the same astrophysical source (e.g. the diffuse ISM, PNe, RCBs, etc.) and/or astronomical spectra. By searching in the literature, we found diffuse ISM abundance-limit estimates of both anthracene and C60 towards the same sources (the reddened stars HD 169454 and HD 183143) and using the same astronomical spectra (Gredel et al., 2011; Rouillé et al., 2021). The upper limits for column densities of interstellar C60 and anthracene are in the ranges \sim1-8 ×\times 1012 cm-2 and 1-3 ×\times 1012 cm-2, respectively. Other small PAHs like pyrene and 2,3-benzofluorene studied by Gredel et al. (2011) display upper limits for column densities on the same order (\sim2-8 ×\times 1012 cm-2). The upper limits for C60 column densities in the ISM are very similar to those estimated towards the C60-rich circumstellar environments around PNe and RCB stars (\sim1-4 ×\times 1012 cm-2, García-Hernández et al., 2012a; García-Hernández & Díaz-Luis, 2013; Díaz-Luis et al., 2015). In short, with the UV/Visible astronomical data at hand, it is not possible to know the relative abundances of C60 and small PAHs in astrophysical environments like the diffuse ISM and C60-rich circumstellar envelopes. However, recent radioastronomy observations have detected indene (In) towards the cold dark cloud TMC-1 (e.g. Cernicharo et al., 2021), with a surprisingly high column density of \sim1.6 ×\times 1013 cm-2 (a slightly lower value of \sim0.96 ×\times 1013 cm-2 was obtained by Burkhardt et al., 2021). Although it is not known if TMC-1 harbor both C60 and indene, such observations suggest that some small PAHs can be found in very high abundance (even higher than C60) in space, as long as they are shielded from strong (inter)stellar UV photons (e.g. Stockett et al., 2025). We note that such a possible UV-shielding, permitting small PAHs to survive, may not be the case for all interstellar and circumstellar environments where fullerenes (C60, C60+ and C70) have been detected. However, these shielding conditions (e.g. UV-shielding by circumstellar dust grains) are usually met in evolved stars (sources between AGB stars and PNe) and neutral molecules, both simple (H2, CO) and complex (PAHs, C60) ones, may be formed in the outer parts or within clumps (see e.g. Manchado et al., 2015; Wesson et al., 2024; Gold et al., 2024; Clark et al., 2025); indeed neutral C60 seems to be distributed in rings and/or clumps in young PNe (Díaz-Luis et al., 2018; Cami et al., 2018). In addition, for the particular case of the circumstellar environments of evolved stars, there is a rapidly changing UV radiation field (from no or very little UV photons in the AGB phase to ionizing UV photons in the PNe phase). At present, it is not known at which exact AGB-PNe evolutionary stage PAHs and fullerenes are formed. It is only known that fullerenes are detected for a very short time when the UV photons from the central star are strong enough to photoionize H (central stars with Teff\sim30-40 kK) but the excitation mechanism (e.g. fluorescence vs thermal) of the C60 IR emission is unknown (e.g. Brieva et al., 2016). In other words, fullerenes could be formed at shielding conditions (and thus coexist with small PAHs for a certain time) but they only could be detected when excited by relatively strong UV photons and emit in the mid-IR; although this is not always the case and C60 emission is also detected in astrophysical environments with rather little UV radiation like proto-PNe and RCB stars.

Therefore, we tentatively estimate the possible order of magnitude of the expected abundances of PAH-C60 adducts by using the thermodynamics131313In this case, we have used the enthalpy energy, which is directly obtained from our quantum-chemistry calculations as well as the free Gibb’s energies listed in Table Thermodynamic stability. and kinetics data presented in the previous Sections (see the Appendix D for more details). We note that our abundance estimation assumes the average thermal conditions of the C60-rich circumstellar envelopes around evolved stars (T\sim300 K) but uses the estimated column densities of indene towards TMC-1 (as obtained from extremely sensitive radio observations, Burkhardt et al., 2021; Cernicharo et al., 2021). This means that our abundance estimates may be useful as rough illustrations only, because the physical conditions in the dense molecular cloud TMC-1 can be very different to those in C60-rich circumstellar environments. This is due to the lack of indene abundance estimates toward fullerene-rich circumstellar environments; e.g. similar extremely sensitive radio astronomical observations towards C60-rich environments are not available in the literature. The DFT simulated mixture spectra in Figure 6 have been constructed following a procedure similar to that in Barzaga et al. (2023a, b); see Appendix D141414In the case of the total abundance-weighted spectra, they both have been normalized with respect to the maximum intensity peak obtained from the Cernicharo et al. (2021) indene abundance data (spectrum C in Figure 6). This way, it is noticeable the change in IR intensity with respect to the lower indene abundance from Burkhardt et al. (2021) (spectrum B in Figure 6). for more details.

From Figure 6, it is very clear the important effect of considering the PAH-C60 abundances in the resulting IR spectra. The pure summed PAH-C60 mixture in Figure 6a exhibits a spectrum, mainly characterized by strong emission features at \sim3.4-3.6, 14.2, 14.7 and 15.1 μ\mum. These features are characteristic of the Iyl bis-adduct under configuration g (Figure 7g), which is highly unstable. This landscape drastically changes when our estimated PAH-C60 abundances are considered in the simulations. In particular, there is a strong modification in the 3.4-3.6 μ\mum feature, which becomes splitted in the abundance-weighted PAH-C60 mixtures, with distinctive peaks at 3.43, 3.51 and 3.57 μ\mum (see Figures 6b, c and Appendix F). Furthermore, the abundance-weighted spectra are spectrally richer along the full \sim5-20 μ\mum spectral range in terms of distinguishable peaks, showing numerous features weaker than the 3.4-3.6 μ\mum one, but still noticeable. The most intense IR features in the latter spectral range are those at \sim6.7-7.4, 12.2-14.0 and 18.6-19.2 μ\mum. In short, the strong differences seen in the DFT spectra in Figure 6 highlight the importance of considering, at least, the tentative abundances of the several species in the theoretical simulation of IR spectra from PAH-C60 adducts.

Interestingly, there is a common factor in the spectra of Figure 6, which is the general absence of IR contribution at the wavelengths of the astronomical PAHs (at \sim6.2, 7.7, 8.6, 11.2, and 12.7 μ\mum) and the C60 features at \sim8.5 and 17.4 μ\mum. Although both species are the progenitors of PAH-C60 adducts, the resulting theoretical spectra with (by far) the most intense IR bands at \sim3.4-3.6 μ\mum don’t indicate their presence. From the astronomical point of view, this suggests that PAH-C60 adducts could contribute to the \sim3.4-3.6 μ\mum IR emission and usually attributed to aliphatic carbon species (see below). The astronomical PAH and C60 features are well characterized and the astrophysical community has generally accepted that their detection indicates the presence of these species. However, according to our tentative predictions, the absence of IR emission from PAHs and/or C60 does not necessarily mean that they are not present; they could be present, instead, forming more complex hybrid species such as PAH-C60 adducts. A recent theoretical work by Xu et al. (2023) finds that much simpler C60 adducts, like C60H+ and C60O+, exhibit most of the pristine C60 features in their predicted IR spectra. Apparently, such a difference is due to the large diversity in the interactions between the PAHs and C60. Thus, the predictions presented here suggest that one of the possible explanations for the lack of PAH and C60 features in astronomical environments could be the formation of their adducts.

Our theoretical calculations may also give some insights about the possible origin of the \sim3.4-3.6 μ\mum IR features, sometimes observed in astronomical sources (see e.g. Li, 2020, for a recent review); e.g. in proto-PNe, where satellite features from 3.4 to 3.6 μ\mum, accompanying the 3.3 μ\mum feature, have been observed (Geballe et al., 1992; Hrivnak et al., 2007; Materese et al., 2017). Such 3-5 μ\mum observations, however, have been carried out with space- and ground-based telescopes less sensitive than Spitzer, which otherwise didn’t cover this spectral range. In particular, Geballe et al. (1992) analyzed several objects in the short transition phase from AGB stars to PNe (i.e. proto-PNe) and observed remarkably strong 3.4-3.5 μm\mu m emission relative to the usually dominant 3.3 μm\mu m feature. They concluded that the observed 3.4-3.5 μm\mu m features cannot be explained by the presence of stretch overtones of CH produced by PAHs and they suggested that they should be assigned to the stretching modes of aliphatic CH2 and CH3 groups. On the other hand, Materese et al. (2017) suggested that the correlation between the abnormally intense 3.4-3.5 μm\mu m and 6.9 μm\mu m features observed in this kind of post-AGB objects could be due to the presence of superhydrogenated PAHs (Hn-PAHs), containing aliphatic CH2 and CH3 moieties. The carbon chemistry in space is dominated by the 3.3 and 3.4 μ\mum features or aromatic-aliphatic dichotomy, with aromatics (PAHs; 3.3 μ\mum) generally believed to be a main component of carbon in space. This popular perception, however, could be challenged in the JWST era. Thanks to the exceptional sensitivity and spectral resolution of recent JWST observations, the 3.4 μ\mum features are routinely detected in very different types of astronomical sources; many of them with a no significant UV radiation field like proto-PNe, trojans, active galactic nuclei, etc. (see e.g. Lai et al., 2023; Wong et al., 2024) but also in more evolved and old PNe with a much stronger UV radiation field (Clark et al., 2025).

Considerable efforts have been previously made to understand the origin (aliphatics, PAHs, etc.) of these unusual \sim3.4-3.6 μ\mum features. The PAH molecules can barely fit these features by using a very large family of irregular structures, even assuming the effects of anharmonic couplings due to high temperatures (see e.g. Joblin et al., 1995; Candian et al., 2012; Bauschlicher et al., 2009). A more recent theoretical work suggests that IR features beyond 3.4 μ\mum, which is the archetypal aliphatic CH stretching mode, cannot be assigned to pure CH bonds (Sadjadi et al., 2017). According to Sadjadi et al. (2017) those features should probably involve another element (e.g. N, S, O) besides C and H, contaminating the CH stretching. On the other hand, hydrogenated fullerenes (fulleranes) have also been suggested to contribute in the \sim3.4-3.6 μ\mum region via their CH bonds (Iglesias-Groth et al., 2012; Zhang et al., 2017). This is reinforced by recent abundance predictions, suggesting that C60H+ should be abundant in the diffuse ISM (Abbink et al., 2024). Also, accurate experimental studies on fulleranes (C60H+\rm C_{60}H^{+}, C70H+\rm C_{70}H^{+}; Palotás et al., 2020; Finazzi et al., 2024) and oxidized fullerenes (C60O+\rm C_{60}O^{+}, C60OH+\rm C_{60}OH^{+}; Palotás et al., 2022) demonstrate that these species have a common IR emission region around \sim6.2-8.3 μ\mum, corresponding to CC vibrations. Thus, since all these fullerene compounds contain CH and CC stretching, their IR spectra could be entangled with those from PAH-C60 adducts in several spectral regions. However, a possible way to tackle this, in order to resolve or distinguish their contribution in astronomical IR spectra, is via the CH/CC stretching band ratio. For example, Finazzi et al. (2024) find that the CH/CC band ratio for fulleranes is quite small due to their low IR intensity of CH stretching; note that CH stretching is not present for oxidized fullerenes (Palotás et al., 2022). In contrast, the PAH-C60\rm C_{60} adducts studied here display a CH/CC band ratio larger than fulleranes (see Figure 3). The potential use of the CH/CC band ratio to distinguish between different H-containing fullerene derivatives would require sensitive astronomical data like those from the JWST; i.e. covering the entire IR range (\sim3-25 μ\mum) of molecular vibrations.

Clearly, there is no consensus yet about the carrier/s of the \sim3.4-3.6 μ\mum features (and similarly occurs for the 6.9 μ\mum feature generally observed in conjunction with them). Our experimental-theoretical study of PAH-C60 adducts suggests other species that may contribute at these wavelengths (see Figures 3, 6 and 13). Remarkably, the PAH-C60 adducts display strong \sim3.4-3.6 μ\mum features without the presence of aliphatic CH bonds (with the exception of the indene derivatives), due to the interplay between C60 and small PAHs (see Figures in the Appendix F). In contrast to Sadjadi et al. (2017) and Zhang et al. (2017), which use a diverse set of model structures, our set of PAH-C60 adducts are rather structurally simple. It is thus likely that the contamination of PAH-C60 adducts by aliphatic chains within the PAH structure would increase, even more, the intensity of the \sim3.4-3.6 μ\mum IR features. Future efforts should enlarge the PAHs diversity, from linear to branched or even helical structures.

Summary

In summary, we have computed IR simulated spectra of multiple PAH-C60 adducts based on experimental data combined with quantum-chemical calculations. According to the experimental data and theoretical results the formation of these new species would almost erase the presence of most of the characteristic IR features from the pristine PAHs and some from isolated C60. Furthermore, we have very roughly estimated the possible abundances of PAH-C60 adducts in C60-rich astrophysical environments by applying the state-of-the-art kinetics data, hopefully building more reliable and global (or weighted) simulated IR spectra. Consequently, the abundance-weighted spectra display a series of new relevant IR features like a broad \sim3.4-3.6 μ\mum feature, a richer 6-10 µm spectral region, a strong modification within the 12-16 μ\mum range together with emission contribution to the 7.0 and 18.9 μ\mum C60 features. Surprisingly, the PAH-C60 adducts presented here, with almost no CH aliphatic in their structure, display strong \sim3.4-3.6 μ\mum features and could be potential carriers of this kind of emission in astronomical environments; especially in those where C60 is known to be abundant. Unfortunately, previous IR observations of C60-rich astrophysical environments were done by Spitzer with no access to the 3-5 μ\mum spectral range. However, such observations could be made at high-sensitivity by the JWST; with its spectral coverage well below 5 μ\mum. This along with high-quality experimental IR spectra of diverse PAH-C60 adducts (e.g., including those formed by C60 and PAHs of different size and chemical nature) will give more insights on the possible IR emission contribution of PAH-C60 adducts to the 3-4 μ\mum spectral region.

Data availability

All the scaled harmonic DFT spectra presented herein, as well as the IR cross sections used to construct them, are available at https://doi.org/10.5281/zenodo.17485506.

Acknowledgements.
We acknowledge the support from the State Research Agency (AEI) of the Ministry of Science, Innovation and Universities (MICIU) of the Government of Spain, and the European Regional Development Fund (ERDF), under grants PID2020-115758GB-I00/AEI/10.13039/501100011033, PID2022-136970NB-I00/AEI/10.13039/501100011033 and PID2023-147325NB-I00/AEI/10.13039/501100011033. The authors also express their gratitude to the Deanship of Graduate Studies and Scientific Research at Taif University for their financial support of this work. This publication is based upon work from COST Action CA21126 - Carbon molecular nanostructures in space (NanoSpace), supported by COST (European Cooperation in Science and Technology). R.B and B.K also acknowledge the generous allocation of computer time at LaPalma-IAC Supercomputer and CITIC Universidade da Coruña. This article used flash storage and CPU/GPU computing resources as Indefeasible Computer Rights (ICRs) being commissioned at the ASTRO POC project that Light Bridges will operate in the Island of Tenerife, Canary Islands (Spain). The ICRs used for this research were provided by Light Bridges in cooperation with Hewlett Packard Enterprise (HPE) and VAST DATA. The authors wish to acknowledge the contribution of the IAC High-Performance Computing support team and hardware facilities to the results of this research. BK is grateful to Science by Women program for funding a six-months visiting senior research fellowship to IAC. The authors acknowledge Deanship of Scientific Research, Taif University for funding this work.

References

  • Abbink et al. (2024) Abbink, D., Foing, Bernard, & Ehrenfreund, Pascale. 2024, A&A, 684, A165
  • Allamandola et al. (1989) Allamandola, L., Tielens, A., & Barker, J. 1989, ApJS, 71, 733
  • Arun et al. (2023) Arun, R., Mathew, B., Manoj, P., et al. 2023, MNRAS, 523, 1601
  • Avery et al. (1982) Avery, L. W., MacLeod, J. M., & Broten, N. W. 1982, ApJ, 254, 116
  • Barham et al. (2018) Barham, J. P., Tanaka, S., Koyama, E., et al. 2018, J. Org. Chem., 83, 4348
  • Barzaga et al. (2023a) Barzaga, R., García-Hernández, D. A., Díaz-Tendero, S., et al. 2023a, ApJ, 942, 5
  • Barzaga et al. (2023b) Barzaga, R., García-Hernández, D. A., Díaz-Tendero, S., et al. 2023b, ApJS, 269, 26
  • Bauschlicher et al. (2009) Bauschlicher, Jr., C. W., Peeters, E., & Allamandola, L. J. 2009, ApJ, 697, 311
  • Bauschlicher et al. (2018) Bauschlicher, Charles W., J., Ricca, A., Boersma, C., & Allamandola, L. J. 2018, ApJS, 234, 32
  • Beheshtian et al. (2012) Beheshtian, J., Peyghan, A. A., & Bagheri, Z. 2012, ApSS, 258, 8980
  • Berné et al. (2017) Berné, O., Cox, N. L. J., Mulas, G., & Joblin, C. 2017, A&A, 605, L1
  • Boersma et al. (2014) Boersma, C., Bauschlicher, C. W., Ricca, A., et al. 2014, ApJS, 211, 8
  • Brieva et al. (2016) Brieva, A. C., Gredel, R., Jäger, C., Huisken, F., & Henning, T. 2016, ApJ, 826, 122
  • Burkhardt et al. (2021) Burkhardt, A. M., Lee, K. L. K., Changala, P. B., et al. 2021, ApJL, 913, L18
  • Cami et al. (2010) Cami, J., Bernard-Salas, J., Peeters, E., & Malek, S. E. 2010, Science, 329, 1180
  • Cami et al. (2011) Cami, J., Bernard-Salas, J., Peeters, E., & Malek, S. E. 2011, Proc. IAU, 7, 216–227
  • Cami et al. (2018) Cami, J., Peeters, E., Bernard-Salas, J., Doppmann, G., & De Buizer, J. 2018, Galaxies, 6, 101
  • Campbell et al. (2015) Campbell, E., Holz, M., Gerlich, D., & Maier, J. 2015, Nature, 523, 322
  • Candian et al. (2012) Candian, A., Kerr, T. H., Song, I.-O., McCombie, J., & Sarre, P. J. 2012, MNRAS, 426, 389
  • Cataldo et al. (2014) Cataldo, F., García-Hernández, D. A., & Manchado, A. 2014, FNCN, 22, 565
  • Cernicharo et al. (2021) Cernicharo, J., Agúndez, M., Cabezas, C., et al. 2021, A&A, 649, L15
  • Cernicharo et al. (2024) Cernicharo, J., Cabezas, C., Fuentetaja, R., et al. 2024, A&A, 690, L13
  • Chown et al. (2024) Chown, R., Sidhu, Ameek, Peeters, Els, et al. 2024, A&A, 685, A75
  • Clark et al. (2025) Clark, N., Peeters, E., Cox, N. L. J., et al. 2025, MNRAS, 540, 1984
  • Díaz-Luis et al. (2015) Díaz-Luis, J. J., García-Hernández, D. A., Kameswara Rao, N., Manchado, A., & Cataldo, F. 2015, A&A, 573, A97
  • Díaz-Luis et al. (2018) Díaz-Luis, J. J., García-Hernández, D. A., Manchado, A., et al. 2018, AJ, 155, 105
  • Dunk et al. (2013) Dunk, P. W., Adjizian, J.-J., Kaiser, N. K., et al. 2013, PNAS, 110, 18081
  • Dunn & Rashed (2018) Dunn, J. L. & Rashed, E. 2018, J. Phys. Conf. Ser., 1148, 012003
  • Ehrlich et al. (2013) Ehrlich, S., Moellmann, J., & Grimme, S. 2013, AcChR, 46, 916
  • El Bakouri et al. (2018) El Bakouri, O., Garcia-Borràs, M., Girón, R. M., et al. 2018, PCCP, 20, 11577
  • Finazzi et al. (2024) Finazzi, L., Esposito, V. J., Palotás, J., et al. 2024, ApJ, 971, 168
  • Frisch et al. (2016) Frisch, M. J., Trucks, G. W., Schlegel, H. B., et al. 2016, Gaussian 16 Revision C.01, Gaussian Inc. Wallingford CT
  • Fu & Hopkins (2018) Fu, W. & Hopkins, W. S. 2018, JPCA, 122, 167
  • García-Hernández & Díaz-Luis (2013) García-Hernández, D. A. & Díaz-Luis, J. J. 2013, A&A, 550, L6
  • García-Hernández et al. (2011a) García-Hernández, D. A., Kameswara Rao, N., & Lambert, D. L. 2011a, ApJ, 729, 126
  • García-Hernández et al. (2012a) García-Hernández, D. A., Kameswara Rao, N., & Lambert, D. L. 2012a, ApJ, 759, L21
  • García-Hernández et al. (2010) García-Hernández, D. A., Manchado, A., García-Lario, P., et al. 2010, ApJL, 724, L39
  • García-Hernández et al. (2012b) García-Hernández, D. A., Villaver, E., García-Lario, P., et al. 2012b, ApJ, 760, 107
  • García-Hernández et al. (2013) García-Hernández, D. A., Cataldo, F., & Manchado, A. 2013, MNRAS, 434, 415
  • García-Hernández et al. (2011b) García-Hernández, D. A., Iglesias-Groth, S., Acosta-Pulido, J. A., et al. 2011b, ApJL, 737, L30
  • Geballe et al. (1992) Geballe, T., Tielens, A., Kwok, S., & Hrivnak, B. 1992, ApJL, 387, L89
  • Gerlich et al. (2018) Gerlich, D., Jašík, J., Strelnikov, D. V., & Roithová, J. 2018, ApJ, 864, 62
  • Gold et al. (2024) Gold, K. R., Schmidt, D. R., & Ziurys, L. M. 2024, ApJ, 976, 196
  • Gredel et al. (2011) Gredel, R., Carpentier, Y., Rouillé, G., et al. 2011, A&A, 530, A26
  • Grimme et al. (2010) Grimme, S., Antony, J., Ehrlich, S., & Krieg, H. 2010, JChPh, 132, 154104
  • He et al. (2010) He, Y., Chen, H.-Y., Hou, J., & Li, Y. 2010, J. Am. Chem. Soc., 132, 1377
  • Hirsch (1995) Hirsch, A. 1995, Synthesis, 1995, 895
  • Hrivnak et al. (2007) Hrivnak, B. J., Geballe, T. R., & Kwok, S. 2007, ApJ, 662, 1059
  • Iglesias-Groth (2019) Iglesias-Groth, S. 2019, MNRAS, 489, 1509
  • Iglesias-Groth et al. (2012) Iglesias-Groth, S., García-Hernández, D. A., Cataldo, F., & Manchado, A. 2012, MNRAS, 423, 2868
  • Jin et al. (2019) Jin, H., Xing, L., Hao, J., et al. 2019, CoFl, 206, 1
  • Joblin et al. (1995) Joblin, C., Boissel, P., Leger, A., d’Hendecourt, L., & Defourneau, D. 1995, A&A, 299, 835
  • Katsyuba et al. (2013) Katsyuba, S. A., Zvereva, E. E., & Burganov, T. I. 2013, JPCA, 117, 6664
  • Kempkes et al. (2019) Kempkes, L. J. M., Martens, J., Berden, G., Houthuijs, K. J., & Oomens, J. 2019, FaDi, 217, 434
  • Khodam Hazrati & L. Hadipour (2016) Khodam Hazrati, M. & L. Hadipour, N. 2016, Comput. Theor. Chem., 1098, 63
  • Krusic et al. (1991) Krusic, P. J., Wasserman, E., Keizer, P. N., Morton, J. R., & Preston, K. F. 1991, Science, 254, 1183
  • Lai et al. (2023) Lai, T. S.-Y., Armus, L., Bianchin, M., et al. 2023, ApJL, 957, L26
  • Li (2020) Li, A. 2020, NatAs, 4, 339
  • Li (2022) Li, F.-B. 2022, in Handbook of Fullerene Science and Technology (Springer), 273–312
  • Lifshitz et al. (2004) Lifshitz, A., Tamburu, C., Suslensky, A., & Dubnikova, F. 2004, JPCA, 108, 3430
  • Manchado et al. (2015) Manchado, A., Stanghellini, L., Villaver, E., et al. 2015, ApJ, 808, 115
  • Mas-Torrent et al. (2002) Mas-Torrent, M., Rodríguez-Mias, R. A., Solà, M., et al. 2002, J. Org. Chem., 67, 566
  • Materese et al. (2017) Materese, C. K., Bregman, J. D., & Sandford, S. A. 2017, ApJ, 850, 165
  • Mattioda et al. (2020) Mattioda, A. L., Hudgins, D. M., Boersma, C., et al. 2020, ApJS, 251, 22
  • McGuire et al. (2021) McGuire, B. A., Loomis, R. A., Burkhardt, A. M., et al. 2021, Science, 371, 1265
  • Mulholland et al. (2000) Mulholland, J. A., Lu, M., & Kim, D.-H. 2000, PComI, 28, 2593
  • Müller et al. (2020) Müller, F., Stückrath, J. B., Bischoff, F. A., et al. 2020, J. Am. Chem. Soc., 142, 18050
  • Otsuka et al. (2013) Otsuka, M., Kemper, F., Hyung, S., et al. 2013, ApJ, 764, 77
  • Palotás et al. (2020) Palotás, J., Martens, J., Berden, G., & Oomens, J. 2020, NatAs, 4, 240
  • Palotás et al. (2022) Palotás, J., Martens, J., Berden, G., & Oomens, J. 2022, JPCA, 126, 2928, pMID: 35533303
  • Peeters et al. (2024) Peeters, E., Habart, Emilie, Berné, Olivier, et al. 2024, A&A, 685, A74
  • Peeters et al. (2021) Peeters, E., Mackie, C., Candian, A., & Tielens, A. G. G. M. 2021, AcChR, 54, 1921
  • Petersson & Al‐Laham (1991) Petersson, G. A. & Al‐Laham, M. A. 1991, JChPh, 94, 6081
  • Petersson et al. (1988) Petersson, G. A., Bennett, A., Tensfeldt, T. G., et al. 1988, JChPh, 89, 2193
  • Polfer (2011) Polfer, N. C. 2011, CSRev, 40, 2211
  • Pousse et al. (2010) Pousse, E., Tian, Z., Glaude, P., Fournet, R., & Battin-Leclerc, F. 2010, CoFl, 157, 1236
  • Puget & Léger (1989) Puget, J. L. & Léger, A. 1989, ARA&A, 27, 161
  • Pujals et al. (2022) Pujals, M., Pèlachs, T., Fuertes-Espinosa, C., et al. 2022, CRPS, 3, 100992
  • Roberts et al. (2012) Roberts, K. R. G., Smith, K. T., & Sarre, P. J. 2012, MNRAS, 421, 3277
  • Rodrigues et al. (2023) Rodrigues, D. J., Pina, I. B., Santos, L. M., & Lima, C. F. 2023, DRM, 136, 110031
  • Rodrigues et al. (2022) Rodrigues, D. J. L., Santos, L. M. N. B. F., Melo, A., & Lima, C. F. R. A. C. 2022, Organics, 3, 364
  • Roithová et al. (2016) Roithová, J., Gray, A., Andris, E., Jašík, J., & Gerlich, D. 2016, AcChR, 49, 223
  • Rosenberg et al. (2014) Rosenberg, M. J. F., Berné, O., & Boersma, C. 2014, A&A, 566, L4
  • Rouillé et al. (2021) Rouillé, G., Krasnokutski, S. A., & Carpentier, Y. 2021, A&A, 656, A100
  • Sabirov et al. (2016) Sabirov, D. S., Terentyev, A. O., & Cataldo, F. 2016, Comput. Theor. Chem., 1081, 44
  • Sadjadi et al. (2017) Sadjadi, S., Zhang, Y., & Kwok, S. 2017, ApJ, 845, 123
  • Sarova & Berberan-Santos (2004) Sarova, G. H. & Berberan-Santos, M. N. 2004, CPL, 397, 402
  • Sellgren et al. (2010) Sellgren, K., Werner, M. W., Ingalls, J. G., et al. 2010, ApJL, 722, L54
  • Stephens et al. (1994) Stephens, P. J., Devlin, F. J., Chabalowski, C. F., & Frisch, M. J. 1994, JPhCh, 98, 11623
  • Stockett et al. (2025) Stockett, M. H., Subramani, A., Liu, C., et al. 2025, The Journal of Chemical Physics, 162, 184306
  • Wenzel et al. (2025a) Wenzel, G., Gong, S., Xue, C., et al. 2025a, ApJ, 984, L36
  • Wenzel et al. (2025b) Wenzel, G., Speak, T. H., Changala, P. B., et al. 2025b, Nature Astronomy, 9, 262
  • Wesson et al. (2024) Wesson, R., Matsuura, M., Zijlstra, A. A., et al. 2024, MNRAS, 528, 3392
  • Wong et al. (2024) Wong, I., Brown, M. E., Emery, J. P., et al. 2024, PSJ, 5, 87
  • Wu et al. (2024) Wu, Y., Hu, X., Zhen, J., & Yang, X. 2024, MNRAS, 531, 682
  • Xu et al. (2023) Xu, J., Li, A., Li, X., & Hou, G.-L. 2023, MNRAS, 525, 3061
  • Zapata Trujillo & McKemmish (2022a) Zapata Trujillo, J. C. & McKemmish, L. K. 2022a, WIREs Comput. Mol. Sci., 12, e1584
  • Zapata Trujillo & McKemmish (2022b) Zapata Trujillo, J. C. & McKemmish, L. K. 2022b, JPCA, 126, 4100
  • Zapata Trujillo & McKemmish (2023) Zapata Trujillo, J. C. & McKemmish, L. K. 2023, JPCA, 127, 1715
  • Zhang & Kwok (2015) Zhang, Y. & Kwok, S. 2015, ApJ, 798, 37
  • Zhang et al. (2017) Zhang, Y., Sadjadi, S., Hsia, C.-H., & Kwok, S. 2017, ApJ, 845, 76
  • Zhen et al. (2019a) Zhen, J., Zhang, W., Yang, Y., & Zhu, Q. 2019a, MNRAS, 490, 3498
  • Zhen et al. (2019b) Zhen, J., Zhang, W., Yang, Y., Zhu, Q., & Tielens, A. G. G. M. 2019b, ApJ, 887, 70

Appendix A Empirical-corrected scaling factors

The scaling factors reported in Table Scaling factors were empirically corrected in the \sim3.3 μ\mum region to improve the accuracy with respect to the experimental spectra (see the corresponding Section and related discussion); Table A shows the final values.

Table 3: Frequency-range-specific scaling factors for the PAHs in Figure 1.151515The scaling factors from our 6-31+G(d)/B3LYP+GD3 (SFD33{}_{3}^{\mathrm{D3}}) method and their corresponding empirical-corrected (SF3empD3\mathrm{SF_{3\,{\mathrm{emp}}}^{D3}}) are reported. The legends remain as in Table Scaling factors.
PAH SF3D3\mathrm{SF_{3}^{D3}} SF3empD3\mathrm{SF_{3\,{\mathrm{emp}}}^{D3}}
0.9852(1) 0.9852
In 0.9675(2) 0.9675
0.9614(3) 0.9090
\arrayrulecolorgray 0.9838 0.9838
Iyl 0.9664 0.9664
0.9627 0.9103
0.9807 0.9807
An 0.9711 0.9711
0.9610 0.9086
0.9818 0.9818
Tn 0.9679 0.9679
0.9614 0.9090
0.9848 0.9848
Pn 0.9648 0.9648
0.9618 0.9094
\arrayrulecolorblack

Appendix B Bis-adducts

Refer to caption
Figure 7: Left: Models of the PAH-C60 bis-adducts after geometry optimization for the PAHs 2H-indene (2InC60, a-c), indenyl (2IylC60, d-g), anthracene (2AnC60, h), and tetracene (2TnC60, i). In the case that the bis-adduct has different configurations (isomers), the most stable structure is highlighted with a red asterisk. Right: The theoretical IR spectra corresponding to the bis-adducts structures a-i presented in the left panel. In all panels, the red dashed spectra correspond to the mono-adduct analogous, while the blue dashed lines mark the four strongest C60 features (\sim7.0, 8.5, 17.4 and 18.9 μ\mum). Note that the convolution parameters are the same as in Figures 4 and 5.

Multiple configurations (isomers) can be generated for the bis-adducts according to the arrangement of the PAHs bonded to the C60 cage. The isomers are molecular structures defined by the same chemical formula but with different spatial distribution. The number of possible isomers (or structures) is extremely large because it depends on the degrees of freedoms for the PAHs and C60. Thus, in order to reduce the number of possible models we have used the stability and isomerization information available in the literature. Previous theoretical works have studied the stability and isomerization of anthracene bis-adducts, allowing us to discriminate between the different models and select the most stable structures (Sabirov et al. 2016; Rodrigues et al. 2022). This isomerization information was also used to built the rest of the bis-adduct models, facilitating our theoretical predictions by significantly reducing the amount of models under consideration. Figure 7 displays the bis-adducts structures considered together with their corresponding most stable isomers. For the anthracene (An) case, only one structure is shown; i.e. the most stable configuration as demonstrated by previous and extensive theoretical studies (Sabirov et al. 2016; Rodrigues et al. 2022). Furthermore, as already mentioned in Sect. Scaling factors validation, the anthracene bis-adduct (2AnC60) has been obtained precisely using a novel experimental synthetic method (Pujals et al. 2022) and used for the validation of scaling factors.

The spectra displayed in Figures 7a-c, which correspond to the bis-adducts of 2H-indene, do not exhibit significant differences between them; the only difference is the isomer c (Figure 7c). Actually, it is quite difficult to distinguish the IR emission produced by these isomers of 2H-indene bis-adducts. The different configuration of the PAHs with respect to the C60 cage neither the increment of PAH units implicate a noticeable change of the IR spectra, in comparison to the mono-adduct analogous (1InC60{}^{{}^{\prime}}, see Figure 4). The most prominent change is seen for the C-H stretching feature at \sim3.3-3.6 μ\mum, which increases its intensity, reflecting the increment of CH units from mono- to bis-adducts.

On the contrary, the IR spectra produced by the isomers of indenyl bis-adducts display noticeable differences between them and also when compared to their mono-adduct analogous (Figures 7d-g and 11); with the clear presence of distinctive IR spectral features depending on the isomer. The indenyl isomers have these specific features mainly due to the change of the PAH binding to the C60 cage. The PAH can bind C60 with one or two C-C bonds depending on the structure, while for 2H-indene the binding is always the same. Models like those in Figures 7d-e denote single C-C bonds between the PAH and C60. Both exhibit subtle differences in terms of intensity, but they possess four distinctive features: (i) a broad plateau feature from 13 to 14.7 μ\mum; (ii) a clear discrete \sim12.4 μ\mum feature accompanied with well-defined peaks around it; (iii) a broad 6.7-7.5 μ\mum band; and (iv) a strong \sim3.3-3.5 μ\mum feature accompanied by a satellite feature at \sim3.6 μ\mum.

A higher stability is not always an indication of more detectable features, and the complexity of the IR spectra depends more on the charge reordering and change in the dipole moment induced by the PAH-C60 bonding. The former is clearly understandable observing the spectra in Figures 7f-g, corresponding to the less stable isomers of indenyl bis-adducts. The specific features observed in these spectra are quite noticeable in terms of intensity; for instance, Figure 7f shows an intense narrow signal at \sim14.2 μ\mum accompanied by broader feature centered at \sim13.3 μ\mum. Under this configuration f, the IR spectrum shows a resemblance with its mono-adduct analogous (1IylC60{}^{{}^{\prime}}), maintaining almost the same features. Such behavior indicates that binding two indenyl molecules to C60 following this geometry almost does not destroy the symmetry of the vibrations. Interestingly, the model in Figure 7g results in weakly bonded indenyl molecules, but it shows the richer IR spectrum of all of the bis-adduct models. Multiple features (e.g. at \sim13.0, 14.6, 15.0 and 20.6 μ\mum) surpass in intensity the \sim3.3-3.6 μ\mum band, which is usually the most intense signal in the PAH-C60 adducts described so far; where most of the vibrations, implying C-C stretching and C-H out-of-plane, close to the PAH-C60 bond create a distortion similarly to a Jahn-Teller effect (Dunn & Rashed 2018)161616In the case of the model in Figure 7g the point group is defined as C2h, which is one of the symmetry subgroups producing a Jahn-Teller distortion in C60 derivatives.. Seemingly, these C-H out-of-plane vibrations are IR active in the bis-adduct due to the higher planarity and aromaticity of the indenyl structure. This combined to the configuration in Figure 7g produce the strong features at \sim14.6, 15.0 and 20.6 μ\mum. However, the analogous bis-adduct for 2H-indene with a similar configuration (Figure 7c) does not display such C-H out-of-plane active vibrations because the presence of CH2 in 2H-indene reduces the planarity and aromaticity of the molecular structure. Nevertheless, it is likely that the indenyl bis-adduct (Figure 7g) should be a species with a rather short-lifetime as a consequence of the weak bonds between the PAHs and C60. Finally, the IR features seen in the simulated spectra of the models for the An and Tn bis-adducts (Figures 7h-i) are characterized only by an increment of CH stretching (3.3-3.6 μ\mum) due to the amount of added bonds.

Appendix C Tris-adducts

Refer to caption
Figure 8: Left: Models of the PAH-C60 tris-adducts after geometry optimization for the PAHs 2H-indene (3InC60, a-b), indenyl (3IylC60, c-e) and anthracene (3AnC60, f-g). In the case that the tris-adduct has different configurations (isomers), the most stable structure is highlighted with a red asterisk. Right: The theoretical IR spectra corresponding to the tris-adducts structures a-g presented in the left panel. In all panels, the red dashed spectra correspond to the bis-adduct analogous, while the blue dashed lines mark the four strongest C60 features (\sim7.0, 8.5, 17.4 and 18.9 μ\mum). Note that the convolution parameters are the same as in Figures 4 and 5.

In order to build the models of the PAH-C60 tris-adducts, we have selected the corresponding bis-adduct from the previous Appendix B and added another PAH unit to the structure. Arguably, this could bias the final structures since we did not explore multiple possible isomers, but it was our preferred way to screening the number of structures; which again, would be huge due to the high number of degrees of freedom. Figure 8 shows the models and IR spectra of the tris-adducts for In, Ilyl and An. As we have mentioned before, the C60 tris-adducts models for tetracene and pentacene were not built because of the computational limitations.

Figure 8 shows non-significant differences in the IR spectra of the same PAH-C60 tris-adduct; except, again for the C60 tris-adducts with indenyl (3IylC60, panels c-e in Figure 8). In particular, for configuration e where the spectral change is due to the single C-C bond between the PAH and C60. Furthermore, the 3IylC60 spectra are still very similar to the ones from their mono-adducts analogous; although, the intensity varies for the most important vibrations that depend on the number of CH bonds (see also Figure 4). The forces exerted over the C60 cage by the binding of three PAH units create a highly symmetric carbon cage, which is equivalent to the mono-adducts and that is reflected in the spectra similarities. Thus, according to our theoretical predictions it should be very difficult to distinguish, in terms of pure IR emission features, mono- from tris-adducts.

In summary, the main characteristics of the C60 tris-adducts spectra are: (i) the more noticeable broad feature, centered at \sim3.5 μ\mum, for 2H-indene (Figures 8a-b); and (ii) the presence of a weaker red-shifted feature at \sim3.6 μ\mum (Figures 8e-g), which denotes the increment of non-equivalent CH bonds by symmetry directly connected to C60. Such features were difficult to distinguish in the case of mono- and bis-adducts due to the lower number of CH bonds implicated. Both features are related to CH bonds but imply different chemical environments; for 2H-indene tris-adducts, they are a consequence of the contribution of more CH2 bonds, while for the indenyl and anthracene cases, they are related to the CH closest to the C-C bonded to C60. Finally, it is worth noting that in the case of indenyl, the feature at \sim3.6 μ\mum only appears in the model e with a C-C single bond to C60 (see Figures 8c-e and Appendix F for a comparison).

Appendix D Abundance estimation from kinetics data for the construction of abundance-weighted IR spectra of PAH-C60 adducts

We have estimated the expected abundances of PAH-C60 adducts from the kinetics data previously reported. For the case of indenyl (C9H7\rm C_{9}H_{7}), the data extracted from indene (C9H8\rm C_{9}H_{8}) pyrolysis (Jin et al. 2019; Lifshitz et al. 2004; Pousse et al. 2010) has been used. Note that indene pyrolysis is also connected to the anthracene formation, through the following mechanism in the presence of cyclopentadiene (C5H6\rm C_{5}H_{6}) (Mulholland et al. 2000):

C9H8+C5H62HC14H12(C9H7C5H5)2HC14H10\rm C_{9}H_{8}+C_{5}H_{6}\xrightarrow{-2H}C_{14}H_{12}\mkern 3.0mu(C_{9}H_{7}-C_{5}H_{5})\xrightarrow{-2H}C_{14}H_{10}

Therefore, with the exception of pentacene for which not kinetic data is available, it is possible to roughly estimate the expected abundances of PAH-C60 adducts by using the corresponding equations, knowing that they follow a second order rate law:

[A]=[A]01+kt[A]0\rm[A]=\frac{[A]_{0}}{1+kt\cdot[A]_{0}} (6)

where k is the rate constant, [A] is the concentration of the reactant at a given time t and [A]0 is the initial concentration of the reactant. In order to apply Equation 6 we have assumed the following assumptions:

  • The PAH-C60 abundances depend only on the PAH abundance.

  • The abundances of the Iyl-C60 and An-C60 adducts depend only on the indene abundance, following the pyrolysis mechanism mentioned above.

  • The initial abundance [A]0 of indene is obtained from the column density estimates towards the cold dark cloud TMC-1 reported by Cernicharo et al. (2021) and Burkhardt et al. (2021), which are 1.6 ×\times 1013 and 0.96 ×\times 1013 cm-2, respectively171717The column density is transformed to abundance (or concentration, in mol/cm3) by assuming the TMC-1 geometrical model of Avery et al. (1982) and a distance of 140 pc; i.e. basically along a path of 0.14 pc..

  • The transition from mono- to bis-adducts follows the same kinetics as from bis- to tris-adducts, using as reference the solid state kinetics of InC60 by Rodrigues et al. (2023).

According to the above-mentioned principles we have determined the PAH abundances of indenyl, anthracene and tetracene from indene in order to estimate the corresponding PAH-C60 adducts abundances. For this purpose, we refer again to Equation 6, where it can be noticed that the timescale (t) has to be introduced. This timescale can be first computed from the indene pyrolisis (Mulholland et al. 2000) using the half-life time expression:

t1/2=1k[A]0\rm t_{1/2}=\frac{1}{k[A]_{0}} (7)

we remind that in Equation 7 [A]0 is the indene abundance from Cernicharo et al. (2021) and Burkhardt et al. (2021). Using the value in t1/2\rm t_{1/2} and the experimental constant rates (k) of the PAH-C60 adducts we can determine [A] and [A]0 for indene and anthracene. For indenyl, we directly used the result obtained from the indene pyrolisis, while for the case of the tetracene adducts we have used the relation of the constant rates with anthracene (see Tables 1 and 2 in Sarova & Berberan-Santos 2004):

kTnC60=1.94102kAnC60\rm k_{TnC_{60}}=1.94\cdot 10^{2}k_{AnC_{60}} (8)

It should be noted that [A] and [A]0 the reactants, and thus the abundances of products must be derived from [A]prod=[A]0-[A], which describes how much of the reactant has become a product. Table D displays the estimated expected abundance for each PAH-C60 adduct under study. Clearly, from the listed values, the abundances of anthracene and tetracene are negligible compared to the ones of indene and indenyl.

Table 4: Estimated abundances (concentrations) for the different PAH-C60 adducts.181818[A]Cerprod{}_{prod}^{\rm Cer} and [A]Bukprod{}_{prod}^{\rm Buk} correspond to the abundances determined using the initial abundance [A]0 of indene as estimated from the column density values by Cernicharo et al. (2021) and Burkhardt et al. (2021), respectively (see the text for more details).
PAH-C60 Adduct [A]Cerprod{}_{prod}^{\rm Cer} [A]Bukprod{}_{prod}^{\rm Buk}
mol/cm3 mol/cm3
Mono 2H-indene (In) 6.151029\cdot 10^{-29} 3.691029\cdot 10^{-29}
indenyl (Iyl) 3.071029\cdot 10^{-29} 1.841029\cdot 10^{-29}
anthracene (An) 2.60 1032\cdot 10^{-32} 3.11 1032\cdot 10^{-32}
tetracene (Tn) 1.00 1037\cdot 10^{-37} 1.00 1037\cdot 10^{-37}
\arrayrulecolorgray    Bis 2H-indene (In) 3.551029\cdot 10^{-29} 2.131029\cdot 10^{-29}
indenyl (Iyl) 1.301029\cdot 10^{-29} 7.781030\cdot 10^{-30}
anthracene (An) 1.101032\cdot 10^{-32} 1.311032\cdot 10^{-32}
tetracene (Tn) 4.221038\cdot 10^{-38} 4.221038\cdot 10^{-38}
\arrayrulecolorgray    Tris 2H-indene (In) 3.551029\cdot 10^{-29} 2.131029\cdot 10^{-29}
indenyl (Iyl) 1.301029\cdot 10^{-29} 7.781030\cdot 10^{-30}
anthracene (An) 1.101032\cdot 10^{-32} 1.311032\cdot 10^{-32}
\arrayrulecolorblack

On the other hand, we also consider the probability of formation according to the thermodynamic stability of PAH-C60 adducts and its influence on the IR spectra. For this, a Boltzmann distribution (equation 9 below) has been used to weight the different regioisomers for the bis- and tris-adducts (see Table Thermodynamic stability), which is described as:

Pi=eΔHfi/kTj=1MeΔHfj/kT\rm P_{i}=\frac{e^{-\Delta H^{i}_{f}/kT}}{\displaystyle\sum^{M}_{j=1}e^{-\Delta H^{j}_{f}/kT}} (9)

with ΔHfi\rm\Delta H^{i}_{f} being the Enthalpy of formation of the bis- or tris-adduct, M the maximum number of regioisomers, T the temperature and k the Boltzmann constant. Note that Pi=1\rm P_{i}=1 for the 2AnC60 and 2TnC60 adducts because no regioisomers were considered. We also note that ΔGf=ΔHfTΔSf\rm\Delta G_{f}=\Delta H_{f}-T\Delta S_{f} is the expression that connects these ΔHf\rm\Delta H_{f} values with the ΔGf\rm\Delta G_{f} ones in Table Thermodynamic stability; both were calculated at 300 K and are directly extracted from our quantum-chemistry calculations. By using the quantities obtained so far, the IR intensity for each PAH-C60 adduct is thus computed by the following expression:

IPi[A]𝑝𝑟𝑜𝑑μx\rm I\propto P_{i}\cdot[A]_{\mathit{prod}}\cdot\frac{\partial\mu}{\partial x} (10)

It is important to note here that μ/x\partial\mu/\partial x it is the change in dipole moment (μ\partial\mu) as a function of the displacement produced by the vibration or vibrational mode (x\partial x). This μ/x\partial\mu/\partial x is the IR cross section as determined from our quantum-chemistry calculations; i.e. without convolving with any peak profile function (Lorentzian, Gaussian, etc.). Finally, the abundance-weighted PAH-C60 mixture spectra (displayed in Figures 6b and c) were constructed by summing the individual weighted PAH-C60 adduct IR spectra; as obtained for the two initial abundances [A]0 of indene from Cernicharo et al. (2021) and Burkhardt et al. (2021). Evidently, due to [A]prod[A]_{prod} and Pi PAH-C60 adducts have a lower/higher contribution to the abundance-weighted spectra.

Appendix E Details on the 1TnC60 mono-adduct configuration

Structurally, it is possible to build two mono-adduct configurations for 1TnC60, bonding the PAH by a symmetric or asymmetric attachment (see Figure 9). The symmetric model represented in Figure 9 is much more thermodynamically unstable than its asymmetric counterpart. Since we have established the possible formation of 1TnC60 through an exohedral cycloaddition, the binding to C60 must occur inside an hexagon ring instead on the C-C bond bridging the fused rings, which corresponds to the asymmetric model. In addition, the asymmetric model is in agreement with previous experimental studies that suggest that tetracene binds C60 in such a way (Sarova & Berberan-Santos 2004).

Refer to caption
Figure 9: Left: Representation of the 1TnC60 mono-adduct models in their asymmetric (1TnCasym60{}_{60}^{\rm asym}) and symmetric (1TnCsym60{}_{60}^{\rm sym}) binding configurations. The corresponding Gibbs free energies of formation are indicated. Right: The theoretical IR spectra for each 1TnC60 mono-adduct configuration.

Appendix F High-resolution windows of PAH-C60 adduct spectra

Refer to caption
Figure 10: High-resolution windows of the more stable mono-adducts presented in Figure 4. Red dotted lines represent the free PAHs.
Refer to caption
Figure 11: High-resolution windows of the more stable bis-adducts presented in Figure 7.
Refer to caption
Figure 12: High-resolution windows of the more stable tris-adducts presented in Figure 8.
Refer to caption
Figure 13: High-resolution windows of the summed quantum-chemical spectra (top) and the abundance-weighted spectra (bottom) presented in Figure 6.