© ESO, 2014

1. Introduction

The infrared (IR) spectra originating from protoplanetary disks are rich and diverse. Polycyclic aromatic hydrocarbons (PAHs) are observed in disks around Herbig Ae/Be stars (Meeus et al. 2001; Acke & van den Ancker 2004; Keller et al. 2008; Acke et al. 2010) and T Tauri stars (Geers et al. 2007b; Bouwman et al. 2008). The strength of the PAH features decreases with stellar effective temperature. No T Tauri stars of the spectral type later than G8 shows the IR emission features in their spectra (e.g. Tielens 2008; Kamp 2011, and references therein). In addition, the strength of the IR emission features is generally much weaker relative to the total IR emission in the Herbig star spectra as compared to the diffuse interstellar medium (ISM; Sloan et al. 2005; Keller et al. 2008; Acke et al. 2010; Boersma et al. 2008). Direct imaging observations taken in the PAH bands have shown that PAHs can be used as a tracer of the outer disk surface (Doucet et al. 2006; Lagage et al. 2006). Even so, the PAH emission is not dominated by PAHs on the surface of the outer disk in some transitional disks (Geers et al. 2007a; Maaskant et al. 2013).

After being electronically excited by ultraviolet (UV) photons, the PAHs cool by emission in the CH- and CC-stretching and bending modes. Laboratory studies and quantum chemical calculations show that many of the IR features shift in peak position, vary in width, and/or show substructure depending on ionization stage and detailed molecular structure (Hudgins & Allamandola 1999; Allamandola et al. 1999; Bauschlicher et al. 2009; Ricca et al. 2012). An important parameter that influences the relative feature strength of the CH- and CC-modes is the effect of ionization with CC modes being carried predominantly by ions and CH modes by neutrals. Hence, the 6.2-μm (CC-mode) to 11.2-μm (CH-mode) ratio provides a good measure of the ionization balance of the emitting PAHs (Allamandola et al. 1999). This has been demonstrated in a variety of astrophysical environments (Hony et al. 2001; Peeters et al. 2002; Rapacioli et al. 2005; Galliano et al. 2008).

Fig. 1 Left: Spitzer/IRS spectra of the four transitional disks HD 97048, HD 169142, HD 135344 B and Oph IRS 48. The spectra are scaled relative to the 6.2-μm PAH peak flux and show a wide range in the I6.2/I11.3 ratio; most neutral in HD 97048 and most ionized in Oph IRS 48. Middle: VLT/VISIR observations show that FWHM as a function of wavelength. The size of the 11.3-μm feature is respectively larger, similar and smaller than the N-band continuum for HD 97048, HD 169142, HD 135344 B, and Oph IRS 48 (data from Maaskant et al. 2013). The diamond symbols show the sizes of the point spread functions of the observations. The dash-dotted black line is the diffraction limit of the telescope. Right: Sketches showing a qualitative interpretation of the spatial distribution of the dominant contribution of PAH emission, as shown by the colored ovals. Each object is given a similar color in the left, middle, and right plots.

The intense radiation field of the central star may photo-dissociate the PAH molecules. As a result, the average PAH size increases, and the total abundance decreases (Visser et al. 2007; Siebenmorgen & Heymann 2012). As a possible path for PAH-survival, Siebenmorgen & Krügel (2010) suggest that turbulent motions in the disk can replenish or remove PAHs from the reach of hard photons. Several modeling studies have been carried out to understand global trends in the composition and evolution of PAHs in relation to the disk (Habart et al. 2004; Dullemond et al. 2007; Berné et al. 2009). However, none of them are focused on understanding PAH ionization.

In this paper, we investigate the effect of disk environments on PAH charge, spectra, and luminosity. The key question is: can we use the ionization balance of PAHs as a tracer of processes in protoplanetary disks? In Sect. 2, we first identify new observational trends seen in the I_6.2/I_11.3 PAH feature strength ratios. Section 3 describes the PAH model in the Monte Carlo radiative transfer tool MCMax. In Sect. 4, we show how the ionization of PAHs behaves as a function of star and disk properties. We use these results to demonstrate how we can understand the PAH spectra of four transitional disks in Sect. 5. In the discussion (Sect. 6), we discuss what we can learn about disk evolution from PAH ionization. The conclusions are given in Sect. 7.

2. New observational trends in PAH features of Herbig stars

In this section, we discuss new observational trends of PAH ionization in relation to the disk structure of Herbig stars. For a detailed observational overview of the four transitional disks in our sample, we refer to Maaskant et al. (2013). The PAH properties of all other Herbig stars are shown in Table 1 and are adopted from Acke et al. (2010).

2.1. Higher PAH ionization in dust-depleted “gaps”

For HD 97048, HD 169142, HD 135344 B, and Oph IRS 48, the Spitzer/IRS spectra are shown on the left of Fig. 1. The intensity of the spectra is scaled relative to the peak flux of the 6.2-μm feature. While the general PAH profile is highly similar, a large spread is observed in the relative strength of the 11.3-μm to the central 6.2-μm feature. The I_6.2/I_11.3 ratio is lowest for HD 97048 and increases for HD 169142, HD 135344 B, and Oph IRS 48.

In contrast to HD 97048 (Lagage et al. 2006; Doucet et al. 2006), it was found using VLT/VISIR mid-infrared (MIR) imaging that PAH emission in Oph IRS 48 is not co-spatial to the dust continuum (Geers et al. 2007a). While the dust continuum at 18.8 μm shows a ring structure, the PAH emission at 8.6 μm peaks at the center. In a recent paper by Maaskant et al. (2013), observations show that PAH emission is also not dominated by the outer flaring disk for HD 169142 and HD 135344 B. These observations are shown in the middle of Fig. 1, where the spatial full width at half maximum (FWHM) size of the 11.3-μm feature can be seen as imaged by VISIR (Geers et al. 2007a; Maaskant et al. 2013). For HD 97048, it is clear that the 11.3-μm feature originates in larger spatial scales than the continuum. The feature is not spatially resolved relative to the continuum in HD 169142 and HD 135344 B and is smaller than the continuum in Oph IRS 48. Using radiative transfer modeling, Maaskant et al. (2013) concluded that a significant contribution must come from close to the star and possibly from within the dust depleted “gap” if the PAH feature is not larger than the continuum.

On the left of Fig. 1 sketches are shown of the disk geometries which include colored ovals to indicate the origins of the dominant PAH contributions. The sizes of the disks are taken from the models of Maaskant et al. (2013), where the inner edge of the outer disk is set by MIR imaging. The inner disks represent the origins of near-infrared (NIR) emission in the spectral energy distributions (SEDs) though their structures are not well constrained.

The four transitional Herbig stars show a trend between the relative band strength ratio I_6.2/I_11.3 and the spatial origin in the disk. For HD 97048, where the PAH emission is dominated by the outer disk, the I_6.2/I_11.3 ratio is low: 0.87 ± 0.01. In HD 169142 and HD 135344 B, the 11.3-μm PAH features are not extended when compared to the continuum and have I_6.2/I_11.3 ratios of respectively 1.24    ±    0.01 and 1.42 ± 0.08. In Oph IRS 48, the PAHS are emitted from a location that is closer to the star than the outer disk, and have a high I_6.2/I_11.3 ratio of 3.76 ± 0.01. Thus it seems the PAH emission which originate from the outer dusty disk has a lower I_6.2/I_11.3 ratio, while PAH emission from within the “gap” has a high I_6.2/I_11.3 ratio. Since the I_6.2/I_11.3 ratio can be attributed to the charge state of the PAHs (Allamandola et al. 1999), we find that the PAHs in the optically thick dust disk are predominantly neutral, while the PAHs in the disk gap are ionized.

2.2. PAH ionization in Herbig stars

To better understand the behavior of the I_6.2/I_11.3 ratio for a broad range of disk geometries, we compare a large representative sample of Spitzer/IRS observations of Herbig objects (see Table 1). For 47 objects, Spitzer/IRS spectra are adopted from Acke et al. (2010). For 30 stars out of this sample, PAH luminosities were derived from detected features in the Spitzer/IRS spectrum. We discarded objects for which evidence of significant contamination from the surroundings exists in the literature. The objects are sorted by the F₃₀/F_13.5 flux ratio between 30.0 and 13.5 μm. This flux ratio is indicative for the disk geometry as it may trace gap-sizes in the temperature range of ~200–500 K and thus allows us to search for trends as a result of the geometry of the disk. The PAH luminosities and band-strength ratio 6.2/11.3 are included in Table 1. There is no evidence that the objects in Table 1, for which no PAH luminosity could be derived, represent a particular range in F₃₀/F_13.5 or L₁₃₀₀/L_∗.

For all the objects with PAH emission, Fig. 2 shows the relative PAH luminosities (L_PAH/L_∗) compared to the MIR spectral index F₃₀/F_13.5. The three different symbols indicate the SED group classifications from Meeus et al. (2001) and whether the object shows the amorphous 10- and 20-μm silicate features in the spectrum. The blue squares (group IIa) are interpreted as the “flat”, self-shadowed disks with silicate features. The red circles (group Ia) are “flaring/transitional” objects with silicate features. The orange objects (group Ib) are “flaring/transitional” objects without silicate features. The latter typically show prominent PAH spectra due to large gaps in the disk, which lower the continuum emission, thereby, increasing the contrast to the PAH features (Maaskant et al. 2013). Within the groups, L_PAH/L_∗ spans ~2−3 orders of magnitude. The averaged difference between the groups is much smaller; flaring/transitional sources are, relatively speaking, only a factor of ~2−3 more luminous in the fractional PAH luminosity than flat objects.

Fig. 2 PAH luminosity (LPAH/L∗) of a sample of Herbig stars compared to the MIR spectral index (F30/F13.5). The MIR spectral index can be used as a tracer of disk gaps, as indicated by the symbols showing the flaring/transitional disks (group I) with silicate features (red circles), without silicate features (orange triangles), and self-shadowed disks (group II) with silicate features (blue squares). While the PAH luminosities among the objects span ~2–3 orders of magnitude, the averaged luminosities of the disk groups only differ by a factor of ~2–3. The parameters of all the objects are listed in Table 1.

Fig. 3 I6.2/I11.3 ratio compared to MIR spectral index (F30/F13.5). Ionized PAHs are traced by a high I6.2/I11.3 ratio. The wider range in the I6.2/I11.3 ratio in flaring/transitional disks (⟨ I6.2/I11.3 ⟩  = 2.61 ± 1.63) suggests that the PAH emission originates in a wider range of physical environments while the origin of PAH emission in flat disks (⟨ I6.2/I11.3 ⟩  = 2.08 ± 0.56) seems to be more homogeneous. The parameters of all the objects are listed in Table 1.

Figure 3 shows the I_6.2/I_11.3 ratio compared to the F₃₀/F_13.5 MIR spectral index. This figure shows that the averaged ⟨ I_6.2/I_11.3 ⟩ ratio in flaring/transitional objects (2.61 ± 1.63) is higher and has a larger spread than in flat objects (2.08 ± 0.56). This may reflect that the PAH emission in flaring/transitional disks originates in more varied physical environments (with some cases, a strong contribution from ionized PAHs) while flat disks may be more homogeneous. In this paper we show that ionized PAH emission can arise from the gas inside the gap for transitional disks, which may arise in the gas flowing from the outer disk through the gap to the inner disk. For flat disks, there are few observational constraints on their disk geometries. However, for one particular case, there is evidence of extended PAH emission. Verhoeff et al. (2010) has shown for HD 95881 that PAHs have a spatial extent of ~100 AU and are well resolved with respect to the continuum. They, thus, confirmed the presence of an illuminated gas surface. As discussed by Dullemond et al. (2007), the opacity in the upper layers of the disk decreases, causing a weaker far-infrared (FIR) flux when dust decouples from the gas and settles to the mid-plane. As a consequence, PAH emission from flat disks is predicted to be stronger because the gas is still flaring. However, as Fig. 2 shows, this is not what is observed, suggesting that some mechanism is not only removing small grains but also the PAHs from the disk upper layers in flat disks.

2.3. PAH ionization and the millimeter luminosity

For most of the sources, the 1.3 millimeter (mm) photometry is collected from the literature. An overview of all available photometry is listed in Table 1. Since mm emission is optically thin and most prominently produced by cold mm-sized grains in the outer disk, its luminosity compared to the stellar luminosity can be used as a proxy of the total disk mass. Large uncertainties remain in the conversion of the sub-mm flux toward an estimate of the total disk mass, because of the unknown outer disk temperatures and the opacities of mm grains. Nevertheless, during the evolution of protoplanetary disks toward planetary systems, the disk loses most of the mass and the mm emission decreases. This makes it interesting to compare the mm fluxes of disks to other disk characteristics.

We report a new trend seen in the PAH I_6.2/I_11.3 flux ratio as a function of sub-mm brightness. Figure 4 shows the sub-mm luminosity at 1.3 mm (L₁₃₀₀/L_∗) compared to the I_6.2/I_11.3 PAH ratio for all Herbig stars with detected mm photometry. The I_6.2/I_11.3 flux ratio is higher for the sources with a lower mm flux. The gray line gives the best fit to the data. Using least squares minimization we find the following exponential best fit to the data for the disks of Herbig stars with mm luminosities between 10^-7 ≲ L₁₃₀₀/L_∗ ≲ 10^-4$\frac{{\mathit{L}}_{\mathrm{1300}}}{{\mathit{L}}_{\mathrm{\ast }}}\mathrm{=}\mathrm{1.1}\mathrm{\times }{\mathrm{10}}^{-4}\exp \left(\mathrm{-}\frac{{\mathit{I}}_{\mathrm{6.2}}}{{\mathit{I}}_{\mathrm{11.3}}}\right)\mathrm{·}$(1)The PAHs in massive disks are predominantly neutral (low I_6.2/I_11.3 flux ratio), while PAHs in lower mass disks are ionized (high I_6.2/I_11.3 flux ratio).

Fig. 4 An observational trend between the I6.2/I11.3 ratio as a function of the luminosity at 1.3 mm (L1300/L∗). The red dots show the objects studied in this paper. The gray line gives the best fit to the data (Eq. (1)). The parameters of all the objects are listed in Table 1. For high mm luminosities, the I6.2/I11.3 ratio is lowest (tracing neutral PAHs), while the I6.2/I11.3 ratio is highest (tracing ionized PAHs) for low mm values.

3. PAHs in the radiative transfer code MCMax

In this section we present our PAH ionization model and discuss how this is implemented in the radiative transfer code MCMax.

3.1. General description of PAH model

We model the disks using the radiative transfer code MCMax (see Min et al. 2009). MCMax is a code that not only solves for the radiative transfer but also solves the disk structure. For the continuum radiative transfer, it uses the Monte Carlo radiative transfer scheme of Bjorkman & Wood (2001) with direct re-emission and temperature correction. This scheme cannot be directly applied to the PAH molecules, since they are not in thermal equilibrium. Therefore, we need to compute the emission spectra of the thermal equilibrium grains first and iterate the radiative transfer to include the radiative transfer effect on the emissivities of the PAHs.

The emission spectra of large molecules, like PAHs or very small dust grains, can be computed using the continuous thermal approximation (Guhathakurta & Draine 1989; Draine & Li 2001). In this approximation it is assumed that the emission is well described by a single temperature at any time. This temperature fluctuates strongly due to transient heating events and fast cooling. A PAH molecule cools down significantly before it is hit by the next energetic photon, which heats it again. These fluctuation, and the corresponding temperature distribution can be computed from the statistics of the radiation field and the heating and cooling properties of the particles. This method fully accounts for multi-photon heating events, or, the case where the particle is not yet cooled down completely before it is hit again. This is important for larger particles with higher heat capacities (i.e., longer cooling times) and for the intense radiation field close to the star. The temperature distribution of the PAHs at each location in the disk is solved using a Monte Carlo method, simulating the transient heating events. An analytical method, using the solution of a set of linear equations to solve for the temperature distribution, is also implemented in MCMax (also used in Dullemond et al. 2007). The two methods perfectly agree on the emerging PAH spectrum.

Once the spectral shape of the emissivity of the PAHs is computed using the temperature distribution, the PAHs can be included in the Monte Carlo radiative transfer in a very similar way to the other dust components. For ach time a photon package in the Monte Carlo radiative transfer run is absorbed by a PAH molecule, we make sure that it is re-emitted with the same energy and with the spectral signature computed using the method above. In this way, we still strictly conserve energy, ensuring a quick convergence of the method.

The optical properties of the PAHs are computed using the method by Draine & Li (2001, 2007). The adopted opacities for neutral and ionized PAHs are shown in Fig. 5 (with N_C = 100, the number of carbon atoms in the PAH molecule, and PAH radius a_PAH = 6 Å). The PAH optical properties are constructed from a number of laboratory and observational studies. Thus our PAH model represents the emission from various PAH molecules and is designed to conform to observed PAH properties in the ISM. We use a single sized (N_C = 100) PAH molecule with two possible charge states. This PAH model is sufficient for the goal of our paper, which is to provide a tool to study the effect of the star and disk properties on the charge state of PAH molecules and its influence on the resulting IR spectra. PAH processing may also play a role in protoplanetary disks, either at the molecular level or in a disk environment. Specifically, the PAH destruction is predicted to act on the smallest PAH molecules (Visser et al. 2007), but the PAHs can be replenished from deeper layers of the disk by turbulence (Siebenmorgen & Heymann 2012).

Fig. 5 Opacities for neutral PAHs (solid blue line) and ionized PAHs (dashed red line) with NC = 100. Computed using the method described by Draine & Li (2007).

The I_6.2/I_11.3μm band ratio may also be affected by dehydrogenation of the PAH molecule. The effect of ionization versus dehydrogenation on the I_6.2/I_11.3 ratio are reviewed in Tielens (2008). The hydrogen coverage on a PAH is a balance between dehydrogenation of highly excited PAHs following UV absorption and hydrogenation reactions with atomic hydrogen (Tielens et al. 1987; Le Page et al. 2003). Theoretically and experimentally, it has been well established that dehydrogenation is an “on/off” switch; for example, when dehydrogenation becomes important, a PAH quickly loses all its hydrogens; over a small range of G₀/n_H (Ekern et al. 1997; Tielens 2005; Berné & Tielens 2012; Zhen et al. 2014). Hence, for any PAH, there is a critical G₀/n_H ratio and a PAH becomes fully hydrogenated (completely dehydrogenated) if the G₀/n_H ratio is less (larger) than this ratio. As the unimolecular dissociation rate depends exponentially on the internal excitation temperature and, therefore, in a given radiation field, on the size of the PAH, PAHs with a size smaller than a critical size are fully stripped of all their H (and do not contribute emission in any C-H mode), while PAHs larger than this size are fully hydrogenated. Small fully dehydrogenated PAHs are quickly destroyed through C₂ loss (or isomerized to cages and fullerenes) by further UV photon absorption (Joblin 2003; Zhen et al. 2014) and this is expected to be the dominant photo destruction route for small PAHs (Berné & Tielens 2012; Montillaud et al. 2013). Hence, as small dehydrogenated PAHs are quickly destroyed, the I_6.2/I_11.3μm band ratio would show little variation. Moreover, if dehydrogenation were important, one would expect to see the solo/duo/trio ratio change in the direction of a larger solo/trio ratio (Schutte et al. 1993), which is not what is observed (Acke et al. 2010; Hony et al. 2001). For these reasons, the I_6.2/I_11.3 ratio have been widely accepted as variations in the ionization state of the emitting PAHs, and we follow this analysis here.

The PAH-to-dust mass fraction is difficult to constrain as PAHs may take part in various processes. The PAHs may be prone to photo-destruction (e.g. Guhathakurta & Draine 1989; Le Page et al. 2003; Visser et al. 2007). Moreover, it has been argued that PAHs may take part in the dust coagulation process (Dullemond et al. 2007) or replenished by mixing processes in the disk (Siebenmorgen & Krügel 2010). Observational evidence of the role of PAHs in these processes are not conclusive and it is presently unclear how important these routes are. As for PAH destruction, very small, catacondensed PAHs can also lose C₂H₂ through UV photon absorption while the H-loss channel dominates. However, the dominant photochemical destruction channel for the carbon skeleton for large pericondensed PAHs is through complete dehydrogenation followed by C₂ loss (Ekern et al. 1997; Joblin 2003; Zhen et al. 2014). Larger PAHs are notoriously difficult to destroy (Jochims et al. 1994). Here, we treat the PAH/dust abundance ratio as a free parameter to be determined through our model fits of the observed SEDs.

3.2. The ionization balance

The charge of PAHs is determined by the balance between photo-ionization and electron recombination. The description of the ionization balance of PAHs is largely adopted from Bakes & Tielens (1994) and Tielens (2005). Here, we provide a brief summary. The ionization balance of PAHs is controlled by the parameter γ: $\mathit{\gamma }\mathrm{=}{\mathit{G}}_{\mathrm{0}}{\mathit{T}}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathit{/}\mathit{n}e$(2)with G₀ as the intensity of the radiation field between 6.0 and 13.5 eV in units of the Habing field ${\mathit{u}}_{\mathrm{UV}}^{\mathrm{Hab}}$, T the gas temperature, and n_e the electron density.

Fig. 6 Fraction of neutral PAH molecules f(0) given by the UV field, the temperature, and the electron number density and controlled by the ionization parameter G0T1/2/ne. In our models, we assume a PAH size of NC = 100.

The intensity of the UV radiation field inside the disk can be characterized by ${\mathit{G}}_{\mathrm{0}}\mathrm{=}{\mathit{u}}_{\mathrm{UV}}\mathit{/}{\mathit{u}}_{\mathrm{UV}}^{\mathrm{Hab}}$, where ${\mathit{u}}_{\mathrm{UV}}\mathrm{=}{\mathrm{\int }}_{\mathrm{6}\mathrm{=}\mathrm{eV}}^{\mathrm{13.5}\mathrm{eV}}{\mathit{u}}_{\mathit{\nu }}\mathrm{d}\mathit{\nu }\mathrm{=}{\mathrm{\int }}_{\mathrm{6}\mathrm{eV}}^{\mathrm{13.5}\mathrm{eV}}\mathrm{(}\mathit{h}\mathit{/}\mathit{c}\mathrm{)}\mathit{h\nu }{\mathit{N}}_{\mathrm{ph}}\mathrm{d}\mathit{\nu }\mathit{.}$(3)Here ${\mathit{u}}_{\mathrm{UV}}^{\mathrm{Hab}}\mathrm{=}\mathrm{5.33}\mathrm{\times }{\mathrm{10}}^{-14}$ erg cm^-3 is the energy density in the UV mean interstellar radiation field (Habing 1968), and N_ph is the number of photons in units of cm^-2 s^-1 erg^-1. The contribution from the interstellar radiation field is also taken into account. Regions in the disk, where the optical depth is high may not receive enough photons between 6 and 13.6 eV to calculate G₀. In those cases, G₀ is set to a lower limit value of 10^-6. Because of their high temperatures, Herbig stars generate strong photospheric UV fluxes. Additional UV excess due to accretion is relatively unimportant for these stars because of the low accretion rates in most targets. The effective temperature can, therefore, be considered as a direct indicator of the stellar UV radiation field.

The surface layer of the disk is essentially a photodissociated region (Woitke et al. 2009; Kamp et al. 2010) and the electron abundance, x_e = n_e/n_H, can be set equal to the C⁺ abundance, based on an equilibrium between the photoionization of neutral C and the recombination of C⁺ (Tielens & Hollenbach 1985). We adopt an ISM abundance of n_C = 1.5    ×    10^-4n₀ (Cardelli et al. 1996) with n₀ the gas density in cm^-3. The C / C⁺ ratio changes with depth but carbon is largely ionized whenever there are UV photons around (Tielens 2005). Therefore, we assume that ionized carbon with a fixed abundance is the dominant supplier of free electrons in the region where PAH emission originates,.

At higher altitudes in the disk surface, the gas temperature may exceed the dust temperature (e.g., Kamp & Dullemond 2004; Nomura & Millar 2005). However, it is beyond the scope of this paper to fully solve the thermo-chemistry balance. As the temperature only enters through a square root, we adopt T_gas = T_dust in our model studies.

We consider PAH with only two accessible ionization stages, the neutral and the singly charged cation. The neutral fraction, f(0) is given by $\mathit{f}\mathrm{\left(}\mathrm{0}\mathrm{\right)}\mathrm{=}\mathrm{[}\mathrm{1}\mathrm{+}{\mathit{\gamma }}_{\mathrm{0}}{\mathrm{]}}^{-1}$(4)with γ₀ the ratio of the ionization rate over the recombination rate, which, adopting the classical expression, is ${\mathit{\gamma }}_{\mathrm{0}}\mathrm{=}\mathrm{3.5}\mathrm{\times }{\mathrm{10}}^{-6}{\mathit{N}}_{\mathrm{c}}^{\mathrm{1}\mathit{/}\mathrm{2}}\frac{{\mathit{G}}_{\mathrm{0}}{\mathit{T}}^{\mathrm{1}\mathit{/}\mathrm{2}}}{{\mathit{n}}_{\mathrm{e}}}\mathrm{·}$(5)Here, the constant 3.5 × 10^-6 is derived from the probability that a colliding electron sticks (Tielens 2005). The MCMax code calculates the neutral fraction f(0) in each grid cell in the disk. For all grid cells with neutral fractions between f(0) = 0 (ionized) and f(0) = 1 (neutral), MCMax linearly interpolates between the neutral and ionized PAH emission profiles.

Figure 6 shows the relation between the neutral fraction and G₀T^1/2/n_e for different PAH sizes. We find that most PAHs are emitting in the physical range G₀T^1/2/n_e ~ 10³ − 10⁶ in our disk models. Therefore, our two charge model covers the physical range of our interest.

The ultimate goal of our modeling is to understand the range of band strength ratios I_6.2/I_11.3 for the four transitional disks (HD 97048, HD 169142, HD 135344 B, and Oph IRS 48) and to learn what differences in disk properties are responsible for neutral versus ionized PAH emission. In the next section we first describe a benchmark model and explore how PAHs behave as a function of the relevant parameter space. After that, we present PAH models of four transitional disks, which fit the I_6.2/I_11.3 ratio by combining contributions from neutral PAHs in the dusty outer disk and ionized PAHs in an optically thin gas flow through the gap.

4. Benchmark model

We study the behavior of the degree of ionization, the PAH luminosity, and the radial dependence of the PAH charge state throughout the disk. We focus on PAH emission from self-consistent disks models, which include high optical depth media. In this approach, the vertical density structure of the gas is set by hydrostatic equilibrium, and dust properties, such as scattering, grain sizes, and grain composition, are important input parameters for the radiative transfer. This allows us to evaluate how the charge state of PAH molecules in the disk translates to the spectrum. We start with a simple hydrostatic flaring disk model for which we adopt the stellar properties of HD 97048. All relevant star and disk parameters, as well as the results of the PAH characteristics, are given in Table 2.

We assume a dust composition of 80% silicate and 20% amorphous carbon. The dust composition with reference to the optical constants is 12% MgFeSiO₄ (Dorschner et al. 1995), 34% Mg₂SiO₄ (Henning & Stognienko 1996), 32% MgSiO₃ (Dorschner et al. 1995), 2% NaAlSi₂O₆ (Mutschke et al. 1998) and 20% C (Preibisch et al. 1993). The shape of our particles is irregular and approximated using a distribution of hollow spheres (DHS, Min et al. 2005). A vacuum fraction of 0.7 has been used, which breaks the perfect symmetry of the spherical compact grain. Though we do not focus on detailed feature shapes this choice therefore does not influence our results.

The dust size range is from a_min = 1μm up to a_max = 1 mm and follows a power-law distribution with a_pow =  −3.5. The surface density of the gas and dust decreases proportional to r^-1. Typically in the ISM, the fraction of the mass locked up in PAHs relative to the total dust mass is 4 × 10^-2 (Draine & Li 2007). A comparison with the PAH abundance in protoplanetary disks may be misleading because gas and dust in disks are expected to be processed. Given the uncertainties in abundances, we find that a PAH-to-dust mass fraction of 5.0    ×    10^-4 fits the PAH luminosity of HD 97048. We use a single size PAH molecule with N_C = 100, which has a radius of a_PAH = 6 Å. The assumption of PAH-size may be arbitrary because the ionization depends weakly this choice. Though this paper only intends to demonstrate the effect of ionization, and the choice of a typical PAH of N_C = 100 therefore seems justified.

Figure 7 presents six diagnostic plots of the benchmark model. The black dashed-dotted line gives the radial τ_UV = 1 surface at 0.1 μm. This is the location where UV-photons from the central star are scattered inward in the disk. The gray dashed line gives the vertical τ_MIR = 1 surface at 10 μm. Above this line, MIR photons can escape the disk. Therefore, it shows the region of the disk where the PAH emission originates. The upper four plots give, respectively, the strength of the UV field G₀, the distribution of the electron number density n_e, the temperature T in the disk, and the ionization parameter γ given by G₀T^1/2/n_e. The ionization parameter γ is a representation of the physical conditions, which control the ionization balance of PAHs. Thus, at the radial τ_UV = 1 surface between 10 AU and 100 AU, γ ~ 10³ − 10⁴. The bottom left plot shows the PAH intensities (arbitrary scaling) responsible for the observed emission. The PAH intensity shows the relative contributions of PAH emission to the final spectrum. The region where the PAH contribution is highest lies between the radial τ_UV = 1 and vertical τ_MIR = 1 surfaces. The solid red, black and blue lines give the locations, where the neutral fractions f(0) are, respectively 0.25, 0.5 and 0.75. The neutral fractions are computed in each grid cell of the disk using Eqs. (4) and (5). From this, we can see that most of the PAH emission originates in values below the f(0) = 0.75 line and, thus, have a high neutral fraction. The bottom right plot shows the intensity of the integrated PAH emission per radial annulus (arbitrary units). The PAHs between ~20−200 AU contribute most to the final spectrum. The blue and red lines give the decomposed neutral and ionized contributions. Of the PAH molecules, 93% in the benchmark flaring disk model are neutral.

An important insight from the benchmark model is that neutral PAHs dominate the emission at all radii in the disk. There is, thus, no gradient in the neutral fraction of PAHs as a function of radius in an optically thick disk. Only the surface layer of the disk produces ionized PAH emission but that contribution is ~1−2 orders of magnitude smaller than the contribution from neutral PAHs from deeper layers in the disk. The effect of a higher UV field does not translate to a higher degree of ionization because the electron density is also higher closer to the star.

The PAH luminosity of the benchmark model is L_PAH/L_∗ = 8.40 × 10^-3 and band strength ratio I_6.2/I_11.3 = 1.35. The results are shown on the top of Fig. 8, where a direct comparison can be made to the observations of our sample of Herbig stars, as shown on the bottom plot. As a next step, we investigate the effects of varying the star and disk parameters of the benchmark model on L_PAH and I_6.2/I_11.3. We discuss the parameters which have a significant effect on the charge state of the PAHs.

Fig. 7 Diagnostic plots of benchmark model. The dashed black line is the vertical τMIR = 1 surface at 10 μm, the dash-dotted grey line gives the radial τUV = 1 surface at 0.1 μm. Upper left: strength of the UV field G0. Upper right: electron number density ne. Middle left: temperature distribution in the disk. Middle right: ionization parameter γ given by G0T1/2/ne. Bottom left: PAH intensity of a 2D cut through the disk. This diagnostic plot shows the origin of the PAH flux in the disk, i.e. the contributions of PAH emission to the final spectrum. The intensity of the PAH flux is shown in logarithmic scale (arbitrary scaling), where each color spans an order of magnitude. The solid red, black and blue lines give the locations where the neutral fractions f(0) are respectively 0.25, 0.5, 0.75. The neutral fractions are computed by Eqs. (4) and (5) in each grid cell of the disk. PAHs in the surface of the disk are largely ionized, while PAHs in the mid-plane are neutral. Bottom right: intensity of PAH emission per radial annulus (arbitrary units). The blue and red lines give the neutral and ionized contributions. 93% of the PAH molecules in the benchmark flaring disk model are neutral. The contribution from neutral PAHs dominate irrespective of the radial location in the disk. The noise is due to photon statistics. The height θ relates to Z by, θ = tan   (Z/r).

Fig. 8 PAH luminosities (LPAH/L∗) compared to the I6.2/I11.3 ratio for the benchmark modeling grid (left) and the observations of Herbig stars (right). Several parameters of the benchmark model are varied and their effect to the PAH luminosity and I6.2/I11.3 ratio is shown. The main results of the benchmark modeling grid study is that in general the PAH neutral fraction is high and ionized PAHs originate from optically thin, low density environments. The red dots show the objects studied in this paper. The parameters of all the objects shown in this right plot are listed in Table 1.

4.1. Stellar properties

We varied the stellar temperatures and luminosities of the benchmark model. A decrease in stellar temperature and luminosity lowers the neutral fraction because the UV field G₀ decreases and fewer UV photons are available to excite the PAHs. Therefore, the I_6.2/I_11.3 ratio goes down. For a higher luminosity and/or stellar temperature, the UV field G₀ and, consequently, the I_6.2/I_11.3 ratio are higher. Figure 8 shows that our chosen range of stellar temperatures and luminosities result in band strength ratios I_6.2/I_11.3 between ~1−2.2. This range covers 11 of the 18 Herbig stars shown on the right of Fig. 8. None of the other 7 objects with higher I_6.2/I_11.3 ratios have temperatures and/or luminosities higher than 20 000 K and/or 400 L_⊙. Thus, for these objects, the stellar temperature and luminosity cannot explain the high band strength ratio I_6.2/I_11.3.

For stars with a lower temperature and/or luminosity, the UV density decreases and the contribution from PAHs heated by multi-photon events is lower. As PAHs heated by multi-photon events have a higher temperature, they also have a slightly higher I_6.2/I_11.3 ratio and vice versa. In the models with a lower stellar luminosity and temperature all the PAHs in the disk are neutral, although the I_6.2/I_11.3 ratio is slightly lower (~0.1 − 0.2) than the 100% neutral benchmark model (see Fig. 8). Thus, this difference is a result of the fraction of PAH emission heated by multi-photon versus single-photon events. However the difference is small and unimportant compared to the effect of ionization.

4.2. Grain properties

To test the effect of the grain size distribution, we change the dust density power-law to a_pow =  −3.0 and a_pow =  −4.0. For a relatively lower abundance of small grains (a_pow =  −3.0), the average dust opacity is lower and the disk is less opaque. As a result, the PAH luminosity increases, since the optical depth in the disk is lower, and the region where PAHs can get excited is larger. For a higher abundance of small grains (a_pow =  −4.0), the PAH luminosity decreases, since the optical depth in the disk is higher and the region where PAHs can get excited is smaller. The size distribution of dust grains has a strong effect on the PAH luminosity through the dust optical depth in the UV. If small dust grains are depleted, then the UV opacities in the disk are dominated by PAHs and a small amount of PAHs can give rise to a strong PAH luminosity. The effect on the degree of ionization is negligible since the ionization parameter γ and the neutral fraction f(0) change accordingly to the optical depth and the fraction of ionized to neutral PAHs stays roughly similar.

4.3. PAH properties

Since evolutionary processes may affect the PAH composition in the disk, it is difficult to constrain PAH sizes and abundances. In our benchmark model, we assume a PAH mass of 5.0 × 10^-4  M_dust. For an order of magnitude higher and lower PAH mass fraction, the PAH luminosities change by factors of 2.7 and 0.21 respectively. The PAH mass fraction may then explain the wide spread in PAH luminosity among Herbig stars (~3 orders of magnitude, see Table 1)

Figure 8 shows a small shift in the I_6.2/I_11.3 ratio for factor of 10 higher or lower PAH-to-dust mass fraction. This is also because the PAH-to-electron abundance is then altered, which changes the amount of electrons per PAH molecule available for recombination. Because PAH sizes in disks are expected to be in the order of N_C ~ 50 − 150 (Tielens 2008), the difference in PAH sizes have a very small effect on the final spectra because the neutral fraction depends only on the square root of the PAH size in our ionization model (Eq. (5)).

4.4. Disk mass

We test the effect of lower disk masses on the PAH characteristics by decreasing the total disk mass from 5    ×    10^-2  M_⊙, which is four orders of magnitude down to 5    ×    10^-6  M_⊙. For a less massive disk, the disk vertical height and, consequently, the fraction of the stellar luminosity captured by the disk are lower. Therefore, the PAH luminosity drops with decreasing disk mass, even when the disk is still optically thick.

The charge state increases when the mass of the disk becomes so small (M_disk ≲ 5 × 10^-5  M_⊙) that the disk becomes optically thin. At that point, the UV field G₀ is high even in the mid-plane of the disk, while the electron density n_e is lower. As a result, the fraction of ionized PAHs in the benchmark model with a disk mass of 5 × 10^-6  M_⊙ is almost 100% and the I_6.2/I_11.3 ratio is 7.82. This low disk mass model and the Herbig star HD 141569 have a similar L_PAH and I_6.2/I_11.3 (see also Table 1). This is consistent with our modeling result because HD 141569 is in the transition of a gas-rich protoplanetary disk to a gas-poor debris disk and largely optically thin (Li & Lunine 2003; Thi et al. 2014).

4.5. Summary

The primary goal of our study is to identify the most important trends in the behavior of L_PAH and the I_6.2/I_11.3 ratio. We have presented a small parameter study around the benchmark model that investigates the behavior of PAH emission as a function of disk properties. We have covered the disk parameters which have a significant effect on the charge state and/or luminosity of the PAHs; the results are shown in Fig. 8. The main insights are the following: 1) the neutral fraction is high (f(0) ~ 0.8 − 0.9) for optically thick disks with properties typical for Herbig stars; 2) the PAHs in an optically thick flaring disk are predominantly neutral at all radii. Thus the PAH neutral fraction shows no gradient with radius; 3) the total PAH luminosity increases if the UV field is stronger if the PAH-to-gas mass fraction is higher and if the dust density is lower; 4) the PAH ionization increases if the UV field is stronger and if the PAHs originate from low-density optically thin environments, where the electron density is low and the UV field is high.

5. Transitional disks

In this section, we show that fitting the PAH spectra of four transitional disks requires a significant contribution of ionized PAHs in gas flows through the gap. We perform quantitative fits to the PAH luminosities L_PAH/L_∗ and band strength ratios I_6.2/I_11.3 of four transitional disks. We fit I_6.2/I_11.3 and L_PAH/L_∗ as close as possible to the observed values with a maximum discrepancy of a factor of 2. The presented PAH models demonstrate that a physical disk solution predicting an enhanced contribution of ionized PAHs from low-density, optically thin regions is successful in explaining the observed spectra. It is also consistent with the available spatial information shown in Fig. 1.

We adopt the disk structures of HD 97048, HD 169142, HD 135344 B, and Oph IRS 48 derived by Maaskant et al. (2013) (see Table 2). All disks have large dust-depleted gaps and flaring outer disks. These large gaps are inferred from analysis of Q-band direct imaging. The absence of the 10- and 20-μm silicate feature indicates that the gaps are radially very large and empty in small dust grains. All disks show a NIR excess indicative of a hot (~1500 K) component. To fit the SED at NIR wavelengths, an inner dust component has been added to the disk. For HD 97048, this inner disk is optically thick, is in hydrostatic equilibrium and is located from 0.3−2.5 AU. For HD 169142, HD 135344 B, and Oph IRS 48, optically thin spherical dust halos are added between 0.1−0.3 AU to account for an even stronger excess of NIR emission. For the outer disks, we assume that the vertical scale height is set by hydrostatic equilibrium. For the dusty disk, we assume a radial dependence of the surface density Σ ∝ r^-1. Therefore, the density is highest in the mid-plane close to the star. If the inner disk is optically thick (as for HD 97048), the inner disk acts as a UV-shield and reduces the radiation field G₀ further out in in the disk.

We use the spatial information from the VLT/VISIR FWHM profile, as shown in Fig. 1, to understand the distribution of PAHs throughout the disk. However, this information is only taken into account qualitatively, as we refrain from doing a detailed fit of the FWHM profile since there are issues with the point spread function (PSF) of the VISIR measurements. In addition, the origin of the N-band continuum is not well understood since it can be produced by very small grains or dust from the inner and/or outer disk (see Maaskant et al. 2013 for an extended discussion of this topic). A qualitative interpretation of the VLT/VISIR data is that the PAHs clearly come from larger scales than the N-band continuum for HD 97048; thus, PAH emission originates in the outer disk (in agreement with Lagage et al. 2006; Doucet et al. 2006). For HD 169142, HD 135344 B and Oph IRS 48, the FWHM in the PAH feature is equal and even smaller than the continuum. If PAH emission would originate only from the outer disk, than the 11.3-μm PAH feature would be resolved with respect to the FWHM of the continuum. Thus we interpret this as a significant contribution of PAH emission that comes from within the dust-depleted gap of the disk.

5.1. Adding low-density gas flows to the disk structure

Our model study in Sect. 4 has shown that the PAH emission from disks is predominantly neutral in character unless the disk mass is very low. In a more general sense, ionized PAHs are produced by high UV fields (or low dust densities) and low electron densities. The fitting of the observed I_6.2/I_11.3 ratio of PAHs in the integrated spectrum of the disk has been done by adding a low-density gaseous disk in the gap. In this way, PAH emission originates not only in predominantly neutral PAHs in the optically thick inner and outer disk but also in ionized PAHs in the disk gap. To model the gas flows though the gap and to prevent PAH self-shielding effects, we assume a radial dependence of the surface density of Σ ∝ r⁰ in the gap. The gas surface densities Σ_gas as a function of radius are shown in Fig. 9. The required surface densities of the gas within the gap are on the order of Σ_gas ~ 10^-3 − 10^-4 g cm^-2, which are consistent with hydrodynamical predictions of gas flows through planets opening a gap (e.g. Bryden et al. 1999; Dodson-Robinson & Salyk 2011; Zhu et al. 2012; Fung et al. 2014). The basic parameter assumptions and results of our best-fit models are given in Table 2. Figure 10 shows the PAH intensities of our derived PAH models as a function of the location in the disk.

Fig. 9 Gas surface densities Σgas for the PAH models as a function of radius. Gas surface densities of Σgas ~ 10-3 − 10-4 g cm-2 in the gaps are consistent with hydrodynamical planet formation models and represent gas flows through protoplanetary gaps. The PAHs in these low-density gaps have a high ionization fraction.

5.2. PAH fit of the transitional disk HD 97048

The PAH luminosity of HD 97048 can be fitted with a PAH-to-dust mass fraction of 5.0 × 10^-4 by giving a luminosity of L_PAH/L_∗ = 9.13 × 10^-3, which is close to the observed 9.01 ± 0.03 × 10^-3 in the Spitzer/IRS spectrum. The averaged neutral fraction is f(0) = 0.97; thus, only 3% of the PAHs in the disk are ionized. We note that the calculated I_6.2/I_11.3 ratio (1.27) is larger than the observed one (0.87). However, that may well be an “artefact” of our adopted PAH properties as the intrinsic I_6.2/I_11.3 ratio of neutral PAHs in the Draine & Li (2007) model is arbitrarily set. Figure 8 shows that the benchmark model in which 100% of the PAHs are neutral, which still has a higher I_6.2/I_11.3 fraction (0.94) than HD 97048. Thus, a better fit to the the observed I_6.2/I_11.3 ratio of HD 97048 is not possible because the adopted PAH properties do not cover lower I_6.2/I_11.3 ratios.

The gap of HD 97048 may not be completely empty because no gap is seen in polarized scattered light (Quanz et al. 2012). Therefore, we fill the gap of HD 97048 with a gas surface density which of Σ_gas = 2.5 × 10^-4 g cm^-2, which is similar to the solutions we find for the other transitional disks in our sample (in the next Sect. 5.3). Due to the optically thick inner disk, most of the material in the gap is in the shadow of the inner disk. Thus, the PAHs in the disk gap only contribute 5% to the total PAH luminosity. A small fraction of dust may also be present in the gap, which may explain the scattered light observations of Quanz et al. (2012). However, we do not further investigate that scenario since that is not the focus of our paper. The spectra of the PAH components from the optically thick dust disks are shown by the blue lines in Fig. 11. The red lines show the PAH emission components from the low-density optically thin gas flows in the gaps.

There are two possibilities why HD 97048 is largely neutral. One explanation is because less (ionized) PAH emission originates in the gap. This occurs because they are in the shadow of the optically thick inner disk. The second explanation is because the contribution from PAHs in the large flaring outer disk is still large and dominating (compared to the other disks in our sample which have smaller outer disk radii). We have run several test models without optically thick inner disks, test models with a higher surface density in the gap, and test models with smaller outer disk radii. The results suggest that the first explanation is most important, because models without an optically thick inner disk show a strong increase of L_PAH and I_6.2/I_11.3 due to the fact that (ionized) PAHs in the gap are now directly exposed to the UV emission. Since this is not observed, we suggest that the inner disk of HD 97048 is still optically thick.

5.3. PAH fits of the transitional disks HD 169142, HD 135344 B and Oph IRS 48

The model solution of HD 97048 does not work for HD 169142, HD 135344 B, and Oph IRS 48 for two reasons. First, the model is not consistent with the FWHM profiles because the observed PAH features are not extended compared to the continuum. Second, the observed I_6.2/I_11.3 ratios are too high to be fitted by neutral PAHs from the outer disk suggesting that a significant contribution comes from ionized PAHs. To fit the PAH spectra of HD 169142, HD 135344 B, and Oph IRS 48 we increase the contribution from the low density, optically thin gas in the disk gaps. We keep the PAH-to-dust fraction in the outer disk similar to that of HD 97048 (f_PAH,disk = 5.0 × 10^-4).

Fig. 10 Diagnostic plots showing the origin of the PAH flux in the disk, i.e. the contributions of PAH emission to the final spectrum. The intensity of the PAH flux is shown in logarithmic scale (arbitrary scaling), where each color spans an order of magnitude. Top left: HD 97048, top right: HD 169142, bottom left: HD 135344 B, bottom right: Oph IRS 48. The red, black and blue solid lines give the locations where the neutral fractions f(0) are respectively 0.25, 0.5, 0.75. The neutral fractions are computed by Eqs. (4) and (5) in each grid cell of the disk. PAHs in the surface of the disk are largely ionized, while PAHs in the mid-plane are neutral. The dashed black line is the vertical τMIR = 1 surface at 10 μm, the dash-dotted grey line gives the radial τUV = 1 surface at 0.1 μm. The height θ relates to Z by, θ = tan   (Z/r).

To obtain a high PAH luminosity and low dust continuum, we assume that the opacities are dominated by the PAHs in the gap. For that reason, we assume that there are no small (~1  μm) grains and, thus, no dust is visible in the spectrum. To get a lower limit to the gas mass needed in the disk gap, we adopt a PAH mass fraction similar to the ISM of f_PAH,gap = 4.0 × 10^-2. For HD 169142, HD 135344 B and Oph IRS 48, the luminosity L_PAH/L_∗ and the I_6.2/I_11.3 ratio are simultaneously fitted by adjusting the total amount of mass in the gas disk in the gap. The characteristics are given in Table 2. We have aimed to fit L_PAH/L_∗ and I_6.2/I_11.3 as close to the observed values as possible. We find that models could be found, which fit the observed L_PAH/L_∗ and I_6.2/I_11.3 within a factor of 2 for all the objects. Following this fitting procedure, PAHs in the gas flows through the gap produce respectively 33%, 73%, and 91% of the total PAH luminosity of HD 169142, HD 135344 B, and Oph IRS 48.

The diagnostic plots in Fig. 10 show that the disk gap of HD 169142 has a high fraction of ionized PAHs, while the disk gap of HD 135344 B has a much larger fraction of neutral PAHs. This difference can be attributed to the lower stellar temperature of HD 135344 B (T = 6590 K) compared to HD 169142 (T = 8200 K) producing a lower UV-field.

5.4. Summary

We have found a model solution to fit L_PAH and the I_6.2/I_11.3 ratio of the transitional disks HD 97048, HD 169142, HD 135344 B, and Oph IRS 48. Observations constrain a contribution of ionized PAHs that arise from within the dust gap (see Fig. 1). Our PAH models combine PAH contributions from the surface of a flaring dust disk and a gaseous disk component inside the gap where the ionization fraction can be high. The density structure of gas and dust in the disk is key to the PAH charge state, since the UV field and electron density are the main parameters that determine the PAH ionization.

We stress that the scale height and surface densities as a function of radius are poorly constrained for the inner regions. However, the presented PAH model solutions do make a constraint that a particular abundance of ionized PAHs in an optically thin, low density environment is required to fit the spectra. We have chosen to present the models for which the vertical density structure is set by hydrostatic equilibrium. Alternatively we can parameterize the density structure in the disk gap and adopt a higher scale-height of the gas in the gap to obtain similar results. For example, an even better fit to the observed L_PAH/L_∗ can be obtained by artificially increasing the scale height of the gas in the gap for Oph IRS 48 (which may be an effect of a larger gas temperature due to UV irradiation). A fit to the observations is then found using a slightly higher mass M_PAH,gap = 1.2 × 10^-9  M_⊙ and setting the vertical scale height H₀ = 3.0 AU at R = 30 AU where H ∝ r^1.1. For comparison, the hydrostatic disk structure has a pressure scale height of H₀ = 2.1 AU at R = 30 AU.

Fig. 11 Top and bottom spectra represent the benchmark models with respectively 100% neutral and 100% ionized PAHs. The middle four spectra show transitional disks with the dotted black line showing the Spitzer observations, the solid black lines the total fluxes of the best fit models, the solid blue line the contribution from PAHs in the dust disk, and the solid red line the contribution from the low density optically thin disk. Note that for HD 97048, there is no low density optically thin disk component thus the solid black and blue lines are merged into a dashed blue-black line. For Oph IRS 48, the spectrum is dominated by ionized PAHs from the low density optically thin disk.

6. Discussion

In the following section we discuss our results on PAH ionization in relation to the disk structure. First, we discuss observations and models of gas in protoplanetary disk gaps. Thereafter, we discuss our results in the context of disk evolution in Herbig stars.

6.1. Gas flows through protoplanetary gaps

To fit the observed PAH spectra of the transitional disks in our study, we find that a requirement is to add a “gaseous” low-density, optically thin disk in the gap to account for a higher fraction of ionized PAHs. Observations of gas in disk gaps can be understood in the context of planet formation models. The disk should go through a gas-rich/dust-poor stage if planets form by a two-step process where grains first agglomerate into large rocky cores, followed by the accretion of a gaseous envelope (e.g. Lissauer 1993). As PAHs are small molecules, which are thought to travel along with the gas, we first compare PAHs to other gas diagnostics. Thereafter, we discuss our results in the context of theoretical modeling of gas flows through gaps.

Recent Herschel/PACS observations show that Herbig objects with pure rotational CO transitions have a strong UV luminosity and strong PAH emission (Meeus et al. 2013). However, for a number of other molecular lines, no connection to the PAH luminosity has been reported (Fedele et al. 2013). From these observations, it is difficult to deduce how PAHs and the gas are connected. On the other hand, there is ample evidence for (molecular) gas in disk gaps. In high-resolution M-band spectral studies of transitional disks, modeling of the rovibrational CO emission yields a great diversity of inner disk gas content and of gas/dust ratios (Salyk et al. 2007, 2009). Similar studies have been carried out for HD 135344 B (Pontoppidan et al. 2008) and Oph IRS 48 (Brown et al. 2012) where the respective authors find evidence for significant quantities of CO gas within the dust gap. Salyk et al. (2011) show that the detection efficiency of H₂O, OH, HCN, C₂H₂, and CO₂ depends on the disk color, F₃₀/F_13.5, Hα equivalent width, and, tentatively, the mass accretion rate, Ṁ. This may suggest a connection between regions of lower dust/gas ratio in the disk and the strength of these molecular lines. Recent ALMA observations of HD 142527 have shown that reservoirs of gas are present in the dust depleted regions (Casassus et al. 2013). These observations show that gas flows through dust-poor regions are not uncommon in protoplanetary disks, which is consistent with our results.

Hydrodynamical simulations of the disk structure caused by planet formation predict the opening of an annular “gap”, where the structure and depth of the gap is dependent on the mass of the planet and disk viscosity (e.g., Bryden et al. 1999; Varnière et al. 2004). Mass can penetrate through the gap and flow onto the planet, or flow through the gap and replenish the inner disk (Lubow et al. 1999; Paardekooper & Mellema 2004; Edgar & Quillen 2008; Tatulli et al. 2011). Since pressure gradients at the outer edge of the gap, which are cleared by the planet, act as a filter, this process can also produce a very large gas-to-dust ratio in the inner disk, potentially explaining (pre-)transitional disks, which have optically thin inner cavities but still have relatively high accretion rates (Rice et al. 2006; Mendigutía et al. 2011). Recent work by Pinilla et al. (2012) demonstrates how dust particles are trapped further out than the location of the edge of the gap in the gas for a given planet. For example, the outer edge for the gas would be located at 37 AU in the specific case of a 15 M_Jup planet at 20 AU, while the dust ring would be at 54 AU. In these models, the surface densities of the gas in the gap can drop to Σ_gas ~ 10^-3−10^-4 g cm^-2, which is similar to other studies on the gas densities in the gaps of transition disk caused by the presence of (multiple) planets (e.g., Bryden et al. 1999; Dodson-Robinson & Salyk 2011; Zhu et al. 2012; Fung et al. 2014). In a paper by de Juan Ovelar et al. (2013) using the models of Pinilla et al. (2012), it is demonstrated how multi-wavelength observations of the gap size (i.e., the location of the inner edge of the outer disk) may help to constrain the planet mass. The results of this paper may further help to constrain the physical properties of the gap. The required surface densities in the gaps of our models are Σ_gas ~ 10^-3−10^-4 g cm^-2, which are consistent with the surface densities found in hydrodynamical simulations. The ionization balance of PAHs can, therefore, serve as an excellent diagnostic tool to characterize the gas flows through the gaps.

6.2. PAH ionization and Herbig star disk evolution

Our study shows that neutral PAHs originate in optically thick environments where the UV field is low, and the electron density is high. Ionized PAHs originate in low density, optically thin gas flows where the UV field is high and the electron density is low. Figure 12 gives a summary sketch of these two disk environments in a transitional disk.

The two most extreme I_6.2/I_11.3 ratios among Herbig stars come from HD 141569, where the PAHs are likely almost all ionized, and HD 97048, where the PAH spectrum represents the highest neutral fraction (Table 1). In HD 141569, studies have shown that the disk is largely optically thin (Li & Lunine 2003; Thi et al. 2014), which agrees with our result that ionized PAHs originate in optically thin environments. In HD 97048, neutral PAHs originating in the optically thick disk dominate. Even though there may be gas in the gap of HD 97048 due to the higher density of an optically thick inner disk, the PAHs in the low-density gap are not exposed enough to the UV field.

Fig. 12 Sketch of the main result of this paper. Neutral PAH emission originates in an optically thick disk (blue components) and ionized PAH emission from low density, optically thin gas flows through protoplanetary disk gaps (red component).

The effect of the disk geometry on the ionization balance may be visible in the I_6.2/I_11.3 ratio compared to the MIR spectral index F₃₀/F_13.5 (Fig. 3). Flaring/transitional (group I) objects have a wider range in their degree of ionization compared to flat (group II) objects. This can be interpreted, as that PAH emission in flaring/transitional objects originates in more varied physical environments: those that are predominantly neutral in the optically thick disk and predominantly ionized in optically thin gaseous disk gaps. The PAH emission from flat disks show more similar I_6.2/I_11.3 ratios. This suggests that the disks of flat (group II) objects are more homogeneous. Possibly this could be due to the lack of large gaps in these objects. This would be consistent with the scenario that flaring (group I) objects are often (pre-)transitional, while self-shadowed (group II) objects lack any signatures of large dust gaps.

While ionized PAHs do trace low density gas flows, it is not always the case that gas flows have ionized PAHs. For example, the Herbig stars HD 135344 B and HD 142527 have rather low temperatures (respectively T = 6590 K and T = 6260 K), and, therefore, their UV fields are less efficient in ionizing the PAHs. Our models show that PAHs will be highly neutral in the disks of T Tauri stars. Thus, ionized PAHs are most likely found in low density gaps around Herbig stars with temperatures of T ≳ 6000 K. Alternatively, if the gas density in the gas flow is high, then the PAHs may recombine more easily into the neutral state. This scenario may be part of the explanation why the spectrum of HD 97048 is still largely neutral.

A correlation is found between the mm luminosity and the I_6.2/I_11.3 ratio (Fig. 4). Since the disk is optically thin at millimeter wavelengths, the mm luminosity may be a proxy of the disk mass in mm-sized dust grains. Therefore, for lower disk masses, the fraction of ionized PAHs seems to be higher. There are two possible ways to understand this trend. The first is that the inner region of the disk becomes optically thinner when accretion depletes most of the inner disk mass and less material from the outer disk is available to replenish the inner regions. As a consequence, the inner disk becomes optically thin, the low density gas flows become more irradiated by UV emission, and the ionization fraction is higher. Another possibility is that more planets are formed creating larger gaps over the course of time when the disk is loosing its mass. Consequently, more low density gas flows emerge with a high fraction of ionized PAHs.

7. Conclusions

We have modeled the charge state of PAHs in Herbig disks. We find that the neutral fraction of the disk is high for typical star and disk parameters. The key parameter to ionize the PAHs in the disk is to lower the gas density. The PAH luminosity is most influenced by competition from dust grains, either by an increased PAH mass fraction, or the removal of small dust grains. Our results can be summarised as follows:

• Emission from ionized PAHs trace low-density optically thinregions and has a high I_6.2/I_11.3 ratio.

• The PAHs in an optically thick disk are predominantly neutral and produce emission with a low I_6.2/I_11.3 ratio.

• The PAH spectra of transitional disks can be understood as superpositions of neutral PAHs from an optically thick disk and ionized PAHs from that of low density, optically thin gas flows through protoplanetary disk gaps.

• There is no gradient in the neutral fraction of PAHs as a function of radius in an optically thick disk.

• The PAH luminosity is not synonymous to the disk structure. As shown for the transitional disks in our sample, the total PAH mass in the low density, optically thin disk gap can be much smaller compared to the total PAH mass in the outer disk. Nevertheless, they give a dominating contribution to the total PAH luminosity.

• We find a trend in the I_6.2/I_11.3 ratio compared to the 1.3 millimeter luminosity. Since the 1.3 millimeter luminosity can be used as a proxy for the disk mass, the degree of ionization tends to be higher for lower disk masses. We suggest that the contribution to the PAH luminosity from ionized PAHs in low-density optically thin “gas flows” becomes larger for lower mass disks.

• For a sample of 18 Herbig stars, the I_6.2/I_11.3 ratios indicates that PAH emission from flaring/transitional disks (group I) is higher and originates in a wider range of physical environments (⟨ I_6.2/I_11.3 ⟩  = 2.61 ± 1.63), while the origin of PAH emission in self-shadowed flat disks (group II) seems to be more homogeneous (⟨ I_6.2/I_11.3 ⟩  = 2.08 ± 0.56).

With superior mid-IR sensitivity and spatial resolution, observations by the Mid-InfraRed Instrument (MIRI) on the James Webb Space Telescope (JWST) will be able to probe more and fainter targets. This work shows that the ionization balance of PAHs provides an extra diagnostic tool to characterize the gaps of transitional disks and, thus, to study the process of planet formation.

Acknowledgments

The authors thank the referee Ralf Siebenmorgen for useful comments which helped to improve the paper. The authors thank Bram Acke for providing the reduced data from his 2010 paper. The authors thank Carsten Dominik and Alessandra Candian for insightful comments and discussions which improved the analysis presented in this paper. The authors thank Lucas Ellerbroek for proof-reading and providing very useful comments. The authors thank the editor Malcolm Walmsley for providing useful comments on the final version of this paper. K.M. is supported by a grant from the Netherlands Research School for Astronomy (NOVA). M.M. acknowledges funding from the EU FP7-2011 under Grant Agreement No 284405. Studies of interstellar chemistry at Leiden Observatory are supported through advanced-ERC grant 246976 from the European Research Council, through a grant by the Dutch Science Agency, NWO, as part of the Dutch Astrochemistry Network, and through the Spinoza premie from the Dutch Science Agency, NWO.

References

All Tables

All Figures

Fig. 1 Left: Spitzer/IRS spectra of the four transitional disks HD 97048, HD 169142, HD 135344 B and Oph IRS 48. The spectra are scaled relative to the 6.2-μm PAH peak flux and show a wide range in the I6.2/I11.3 ratio; most neutral in HD 97048 and most ionized in Oph IRS 48. Middle: VLT/VISIR observations show that FWHM as a function of wavelength. The size of the 11.3-μm feature is respectively larger, similar and smaller than the N-band continuum for HD 97048, HD 169142, HD 135344 B, and Oph IRS 48 (data from Maaskant et al. 2013). The diamond symbols show the sizes of the point spread functions of the observations. The dash-dotted black line is the diffraction limit of the telescope. Right: Sketches showing a qualitative interpretation of the spatial distribution of the dominant contribution of PAH emission, as shown by the colored ovals. Each object is given a similar color in the left, middle, and right plots.

In the text

Fig. 2 PAH luminosity (LPAH/L∗) of a sample of Herbig stars compared to the MIR spectral index (F30/F13.5). The MIR spectral index can be used as a tracer of disk gaps, as indicated by the symbols showing the flaring/transitional disks (group I) with silicate features (red circles), without silicate features (orange triangles), and self-shadowed disks (group II) with silicate features (blue squares). While the PAH luminosities among the objects span ~2–3 orders of magnitude, the averaged luminosities of the disk groups only differ by a factor of ~2–3. The parameters of all the objects are listed in Table 1.

In the text

Fig. 3 I6.2/I11.3 ratio compared to MIR spectral index (F30/F13.5). Ionized PAHs are traced by a high I6.2/I11.3 ratio. The wider range in the I6.2/I11.3 ratio in flaring/transitional disks (⟨ I6.2/I11.3 ⟩  = 2.61 ± 1.63) suggests that the PAH emission originates in a wider range of physical environments while the origin of PAH emission in flat disks (⟨ I6.2/I11.3 ⟩  = 2.08 ± 0.56) seems to be more homogeneous. The parameters of all the objects are listed in Table 1.

In the text

Fig. 4 An observational trend between the I6.2/I11.3 ratio as a function of the luminosity at 1.3 mm (L1300/L∗). The red dots show the objects studied in this paper. The gray line gives the best fit to the data (Eq. (1)). The parameters of all the objects are listed in Table 1. For high mm luminosities, the I6.2/I11.3 ratio is lowest (tracing neutral PAHs), while the I6.2/I11.3 ratio is highest (tracing ionized PAHs) for low mm values.

In the text

	Fig. 5 Opacities for neutral PAHs (solid blue line) and ionized PAHs (dashed red line) with NC = 100. Computed using the method described by Draine & Li (2007).
In the text

	Fig. 6 Fraction of neutral PAH molecules f(0) given by the UV field, the temperature, and the electron number density and controlled by the ionization parameter G0T1/2/ne. In our models, we assume a PAH size of NC = 100.
In the text

Fig. 7 Diagnostic plots of benchmark model. The dashed black line is the vertical τMIR = 1 surface at 10 μm, the dash-dotted grey line gives the radial τUV = 1 surface at 0.1 μm. Upper left: strength of the UV field G0. Upper right: electron number density ne. Middle left: temperature distribution in the disk. Middle right: ionization parameter γ given by G0T1/2/ne. Bottom left: PAH intensity of a 2D cut through the disk. This diagnostic plot shows the origin of the PAH flux in the disk, i.e. the contributions of PAH emission to the final spectrum. The intensity of the PAH flux is shown in logarithmic scale (arbitrary scaling), where each color spans an order of magnitude. The solid red, black and blue lines give the locations where the neutral fractions f(0) are respectively 0.25, 0.5, 0.75. The neutral fractions are computed by Eqs. (4) and (5) in each grid cell of the disk. PAHs in the surface of the disk are largely ionized, while PAHs in the mid-plane are neutral. Bottom right: intensity of PAH emission per radial annulus (arbitrary units). The blue and red lines give the neutral and ionized contributions. 93% of the PAH molecules in the benchmark flaring disk model are neutral. The contribution from neutral PAHs dominate irrespective of the radial location in the disk. The noise is due to photon statistics. The height θ relates to Z by, θ = tan   (Z/r).

In the text

Fig. 8 PAH luminosities (LPAH/L∗) compared to the I6.2/I11.3 ratio for the benchmark modeling grid (left) and the observations of Herbig stars (right). Several parameters of the benchmark model are varied and their effect to the PAH luminosity and I6.2/I11.3 ratio is shown. The main results of the benchmark modeling grid study is that in general the PAH neutral fraction is high and ionized PAHs originate from optically thin, low density environments. The red dots show the objects studied in this paper. The parameters of all the objects shown in this right plot are listed in Table 1.

In the text

	Fig. 9 Gas surface densities Σgas for the PAH models as a function of radius. Gas surface densities of Σgas ~ 10-3 − 10-4 g cm-2 in the gaps are consistent with hydrodynamical planet formation models and represent gas flows through protoplanetary gaps. The PAHs in these low-density gaps have a high ionization fraction.
In the text

Fig. 10 Diagnostic plots showing the origin of the PAH flux in the disk, i.e. the contributions of PAH emission to the final spectrum. The intensity of the PAH flux is shown in logarithmic scale (arbitrary scaling), where each color spans an order of magnitude. Top left: HD 97048, top right: HD 169142, bottom left: HD 135344 B, bottom right: Oph IRS 48. The red, black and blue solid lines give the locations where the neutral fractions f(0) are respectively 0.25, 0.5, 0.75. The neutral fractions are computed by Eqs. (4) and (5) in each grid cell of the disk. PAHs in the surface of the disk are largely ionized, while PAHs in the mid-plane are neutral. The dashed black line is the vertical τMIR = 1 surface at 10 μm, the dash-dotted grey line gives the radial τUV = 1 surface at 0.1 μm. The height θ relates to Z by, θ = tan   (Z/r).

In the text

Fig. 11 Top and bottom spectra represent the benchmark models with respectively 100% neutral and 100% ionized PAHs. The middle four spectra show transitional disks with the dotted black line showing the Spitzer observations, the solid black lines the total fluxes of the best fit models, the solid blue line the contribution from PAHs in the dust disk, and the solid red line the contribution from the low density optically thin disk. Note that for HD 97048, there is no low density optically thin disk component thus the solid black and blue lines are merged into a dashed blue-black line. For Oph IRS 48, the spectrum is dominated by ionized PAHs from the low density optically thin disk.

In the text

	Fig. 12 Sketch of the main result of this paper. Neutral PAH emission originates in an optically thick disk (blue components) and ionized PAH emission from low density, optically thin gas flows through protoplanetary disk gaps (red component).
In the text