Open Access
Volume 686, June 2024
Article Number A181
Number of page(s) 9
Section Interstellar and circumstellar matter
Published online 11 June 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1 Introduction

Polycyclic aromatic hydrocarbons (PAHs) are a group of organic molecules that contain two or more aromatic carbon rings. These molecules have been observationally confirmed in the interstellar medium (ISM), galaxies, post asymptotic giant branch (post-AGB) stars and protoplanetary discs (e.g. Allamandola et al. 1985, 1989; Meeus et al. 2001; Kennicutt et al. 2003; Geers et al. 2006; Acke et al. 2010) through their distinct infrared emission features between 3 and 20 µm obtained with ground-based and space-based observatories such as VISIR on the Very Large Telescope (VLT), the Spitzer Space Telescope, and the Infrared Observatory (ISO). Future missions such as Ariel and Twinkle will be also capable of detecting PAH emission features in the infrared in protoplanetary discs and in principle also in exoplanetary atmospheres (Ercolano et al. 2022).

Surveys in the past have searched systematically for these molecules in planet-forming discs (i.a. Habart et al. 2004). Geers et al. (2006) searched amongst T Tauri stars for signs of PAHs, finding that only 8% of the T Tauri discs show emission features, excluding tentative detections. Recently, Valegård et al. (2021) found a 27% detection rate of PAHs in intermediate-mass T Tauri star discs (IMMT) through visual inspection of the spectrum, excluding a few tentative detections as well. Acke et al. (2010) searched for PAH features in Herbig Ae/Be discs and found a detection rate of PAH features of about 60%. Despite a clear relation between effective stellar temperature and the PAH detection rate that can most likely be traced back to the weaker UV fields in fainter and colder systems, it is not well understood what causes the difference in PAH emission between two sources of similar stellar type and why some sources are lacking emission at all.

Most of the available data have been obtained by the Spitɀer Space Telescope using the infrared spectrograph (IRS) between 5 µm and 40 µm. Unfortunately, data for the very distinct 3.3 µm PAH feature are only available for a few sources obtained by the ISO, SWS, ISOPHOT or VLT ISAAC (see e.g. Seok & Li 2017, for an overview of available data). One year ago, the James Webb Space Telescope (JWST) launched, carrying the near-infrared spectrograph (NIRSpec), the near infrared camera (NIRcam) and the mid-infrared instrument (MIRI). They allow obtaining high-quality data of the important PAH features at 3.3 µm, 6.2 µm, 7.7 µm, 8.6 µm, and 11.3 µm with a single facility. These instruments are well-suited for characterising the PAH content in protoplanetary discs, especially with the NIR-Spec and MIRI integral field units, which for the first time allow us to obtain spatially resolved spectra of PAH emission especially aimed at T Tauri discs. The JWST is expected to provide an unbiased and complete overview of the spatial and spectral characteristics of PAHs in protoplanetary discs. To guide JWST observations and interpret the amount of newly available data in the coming years, it is necessary to make models so that the upcoming results can be properly interpreted. For this purpose, we created two disc models that contain PAHs and performed radiative transfer calculations to investigate two different aspects of the upcoming observations. First, we wish to analyse the spectra of protoplanetary discs seen under different inclinations and characterise the different PAH feature strengths, their origin, and the possible interpretations of upcoming new PAH detections. Second, both NIRSpec and MIRI on board the JWST can provide integral field unit spectroscopy at unpreceeded resolution and sensitivity, allowing us to extract radial profiles of the important PAH features at 3.3 µm, 6.2 µm, and 11.3 µm. Therefore, we wish to provide calculated profiles and predictions for the variation in the relative feature strengths between the three PAH features as a function of the distance from the host star.

This work is structured in the following way: first, we explain our disc model and our setup for the radiative transfer code RADMC-3D. Then, we present the results we obtained for the SEDs under different inclinations, followed by a discussion of the radial profiles that we expect from a (dust)-gapless disc. Finally, we discuss and summarise our results.

2 Methods

In order to simulate the PAHs in a protoplanetary disc, we require a disc model, a radiative transfer model that determines the radiation field and temperature structure of the disc, and a model for the PAH emission itself. For the radiative transfer model, we used the commonly used code RADMC-3D1 (Dullemond et al. 2012), which can simulate the continuum radiative transfer of dust populations in protoplanetary disc environments, for example. It uses a Monte Carlo method to calculate the radiative equilibrium temperature in each grid cell. Photon packages that originate from heating sources such as a central star are distributed throughout the grid, where they can be absorbed and re-emitted at a different wavelength or scattered into a different direction. Then, if the temperature structure of the disc is known, RADMC-3D uses ray-tracing in order to obtain images and spectra from a certain inclination and distance of the observer to the object. However, RADMC-3D does not self-consistently compute the gas temperature in a disc as it does not feature detailed photochemistry. Hence, as PAHs and the stochastically heating and cooling of PAHs is not included in RADMC-3D itself, we are required to perform these calculations outside RADMC-3D and supply them later for the radiative transfer. Therefore, we separate this section into two parts: the disc model and setup that we used for RADMC-3D to obtain the radiation field in the UV, dust grain temperatures, and the obtained continuum emission. In the second part, we calculate the stochastically heating of the PAHs for different UV fields and the corresponding temperature distributions in order to calculate the PAH emission spectrum at the corresponding UV field in the disc.

2.1 Disc model

We focus on a smooth protoplanetary disc (no dust substructures as seen e.g. with ALMA in the DSHARP survey, Andrews et al. 2018) in order to link our studies to earlier studies on PAHs and because our aim is to define a general framework rather than modelling a specific disc and its properties (for the effect of dust gaps on PAH emission features see Maaskant et al. 2014). We chose to model a disc around a Herbig star, where we take HD 169142 as a reference model for the typical Herbig star. However, these models can be readily extended to T Tauri stars, even though their UV radiation field is much weaker than in a typical Herbig star disc and additional photo-chemical processes of importance such as consideration of highly energetic UV and X-ray photons (Siebenmorgen & Krügel 2010) might be necessary.

The gas in our disc follows the vertical equilibrium, where the density profile is given by (1)

where ρmid(r) marks the midplane density that is linked to the gas pressure scale height H and the surface density Σ(r) through . The gas-pressure scale height is defined as H = csK , with the speed of sound cs and the Keplerian frequency ΩK. The speed of sound is given by , with the Boltzmann constant kB and the mean molecular weight of the gas (µ = 1.37). Our disc follows a gas surface density profile scaled to a minimum mass solar nebula (MMSN), (2)

as described by Weidenschilling (1977), so that the total mass of the disc is Mdisc = 0.01 M*. Our temperature profile is defined through the commonly used Hayashi (1981) profile (3)

For the disc and stellar parameters, we chose the inner disc radius to be rin = 0.5 au, the outer disc radius as rout = 100 au, the stellar mass M* = 2.28 M, the stellar luminosity L* = 15.33 L, and the stellar temperature T* = 8200 K (Honda et al. 2012). For the dust grains, we chose a minimum grain size of amin = 0.05 µm and a maximum grain size of amax = 1 mm following a Mathis-Rumpl-Nordsieck (MRN) size distribution (Mathis et al. 1977), (4)

where γ = 3.5. Our disc is assumed to be a settled disc in which grains are in an equilibrium state between vertical settling and vertical turbulence. We followed the model of Fromang et al. (2007) to describe the vertical distribution of the grains, (5)

where we obtained the steady-state solution by integrating the differential equation to solve for the dust density ρd, (6)

We used the standard relations provided by Fromang et al. (2007) for the diffusion coefficient D and stopping time τs, which were also used and described in Lange et al. (2023).

2.2 Stochastic heating of PAHs

The PAHs are large aromatic molecules that are efficient in absorbing UV photons and emitting the gained energy as spectral features in the near-infrared (NIR) and mid-infrared (MIR). As PAHs are individual molecules, their heat capacity is small, and the absorption of a single UV photon can raise their electronic excitation temperature by several hundred Kelvin (Tielens 2008). As a consequence, the excitation temperature of PAHs is not well described by the radiative equilibrium temperature like for dust grains, but requires a stochastic treatment of the heating and cooling process. Instead of describing the temperature only through the equilibrium value, the temperature probability distribution G(T) needs to be determined to properly treat and model the emission of spectral features (Bakes et al. 2001). For a low photon absorption rate (≤1 s−1), the temperature probability distribution of the PAHs is dominated by the absorption of a single photon heating the PAH to a peak temperature, followed by a phase of cooling down to the background temperature. For a high UV photon flux without consideration of photochemical processes that require energy such as photodissociation, the PAH will not have enough time to cool down to the background equilibrium temperature and will fluctuate around a mean temperature value in a broad Gaussian-like distribution (for details and limits of this Gaussian distribution, see Tielens 2021; Lange et al. 2021). We evaluated the temperature distribution function for PAHs irradiated by a given UV field using a Monte Carlo method that accounts for the stochastic radiative cooling and heating processes.

We modelled the heat capacity of the PAHs by the micro-canonical temperature-energy relation2 (Tielens 2021), (7)

where N describes the number of atoms in a PAH molecule, Tm is the microcanonical temperature, and E is the excitation energy. From this relation, it becomes clear that the smallest PAHs undergo the strongest temperature fluctuations and reach the highest peak temperatures. The absorption rate of photons ra is determined by the stellar spectrum that we assumed for Herbig stars to be a black body at the given stellar effective temperature of the host star and the absorption cross section of the PAH. The absorbed energy flux Φa,E at energy E is then (8)

where FE is the flux density at given E and σE is the absorption cross section of a PAH at a given photon energy. We considered photons between 0.1 and 13.6 eV and bin the photon distribution into 100 bins, so that the absorption rate of photons of a given bin m with average energy Em is (9)

As the ionisation of hydrogen significantly increases the optical depth above 13.6 eV (Ryter 1996), we did not account for higher-energy photons. By assuming an average photon energy in the infrared EIR , we can determine the rate re at which photons are emitted for cooling, (10)

if the cooling rate is known. Together with Eq. (7), we applied the cooling model from Bakes et al. (2001) and Tielens (2021) to determine the cooling rate, (11)

Then, in order to evaluate the statistical heating and cooling, we determined the total rate of events for a Monte Carlo evaluation, (12)

Similar to Zsom & Dullemond (2008), we used an exponential distribution for the time between two events, (13)

with R = {ra,m, re}. Then, the probability of an event j with rate rj ∈ R can be described by (14)

and the event is randomly drawn with the given probabilities. As a result, we received the fluctuating temperature of the PAHs at any time and calculated the temperature probability distribution G(T).

2.3 PAH emission model

As the next step in our calculation, we need to calculate the PAH emission spectra from the temperature probability distribution G(E). We used the NASA AMES PAH IR spectral database, where theoretical PAH spectra computed with density functional theory (DFT) are available (Boersma et al. 2014; Bauschlicher et al. 2018; Mattioda et al. 2020). Because we focused our calculations on the general spectral features of PAHs rather than modelling a specific disc, we use dcoronene C24H12, circu-moronene C54H18, and circumcircumcoronene C96H24 in our calculations in the neutral and cation state to cover the whole size range of astrophysically relevant PAH molecules. The spectra of PAH clusters are available in the database, and we expect that PAHs will be present on dust grains as clusters with different sizes, ranging from a few cluster members up to a few hundred. However, we chose to use the spectrum of individual PAH molecules for both frozen-out (hereafter also called adsorbed) and gas-phase PAHs. When we compare the results for photodissociation and desorption of PAH clusters in protoplane-tary discs (Lange et al. 2021), the timescale for photodissociation of clusters is much shorter than the timescale for desorption of the cluster from the dust grain. Hence, we do not expect a large fraction of clusters to be present in the gas phase, but to be dom-inantly adsorbed. Additionally, the selection of clusters available in the database does not cover all expected cluster sizes.

Using the temperature probability distribution G(T) from our Monte Carlo model, we followed the approach by Bakes et al. (2001), where the full emission profile of a temperature fluctuating PAH in the gas-phase is given by (15)

where nPAH,g defines the number density of gas-phase PAHs, L(υ, Z) describes the emission profile of the PAH as a function of the charge state Z, f gives the fraction of PAHs that have a charge state Z, and B(T, v) is Planck’s function. The integral over the Planck function B and temperature probability distribution G(T) defines the time-averaged black-body emission function. We calculated the PAH emission L(υ, Z) by considering all available transitions as Lorentzians with the standard full width at half maximum (FWHM) value of 15 cm−1 (as in the NASA AMES spectral database) and integrated cross sections that are given in the database3.

For frozen-out PAHs, we considered the temperature fluctuations to be negligible as the heat transfer from the PAH to grain occurs on a much faster second timescale than the IR emission of an IR photon (Lange et al. 2023). Therefore, we assumed that the frozen-out PAHs have the same temperature as their hosting dust grains at radiative equilibrium. Under the assumption that all adsorbed PAHs are neutral, the PAH emission spectrum is calculated by (16)

where nPAH,s is the number density of frozen-out PAHs on dust grains. Because we did not consider the destruction of PAHs, we ensured that nPAH,s + nPAH,g always equaled the expected abundance of PAHs from the ISM.

Because we did not include the statistical heating and cooling of PAHs into RADMC-3D itself, we needed to parametrise our calculations to make them accessible to RADMC-3D. For this purpose, RADMC-3D provides the function userdef_srcalp() to add fluxes to every grid cell without interfering with the initial calculations of the dust temperature and calculation of the local radiation field. Unfortunately, this means that any contributions to the opacity for the whole disc model are not considered or must be taken into account separately in a non-consistent way. However, because we only considered a PAH abundance of 1% of the ISM value, the contribution of PAHs to the opacities is negligible.

We chose the local UV field characterised by G0 to be our parametrisation for all calculations. For this purpose, we sampled the temperature fluctuations G(T) of the gas-phase PAHs in a logarithmic grid between 1 ≤ G0 ≤ 109 with 40 grid points and interpolated then to the nearest neighbour. Opacity effects by the dust were then considered by RADMC-3D during the image and SED calculation itself.

We used two disc models that we wished to compare in our studies. Our smooth-disc model has a gas-phase PAH abundance of 1% of the ISM relative to the gas density (0.01 • CPAH: H = 1.5 × 10−5 Tielens 2008)), where the 99% of the ISM abundance are frozen out on grains. Our PAH-gap disc model has a gas-phase PAH abundance of only 0.1% relative to the ISM in the inner disc (r ≤ 40 au) and a gas-phase PAH abundance 1% relative to the ISM in the outer disc, where again the remaining PAHs are frozen out. We chose this model as a test model for the results of Lange et al. (2023), where coagulation and adsorption can lead to a depletion of PAHs in the inner disc in a Herbig disc model. In this way, we can test the effect on the observable PAH flux ratios and detectability with current facilities such as the JWST.

2.4 lonisation alance

Furthermore, we took the charge in the PAHs into account in our emission model. To do this, we evaluated the PAH charge and followed the ideas of Maaskant et al. (2014); Bakes & Tielens (1994); Tielens (2005), and references therein. In their approach, the ionisation balance is set by the ionisation parameter (17)

where ne is the free electron density. Assuming that the photosphere behaves like a photodissociation region (Woitke et al. 2009), we assumed that the free electrons are supplied by ionised carbon C+, and we assumed an ISM abundance of elemental carbon that is dominantly ionised in an UV environment, nC = ne = 1.5 × 10−4n0. Then, the fraction of neutral PAHs is determined by (18)

where the ratio of ionisation to recombination is given through (19)

where the pre-factor orginates from the electron-sticking probability (Tielens 2005), and NC is the number of carbon atoms of the PAH molecule. This is a simplified approach, and a detailed model would require a complex disc chemistry framework to solve the charge balance in the different molecular species.

2.5 Calculation sequence

The setup and calculation of the models involves calculations that are performed outside RADMC-3D, but require information about the disc setup and the radiative transfer model itself. To calculate the opacities of the dust grains, we used optool4 Dominik et al. (2021), a command line tool for calculating dust opacities for various mixtures of materials and grain size distributions. For our dust, we chose the standard DIANA opacities (Woitke et al. 2016). The calculations were performed in the following order. First we performed Monte Carlo calculations for the PAH temperature fluctuations in the range of the UV field 1 ≤ G0 ≤ 109 that we expect in the disc. Then, we used the PAH emission model to obtain the expected fluxes of each individual PAH molecule in the given UV radiation field. Next, we set up the disc model itself in RADMC-3D and obtain the temperature structure of the disc with mctherm and derived the local radiation field of each grid cell with mcmono in the photon-energy range between 5 eV and 13.6 eV, and we calculated G0 in each grid cell. Based on G0 and the number density of gas-phase PAHs in the grid cell, we calculated the PAH flux that is emitted in each grid cell. Finally, we supplied RADMC-3D with the PAH fluxes in the function userdef_srcalp() and calculated images and spectra. For each of the images, we used 5 × 108 photons, and for the spectra, we used 108 photons at each wavelength.

3 Results

3.1 PAH SEDs from protoplanetary discs

We wish to analyse the expected PAH SEDs from protoplanetary discs in a systematic way. Therefore, we ran radiative transfer models with RADMC-3D to simulate observations of discs at different inclinations and locally different PAH abundances. We aim to understand the general effects and do not wish to reproduce the observational spectra obtained by Acke et al. (2010); Seok & Li (2017) for example, and we used the theoretical emission profiles of coronene (C24H12) from the NASA AMES PAH IR spectral database rather than the commonly used model by Draine & Li (2007). The model by Draine & Li (2007) was built to reproduce the shape of PAH features observed in protoplan-etary discs by an average disc PAH spectrum. It can be scaled to different sizes, including a transition to small carbonaceous grains. However, the strength of the features is built to match general trends from observations and is only guided by calculated PAH spectra. The effect of the radiation environment in the disc (neutral to cation ratio, PAH sizes, and disc structures) on PAHs is therefore already implicitly built into the model by Draine & Li (2007). By using the NASA Ames database, we used only the intrinsic properties of PAHs and avoided these implicit assumptions about discs.

Furthermore, it is not possible to identify individual PAHs or to fit features of specific substituents (e.g. replace N or O with C) with it. Instead, advanced models including many PAH species should be developed that not only reproduce the observations, but also allow conclusions about the underlying chemistry and radiation field. We wish to motivate the usage of such a model by applying single PAH spectra from the database to a radiative transfer code.

In the first part, we consider the contributions of PAHs in the different states (gas-phase and frozen) to the SED, and then, we analyse the effect of disc inclination on the observable PAH feature-to-continuum ratio and the integrated PAH line fluxes.

Figure 1 displays an SED of a protoplanetary disc with PAH emission (dashed black line) observed face-on (0°) that is composed of neutral gas-phase PAH emission (green), cationic gas-phase PAH emission (red), dust continuum emission (blue), and emission from PAHs that are adsorbed on dust grains (orange). Similar to the results of Maaskant et al. (2014), the PAH emission of our disc model without any dust gaps is dominated by neutral PAHs. Because gas-phase PAHs undergo the strongest temperature fluctuations, they contribute most to the PAH typical infrared features between 3 and 13 µm to the spectrum. In contrast, the adsorbed PAHs are thermally coupled to the dust grains at the radiative equilibrium temperature. The non-fluctuating thermal emission of the adsorbed PAHs does not significantly contribute to features between 3 and 9 µm, even though there are 100 times more frozen-out PAHs than gas-phase PAHs in our model and dust grains in the inner disc can reach several hundred Kelvin. However, there is only contribution to the features at 11 µm, 13 µm, and 18 µm. For the smooth-disc model, the PAH clusters contribute 20% of the emission to the 11 µm feature, while the contribution in the PAH-gap model is equal between adsorbed clusters and gas-phase PAH. Consequently, the integrated5 feature ratios are affected by the disc model. For the smooth-disc model, I3/I11 = 0.3, I6/I11 = 0.08, while for the PAH-gap model, I3/I11 = 0.1 and I6/I11 = 0.05.

Therefore, the relative intensity between the long-wavelength features and the short-wavelength features does not only resemble the charge state of PAHs (e.g. 6.2 µm to 11.3 µm ratio Allamandola et al. 1999), but also the relative abundance of PAHs in the disc (hot inner disc versus cold outer disc) and the aggregate form of PAHs (strong temperature fluctuations as monomers in gas-phase versus thermally coupled PAHs and PAH clusters frozen-out on dust grains). We also find that the total integrated PAH emission over all features in the gap model is missing 80% of the smooth-disc model PAH luminosity, even though only 15% of the total gas-phase PAHs is additionally adsorbed in the PAH-gap model. Hence, the total PAH emission might not only reflect the total PAH abundance in the disc, but also strongly depends on the spatial distribution of the PAHs within the disc.

As the next step, we investigated the difference between the smooth-disc model and the PAH-gap model when the disc is observed under different inclinations. Figure 2 shows the SED of our two disc models for the MIR wavelength range that are accessible to the JWST for different selected disc inclinations. For a low and intermediate inclination, the silicate features between 9 and 12 µm contribute significantly to the spectrum, which typically makes the extraction of weak PAH features very difficult because forsterite and enstatite, for example, have peaky features in their opacities (van Boekel et al. 2005). Hence, the lowest contamination from the dust grains exists in edge-on discs, which should provide the best PAH spectra that can be obtained. We wished to analyse only the PAH features, and therefore, we subtracted the continuum from the SED and integrated the residual PAH features in the aforementioned intervals. We then plot the derived integrated PAH flux ratios for the easiest detectable PAH features at 3 µm, 6 µm, and 11 µm as a function of disc inclination in Fig. 3. For the models at low to medium inclinations (face-on to 70° inclination), the ratio of the integrated line fluxes remains constant, and the inclination does not play a role. We note again that the difference in I6/I11 of the two disc models is related to the different temperature structures of the excited PAHs and not to the level of ionisation because only a very small fraction of the PAH emission originates from ionised molecules. Only for very high inclinations (≥75°) do all ratios decrease because the short-wavelength PAH emission originating from the inner disc is blocked by the outer disc and can no longer be observed. As a consequence, the line ratios of the two disc models at high inclinations are barely distinguishable and only provide information on the PAHs of the outer disc. Even though edge-on discs provide the highest feature-to-continuum ratio and suppress the coinciding silicate feature emission, their SED provides only very limited information about the PAHs and their UV environment in most parts of the disc.

thumbnail Fig. 1

Simulated SED contributions of the smooth-disc model (left) and the PAH-gap model (right) seen in a face-on view. The PAH features are dominated by neutral emission, as expected (Maaskant et al. 2014). While the smooth-disc model emission mainly comes from gas-phase PAHs, the PAH-gap model has a significant contribution from the adsorbed PAHs that dominate the longer-wavelength features. The integrated PAH emission in both models does not scale with the total number of PAHs in the disc because five times more PAH flux is missing than in PAH molecules.

thumbnail Fig. 2

Simulated spectra of the two PAH disc models for selected inclinations. The left panel shows the SED of a smooth PAH disc, and the bottom panel shows the disc with a 40 au PAH gap containing only 10% PAHs in the gap compared to the smooth-disc model. Even though the strength of all PAH features decreases in the gap model, the short-wavelength features (3 µm and 6 µm) become much weaker than the 11 µm feature.

thumbnail Fig. 3

Intensity ratio of the 3.3 µm and 6.2 µm feature compared to the 11.3 µm feature. The PAH gap in the inner disc clearly affects the observable feature ratio by a factor 1.6 to 3.

3.2 Radial profiles of the PAH emission

We furthermore used RADMC-3D to produce artifical images to determine the radial profiles of the PAH emission at the most prominent emission features that will become accessible through the JWST. We always produced two artificial images: one including PAH emission, and the other containing only continuum emission produced by the dust grains so that we were able to calculated the contribution from the PAHs by subtracting the continuum for every pixel. We produced our images with a spatial resolution of 2.3 au px−1 at the feature wavelengths for coronene C24H12 at 3.259 µm, 6.277 µm, and 11.59 µm, for circumcoronene C54H18 at 3.26 µm, 6.26 µm, and 11.1 µm, and for circumcircumcoronene C96H24 at 3.26 µm, 6.32 µm and 11.0 µm. As the fraction of ionised emission is minor, we did not need to integrate over a wavelength interval around the feature to take the shift in the line wavelength into account. We then determined the distance of each pixel to the centre of the image and plot the flux in each pixel at the corresponding distance.

Figure 4 shows the radial profiles for all three PAH molecules (coloured lines) and continuum emission (grey). The radial profiles of the PAH emission between the inner and outer disc differ for all molecules, and we find that each region can be described by a simple power law with a shape of f (x) = Axγ. We attribute the difference between inner and outer disc to the temperature fluctuations of the corresponding PAHs. In the inner disc, the temperature of the PAHs fluctuates in a narrow distribution around the radiative equilibrium temperature in a multi-photon regime, and all PAH features between 3 and 11 µm are excited because the temperature is high enough to allow emission the entire time. The shallow power-law slope therefore reflects multiple properties. On the one hand, the radial density gradient in the number of emitting PAHs, as well as the flaring angle of the disc, because more flaring allows more absorption of UV photons at a given distance. On the other hand, the width of the temperature probability distribution, where a broader distribution leads to a slightly steeper radial gradient (the slope is more shallow for larger PAHs). In the outer disc, PAHs are mostly flash-heated by single photons and then cool down to the background temperature. During cooling, the excitation temperature is high enough most of the time to lead to emission of PAH features at 6 µm and 11 µm. However, the 3 µm feature requires at least two photons to be excited for an extended period of time, and the emission of the line depends on the high-energy tail of the temperature fluctuations. Therefore, the larger the PAH, the greater the heat capacity and the shorter the period of time for which the 3 µm feature is excited. As a consequence, the power law becomes steeper for larger PAHs. Additionally, we note that the transition between the two power laws increases with wavelength of the emission feature and size of the PAH, even though the transition at 6µm is very similar for all three PAHs.

Finally, Fig. 5 shows the relative peak intensities Ipeak of the three features normalised to the mean value that is obtained from the whole disc. As expected, the I3/I11 has the strongest relative variation through the disc, showing a variation of a factor 25 between the inner and outer disc. The I3/I6 intensities behave very similarly, but show a weaker variation of only a factor 10. Because both the I6 and I11 intensities are not as sensitive to the temperature fluctuations, the variation throughout the disc is only minor, of a factor of 2. We do not find a relevant change in the variation of intensities for different PAH species and sizes.

thumbnail Fig. 4

Extracted radial profiles of the most prominent PAH emission features from radiative transfer simulations. The inner and outer disc are fitted with power laws. The inner disc is dominated by multi-photon absorptions with a shallow slope, while the PAHs in the outer disc are heated only through single photons.

thumbnail Fig. 5

Ratio of the most prominent PAH line ratios normalised to the mean value obtained from the whole disc. Because the 3 µm line is very sensitive to the UV field of the disc, the line ratio changes very strongly throughout the disc. However, the size of the PAH molecule plays only a minor role.

4 Discussion

Our results show that even for our two very simple radiative transfer models, the information that can be obtained through investigation of the SED is very limited. In the literature, the observed spread in PAH luminosities in Herbig stars ranges from LPAH/L* = 8.3 × 10−5 (HD 58647) to LPAH/L* = 2.6× 10−2 (Oph IRS 48; Acke et al. 2010). This difference has been analysed and interpreted as a difference in PAH mass between 0.052 × 10−6 M and 24.3 × 10−6 M (Seok & Li 2017). This cor-respsonds to an almost linear relation between PAH luminosity and abundance. Because our simple PAH-gap model lacks only 15% of the PAHs compared to the smooth-disc model, but has 80% less PAH emission, it shows that the correlation between PAH abundance and PAH luminosity is not necessarily linear. This might imply that the observed spread in PAH abundance might be smaller than previously calculated and might instead indicate a different spatial distribution of PAHs in the various discs. Especially the various different substructures of simple carbon species seen with ALMA in the MAPS survey (e.g. see Öberg et al. 2021; Law et al. 2021) links the majority, but not all, substructures to the dust continuum and the chemical interaction and evolution between the disc midplane and surface. Hence, given the complexity and amount of substructures seen with ALMA in the dust continuum (e.g. with DSHARP, Andrews et al. 2018), at least some information about the distribution of PAHs and substructures in the disc must be known or assumed in order to derive meaningful PAH abundances because their chemical evolution is likely also linked to dust grains (Lange et al. 2023). Additionally, information about the approximate inclination (i ≥ 70° or not) of the disc is useful to interpret the extracted information from an SED correctly. The JWST facilities and IFU of NIRspec and MIRI will provide this information in the near future.

We furthermore wished to evaluate the impact of the temperature structure on the measured relative line fluxes between 3 µm, 6 µm and 11 µm. Focussing on the example of coronene C24H12, the intrinsic integrated cross section of the 11 µm feature is σ11 = 180 km mol−1 (neutral) and σ11 = 190 km mol−1 (cation), while at 6.2 µm, the integrated cross section is σ6 = 13 km mol−1 (neutral) and σ6 = 455 km mol−1 (cation). Hence, the change in intrinsic cross section through charging of the PAH is a factor of 35 for the σ6/11 ratio. In our radial profiles ,however, we only see a variation of I6/11 throughout the disc of a factor of 2, which arises from the different temperature fluctuations. The exact variation of the intrinsic cross section also depends on the considered molecule and can be much larger when circumcir-cumcoronene C96H24 is considered, for instance. Nevertheless, due to the insensitivity to the temperature structure of the disc, the I6/I11 line ratio is a very good tracer of the ionisation state at all distances in the disc, allowing for more detailed electron chemistry models.

Based on the intrinsic cross sections for coronene, we expect 0.07 ≥ I6/I11 ≥ 2.5. However, when we compare the observed ratios for I6/I11 with Spitɀer (Acke et al. 2010), which typicall lies at 1 ≥ I6/I11 ≥ 3, our simple disc models are not able to reproduce these values because they heavily underestimate the fraction of ionised emission through usage of a smooth disc without dust gaps. ALMA observations show that most of the discs have prominent substructures. This emphasises the need of a detailed model for each individual disc in order to explain the observed feature ratios (see also Maaskant et al. 2014).

When we compare for coronene C24H12 the 3 µm and 11 µm cross sections, the intrinsic integrated cross-section ratio is σ3/11 = 0.8 (neutral) and σ3/11 = 0.3 (cation). Even though the 3 µm feature is suppressed for cationic emission, the observed change in the intensity ratio in our models is a factor of 35. Hence, the I3/I11 is clearly a strong indicator for the UV field rather than the ionisation state because the 3 µm is very sensitive to the temperature fluctuations. Therefore, we propose to use the 3 µm feature to trace back the UV field. The radial profiles of all three PAH features might be used to characterise the UV field by determining the position at which the transition from the multi-photon to single-photon regime occurs. However, detailed modelling of the disc with its structure is likely necessary to obtain reliable results, and further studies of the sensitivity of the transition on the disc and stellar spectrum are required.

5 Conclusions

We have simulated the SED of two protoplanetary discs: one disc with a constant abundance of PAHs, and the other disc where the inner disc has 10% fewer PAHs than the outer disc. We furthermore provided radial profiles for the 3 µm, 6 µm, and 11 µm emission. The strength of the integrated PAH emission barely reflects the PAH abundance, but is rather sensitive to the spatial distribution of the PAHs and their absorbed UV radiation. Considering the complex substructures of protoplanetary discs, it is questionable how conclusive derived PAH abundances from SEDs are without a proper consideration of the gaps, rings, and flaring angle of the disc. The radial power-law slope of the PAH emission provides information about the size of the PAHs in a disc because especially the 3 µm feature is very sensitive to temperature fluctuations. The transition from a multi-photon regime to a single-photon absorption regime provides information about the underlying UV field at the transition location of the radial profiles. However, for quantitative statements about the PAH abundance, the PAH size, and the underlying UV field, a detailed model of the corresponding disc is required, which needs to be supplemented by additional observations such as observations with ALMA and SPHERE.


The authors thank the anonymous referee. The authors further thank Alessandra Candian for very fruitful and useful discussions. K.L. acknowledges funding from the Nederlandse Onderzoekschool Voor Astronomie (NOVA) project number R.2320.0130. C.D. acknowledges funding from the Netherlands Organisation for Scientific Research (NWO) TOP-1 grant as part of the research program “Herbig Ae/Be stars, Rosetta stones for understanding the formation of planetary systems”, project number 614.001.552. Studies of interstellar PAHs at Leiden Observatory are supported through a Spinoza award from the Dutch research council, NWO.


  1. Acke, B., Bouwman, J., Juhász, A., et al. 2010, ApJ, 718, 558 [Google Scholar]
  2. Allamandola, L. J., Tielens, A. G. G. M., & Barker, J. R. 1985, ApJ, 290, L25 [Google Scholar]
  3. Allamandola, L. J., Tielens, A. G. G. M., & Barker, J. R. 1989, ApJS, 71, 733 [Google Scholar]
  4. Allamandola, L. J., Hudgins, D. M., & Sandford, S. A. 1999, ApJ, 511, L115 [Google Scholar]
  5. Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018, ApJ, 869, L41 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bakes, E. L. O., & Tielens, A. G. G. M. 1994, ApJ, 427, 822 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bakes, E. L. O., Tielens, A. G. G. M., & Bauschlicher, Charles W., J. 2001, ApJ, 556, 501 [Google Scholar]
  8. Bauschlicher, C. W., Ricca, A., Boersma, C., & Allamandola, L. J. 2018, ApJS, 234, 32 [NASA ADS] [CrossRef] [Google Scholar]
  9. Boersma, C., Bauschlicher, C. W., Ricca, A., et al. 2014, ApJS, 211, 8 [NASA ADS] [CrossRef] [Google Scholar]
  10. Dominik, C., Min, M., & Tazaki, R. 2021, Astrophysics Source Code Library, [record ascl:2104.010] [Google Scholar]
  11. Draine, B. T., & Li, A. 2007, ApJ, 657, 810 [CrossRef] [Google Scholar]
  12. Dullemond, C. P., Juhasz, A., Pohl, A., et al. 2012, Astrophysics Source Code Library, [record ascl:1202.015] [Google Scholar]
  13. Ercolano, B., Rab, C., Molaverdikhani, K., et al. 2022, MNRAS, 512, 430 [Google Scholar]
  14. Fromang, S., Papaloizou, J., Lesur, G., & Heinemann, T. 2007, A&A, 476, 1123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Geers, V. C., Augereau, J. C., Pontoppidan, K. M., et al. 2006, A&A, 459, 545 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. Habart, E., Natta, A., & Krügel, E. 2004, A&A, 427, 179 [Google Scholar]
  17. Hayashi, C. 1981, Progr. Theor. Phys. Suppl., 70, 35 [Google Scholar]
  18. Honda, M., Maaskant, K., Okamoto, Y. K., et al. 2012, ApJ, 752, 143 [Google Scholar]
  19. Kennicutt, Robert C., J., Armus, L., Bendo, G., et al. 2003, PASP, 115, 928 [Google Scholar]
  20. Lange, K., Dominik, C., & Tielens, A. G. G. M. 2021, A&A, 653, A21 [Google Scholar]
  21. Lange, K., Dominik, C., & Tielens, A. G. G. M. 2023, A&A, 674, A200 [Google Scholar]
  22. Law, C. J., Loomis, R. A., Teague, R., et al. 2021, ApJS, 257, 3 [NASA ADS] [CrossRef] [Google Scholar]
  23. Maaskant, K. M., Min, M., Waters, L. B. F. M., & Tielens, A. G. G. M. 2014, A&A, 563, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425 [Google Scholar]
  25. Mattioda, A. L., Hudgins, D. M., Boersma, C., et al. 2020, ApJS, 251, 22 [Google Scholar]
  26. Meeus, G., Waters, L. B. F. M., Bouwman, J., et al. 2001, A&A, 365, 476 [CrossRef] [EDP Sciences] [Google Scholar]
  27. Öberg, K. I., Guzmán, V. V., Walsh, C., et al. 2021, ApJS, 257, 1 [CrossRef] [Google Scholar]
  28. Ryter, C. E. 1996, Ap&SS, 236, 285 [Google Scholar]
  29. Seok, J. Y., & Li, A. 2017, ApJ, 835, 291 [Google Scholar]
  30. Siebenmorgen, R., & Krügel, E. 2010, A&A, 511, A6 [Google Scholar]
  31. Tielens, A. G. G. M. 2005, The Physics and Chemistry of the Interstellar Medium (Cambridge, UK: Cambridge University Press) [Google Scholar]
  32. Tielens, A. G. G. M. 2008, ARA&A, 46, 289 [NASA ADS] [CrossRef] [Google Scholar]
  33. Tielens, A. G. G. M. 2021, Molecular Astrophysics (Cambridge University Press) [Google Scholar]
  34. Valegård, P. G., Waters, L. B. F. M., & Dominik, C. 2021, A&A, 652, A133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. van Boekel, R., Min, M., Waters, L. B. F. M., et al. 2005, A&A, 437, 189 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Weidenschilling, S. J. 1977, Ap&SS, 51, 153 [Google Scholar]
  37. Woitke, P., Kamp, I., & Thi, W. F. 2009, A&A, 501, 383 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Woitke, P., Min, M., Pinte, C., et al. 2016, A&A, 586, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. Zsom, A., & Dullemond, C. P. 2008, A&A, 489, 931 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]


As in Lange et al. (2021), we smoothed the transition between the low- and high-temperature regime to preserve differentiability, .


Note the unit conversion from km mol−1 to cm2 Hz through the conversion factor 105c/NA using the speed of light and Avogadro’s constant.


To account for the shift of lines between neutral/cation and broadness of the lines, our integration intervals are 3 µm: 3.1–3.4 µm, 6 µm: 6–7 µm, and 11 µm: 11–12 µm.

All Figures

thumbnail Fig. 1

Simulated SED contributions of the smooth-disc model (left) and the PAH-gap model (right) seen in a face-on view. The PAH features are dominated by neutral emission, as expected (Maaskant et al. 2014). While the smooth-disc model emission mainly comes from gas-phase PAHs, the PAH-gap model has a significant contribution from the adsorbed PAHs that dominate the longer-wavelength features. The integrated PAH emission in both models does not scale with the total number of PAHs in the disc because five times more PAH flux is missing than in PAH molecules.

In the text
thumbnail Fig. 2

Simulated spectra of the two PAH disc models for selected inclinations. The left panel shows the SED of a smooth PAH disc, and the bottom panel shows the disc with a 40 au PAH gap containing only 10% PAHs in the gap compared to the smooth-disc model. Even though the strength of all PAH features decreases in the gap model, the short-wavelength features (3 µm and 6 µm) become much weaker than the 11 µm feature.

In the text
thumbnail Fig. 3

Intensity ratio of the 3.3 µm and 6.2 µm feature compared to the 11.3 µm feature. The PAH gap in the inner disc clearly affects the observable feature ratio by a factor 1.6 to 3.

In the text
thumbnail Fig. 4

Extracted radial profiles of the most prominent PAH emission features from radiative transfer simulations. The inner and outer disc are fitted with power laws. The inner disc is dominated by multi-photon absorptions with a shallow slope, while the PAHs in the outer disc are heated only through single photons.

In the text
thumbnail Fig. 5

Ratio of the most prominent PAH line ratios normalised to the mean value obtained from the whole disc. Because the 3 µm line is very sensitive to the UV field of the disc, the line ratio changes very strongly throughout the disc. However, the size of the PAH molecule plays only a minor role.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.