A&A 432, 515-529 (2005)
DOI: 10.1051/0004-6361:20040331

PAH charge state distribution and DIB carriers: Implications from the line of sight toward HD 147889[*],[*]

R. Ruiterkamp 1 - N. L. J. Cox 2 - M. Spaans 3 - L. Kaper 2 - B. H. Foing 4 - F. Salama 5 - P. Ehrenfreund 2


1 - Leiden Observatory, PO Box 9513, 2300 RA Leiden, The Netherlands
2 - Institute Anton Pannekoek, Kruislaan 403, 1098 SJ Amsterdam, The Netherlands
3 - Kapteyn Astronomical Institute, PO Box 800, 9700 AV Groningen, The Netherlands
4 - ESA Research and Scientific Support Department, PO Box 299, 2200 AG Noordwijk, The Netherlands
5 - Space Science Division, NASA Ames Research Center, Mail Stop 245-6, Moffett Field, California 94035, USA

Received 25 February 2004 / Accepted 29 September 2004

Abstract
We have computed physical parameters such as density, degree of ionization and temperature, constrained by a large observational data set on atomic and molecular species, for the line of sight toward the single cloud HD 147889. Diffuse interstellar bands (DIBs) produced along this line of sight are well documented and can be used to test the PAH hypothesis. To this effect, the charge state fractions of different polycyclic aromatic hydrocarbons (PAHs) are calculated in HD 147889 as a function of depth for the derived density, electron abundance and temperature profile. As input for the construction of these charge state distributions, the microscopic properties of the PAHs, e.g., ionization potential and electron affinity, are determined for a series of symmetry groups. The combination of a physical model for the chemical and thermal balance of the gas toward HD 147889 with a detailed treatment of the PAH charge state distribution, and laboratory and theoretical data on specific PAHs, allow us to compute electronic spectra of gas phase PAH molecules and to draw conclusions about the required properties of PAHs as DIB carriers. We find the following. 1) The variation of the total charge state distribution of each specific class (series) of PAH in the translucent cloud toward  HD 147889 (and also of course for any other diffuse/translucent cloud) depends strongly on the molecular symmetry and size (number of $\pi $ electrons). This is due to the strong effects of these parameters on the ionization potential of a PAH. 2) Different wavelength regions in the DIB spectrum are populated preferentially by different PAH charge states depending on the underlying PAH size distribution. 3) The PAH size distribution for HD 147889 is constrained by the observed DIB spectrum to be Gaussian with a mean of 50 carbon atoms. 4) For the given PAH size distribution it is possible to constrain the total small catacondensed PAH column density along the line of sight to HD 147889 to $2.4\times10^{14}$ cm-2 by comparing the total observed UV extinction to the strong UV absorptions of neutral PAHs in the 2000-3000 Å region. 5) Catacondensed PAHs with sizes above some 40 C-atoms are expected to show strong DIBS longward of 10 000 Å. Large condensed PAHs in the series, pyrene, coronene, ovalene, .... , on the other hand, mainly absorb between 4000 and 10 000 Å but extrapolation to even larger pericondensed PAHs in this series also shows strong absorptions longward of 10 000 Å. 6) Only the weak DIBs in HD 147889 could be reproduced by a mix of small catacondensed PAHs (<50 C atoms) while for large pericondensed PAHs (50 < C atoms < 100) the intermediate DIBs are well reproduced. Small catacondensed PAHs cannot contribute more than 50% of the total observed equivalent width toward HD 147889. Strong DIBs can only be reproduced by addition of very specific PAH molecules or homologue series to the sample set (i.e., a small number of PAHs with high oscillator strength or a large number of PAHs with a low oscillator strength). An outline is provided for a more general application of this method to other lines of sight, which can be used as a pipeline to compute the spectroscopic response of a PAH or group of PAHs in a physical environment constrained by independent (non-DIB) observations.

Key words: ISM: clouds - ISM: abundances - ISM: individual objects: HD 147889 - ISM: molecules

1 Introduction

1.1 General considerations

The properties and nature of the carrier(s) of diffuse interstellar bands (DIBs) have so far resisted unambiguous identification (Herbig 1995). An important step forward has been the PAH hypothesis (Leger & d'Hendecourt 1985; van der Zwet & Allamandola 1985; Salama et al. 1996), in which large polycyclic aromatic hydrocarbons (PAHs) are responsible for (at least some) of the many DIB features that are observed toward many diffuse and translucent clouds. In this paper we try to test this hypothesis for the line of sight toward HD 147889.

The DIB phenomenon can be roughly divided into two separate lines of research. First, laboratory experiments provide ionization potentials, electron affinities, absorption band strengths and wavelengths, etc. for possible carrier molecules, such as PAHs, fullerenes or carbon chains. Second, any given candidate will contribute to the global DIB spectrum depending on its abundance as well as ionization and hydrogenation state which in turn is determined by the overall chemical and thermal balance of the interstellar cloud under consideration.

Hence, in order to apply the laboratory results, one needs accurate information on the physical characteristics of interstellar gas. In this work we try to combine experimental data on a number of PAHs with an accurate description of the physical state of an interstellar cloud, e.g., density, electron abundance, and temperature. By constraining the physical parameters of the line of sight as well as possible we can draw conclusions on the required microscopic properties of DIB carriers, such as ionization potential, molecular structure, and number of $\pi $ electrons.

1.2 Model constraints on DIB properties

There are thus two ways in which we can constrain the properties of DIB carriers. In the specific case of HD 147889 a large observational data set is modeled. Computed physical properties of the cloud, e.g., molecular content, temperature, density, internal UV flux, and degree of ionization; all as a function of extinction, can be used to determine the charge balance of PAH molecules along the line of sight. If PAH cations are responsible for (some of) the DIB features (van der Zwet & Allamandola 1985; Leger & d'Hendecourt 1985) then we would expect their ionization balance to be crucial in the explanation of DIB properties.

A more general approach is to model the strength of the internal UV field for a wide range of radiative models that encompass all diffuse and translucent clouds in terms of extinction and mean density, temperature, and electron abundance. Subsequently, we can compute the ionization stages of PAHs and other putative carriers for which ionization potentials and electron affinities are known (based on group symmetric properties) in a statistical manner. In this work we will concentrate exclusively on the former approach since an excellent data set is available for the line of sight toward HD 147889 that allows us to constrain the physical properties of the cloud quite well and to compute electronic spectra.

1.3 Ionization fraction in diffuse clouds

Over the last decades many authors have worked on the challenges involved in calculating the reaction rates for grain and PAH chemistry in the ISM. It has become evident that PAHs constitute an important part of the complex ensemble of molecules in the many different astronomical environments observed, like HII regions, reflection nebulae, molecular, and diffuse clouds. Omont (1986), via infra-red emission, estimated the PAH abundance in the ISM to be of the order of 10-7 (per H atom). PAHs are estimated to be the most abundant molecules in the ISM after H2 and CO.

From previous studies (e.g., Le Page et al. 2003; Allain et al. 1996a,b; Bakes & Tielens 1994; Omont 1986, and references therein) it is clear that hydrogenation and ionization are important for interstellar PAH chemistry. The PAH ionization fraction, or more extensively, the PAH charge state (cation, neutral and anion) distribution is important for the proposition put forward by van der Zwet & Allamandola (1985) and Leger & d'Hendecourt (1985) that the positively charged PAHs could be the carriers of the visible DIBs. Omont (1986) noted that for diffuse and translucent clouds the typical values for the UV radiation field ($\chi$ in units of 108 photons cm-2 ster-1 s-1 or $I_{\rm UV}$ in units of the Draine field), the electron density $n_{\rm e}$, and the gas temperature T yield environmental conditions that are favorable to the ionization of PAHs (i.e., if $\chi$/ $n_{\rm e} \leq 10^4/$T1/2). Therefore, in diffuse and translucent clouds, small effects in the environmental parameters can have significant consequences for the specific PAH charge state distribution. In these photoionization-dominated regions, the PAH charge is mostly positive, and also limited to a single elementary charge e. Nevertheless, we feel that it is prudent to include the negatively charged PAHs. In Sect. 4 we will derive the 3-charge state distribution for PAHs in the translucent cloud HD 147889.

In the remainder of this paper we will only discuss the effects of cloud parameters and PAH characteristics on the charge state distribution. Therefore, for simplicity and clarity, we refrain from implementing the hydrogenation processes. Note that, although hydrogenation is important in the carbon chemistry of the ISM, it does not greatly alter the calculations for the charge state distribution.

This paper is organized as follows. Section 2 summarizes the observations available for HD 147889. In Sect. 3 we present a model that reproduces most of the atomic and molecular line data for HD 147889 and thus constrains the ambient physical parameters. The PAH charge balance in HD 147889 for various point groups is discussed in Sect. 4. We refer the reader to the Appendix for the details on the adopted reaction rates and the ionization potential derivations. Section 5 presents electronic PAH spectra for HD 147889 based on the results of the previous sections and for various PAH size distributions. The basic PAH data set and the condensed homologue series for the cations of pyrene, ovalene and coronene are considered. These model spectra are compared to the available DIB observations of HD 147889 in Sect. 6. Section 7 summarizes and discusses our results.

Table 1: Properties of HD 147889 and telluric standard.

Table 2: LSR velocities ( $V_{\rm LSR}$), temperatures ($T_{\rm A}$) and integrated temperatures ($\int$ $T_{\rm A}$ dV) of observed radio lines toward HD 147889. For references see footnote of Table 3. The rightmost column gives the values from our model as described in Sect. 3.

2 The translucent cloud toward HD 147889

HD 147889 is a B2V star, with a luminosity $L_{\star}\approx 5300~ L_{\odot}$, V = 7.9, B-V = 0.83, distance d = 136 pc, and an effective temperature of 22 000 K (Liseau et al. 1999). HD 147889 illuminates the western edge of the $\rho$ Oph cloud. It is situated in a cavity in its surroundings created by its own eroding UV radiation. The translucent cloud toward HD 147889 has a reddening of E(B-V) = 1.09 mag Federman & Lambert (1988) quote for HD 147889: $n_{\rm H} = 1000$ cm-3 and T = 25 K. This renders this line of sight quite opaque so that we expect the internal UV radiation field and hence also the abundances of the PAH charge states to vary significantly with depth . The gas phase carbon abundance toward HD 147889 is assumed to be ${\approx} 1.4\times 10^{-4}$, the typical diffuse ISM value (Cardelli et al. 1996). Our own cloud model for HD 147889 (Sect. 3) gives slightly different values for the physical parameters. This star is located such that only a single cloud component is in the line of sight, giving the unique opportunity to easily model the density and radiation environment in the cloud.

2.1 Observational overview

The region toward HD 147889 has been observed ubiquitously by many authors (cf. Tables 3 and 6). Recently, in September 2001, HD 147889 was also observed with the Very Large Telescope (VLT) using the ultra-violet echelle spectrograph (UVES) as well as, in April 1999, with the ESO 1.52-m telescope with the fiber-fed extended range optical spectrograph (FEROS). In Table 1 we list some basic properties of the target (HD 147889) and the unreddened standard star that is used to correct for telluric lines. With the FEROS instrument we obtained spectra with a signal-to-noise (S/N) up to 160, and a resolution R of about 75 000. Complementary and partly overlapping UVES spectra have ${\it S/N}\sim$ 200-300, and $R\sim115~000.$ Unfortunately, for the spectral range longwards of 6500 Å the UVES spectra suffer significantly from small scale fringing. This fringing effect is known to appear in UVES red arm exposures of bright objects. Most of the region degraded by fringing was also covered with FEROS, for which the spectra do not show any such fringing. In addition, there are numerous telluric lines longwards of 7000 Å. For most of the spectral range the telluric lines could be removed successfully through division by the unreddened reference star. Lines that were obscured before suddenly appear, although the achieved S/N is lower.

In Tables 2 and 3 we list literature values for abundances (and/or equivalent widths) of a large number of species observed in the translucent cloud toward HD 147889. We re-measured the strongest DIBs in the spectrum and listed the equivalent widths (including those from the literature) in Table 6.

Table 3: Abundances of species observed toward HD 147889. For each of the observed species we list the central wavelength of the line in air ($\lambda $ line), equivalent width ( EW) and the observed and predicted column densities (N). The rightmost column gives the values from our model as described in Sect. 3 (see footnote a).

  
3 Physical and chemical models for the line of sight toward HD 147889

Models for the chemical and thermal balance of HD 147889 were constructed using the Monte Carlo code as described in Spaans (1996) and Spaans & van Dishoeck (1997). The reader is referred to these papers for a detailed description. Here a summary is provided of the physical processes included. The chemical network is based on the UMIST compilation (cf. Le Teuff et al. 2000) and is driven by ion-molecule reactions, typical of diffuse and translucent clouds. The self-shielding transitions of H2 and CO (as well as their isotopes) have been included in detail (van Dishoeck & Black 1988). Heating is provided by photo-electric emission from dust grains (Bakes & Tielens 1994) and cosmic rays. Cooling by fine-structure lines of [C  I] 609, [C  II] 158 and [O  I] 63 $\mu$m as well as rotational lines of CO has been included under statistical equilibrium with optical depth effects treated according to an escape probability formalism (Spaans et al. 1994). The most recent, but historically always uncertain, H3+ recombination rate of McCall et al. (2003) has been used.

The parameters of the models considered here are the total hydrogen density, assumed to be constant and given by $n_{\rm H}=n(H)+2n({\rm H}_2)$, the strength of the impinging UV radiation field in units of the Draine radiation field, $I_{\rm UV}$, the abundances of metals and PAHs, the total visual extinction, $A_{\rm V}$, dust properties (opacity, albedo en phase function) and the geometry of the cloud. In all cases the Henyey-Greenstein phase function and grain parameters have been adopted (Roberge et al. 1991). For HD 147889 an observed extinction curve (Karl Gordon, private communication) is available with $R_{\rm V} = 4.22$. Note that the observed $R_{\rm V}$ value is substantially different from the generic value of 3.1. The observed dust properties are used for the detailed radiative transfer model. The model clouds are assumed to be homogeneous and symmetric, and have a total visual extinction of 4.6 mag, i.e., 2.3 mag to the center.

Observational data suggest that the line of sight toward HD 147889 is a single cloud that has a fairly simple kinematic structure, i.e., it has one velocity component, based on Na I and K I (see Fig. 1) and CO measurements (van Dishoeck et al. 1991; Gredel et al. 1994). Therefore, the assumed model properties of the cloud are justified.


  \begin{figure}
\par\includegraphics[angle=-90,width=8.3cm,clip]{0331fig1.eps}\end{figure} Figure 1: Normalized UVES spectra of the Na I (D1 and D2) and K I (7698.98 Å) interstellar lines toward HD 147889. The Na I line is saturated (from the D1 and D2 line ratio) but the K I line is not. Both lines show only one narrow velocity component, indicative of a single cloud structure in this line of sight.
Open with DEXTER

In order to assess the goodness of the model fit a grid of models is constructed for different values of $n_{\rm H}$ and $I_{\rm UV}$. The constraints for the fit are provided by the observed intensities of 12CO, 13CO and C18O, as well as CH and CN (see Tables 3 and 6, and references therein). The quality of the fit is computed by normalization of the sum of the squared deviations. It is found that the UV radiation field and the ambient density are constrained quite well, within 20%. In all, given the large body of data that can be reproduced quite well, with the exception of the generally anomalous CH+ abundance, our characterization of the physical conditions toward HD 147889 is reliable.

In Tables 3 and 6 we compare the observed and modeled abundances for some key species present in the cloud. The agreement is quite good for HD 147889 and our model gives: $n_{\rm H} = 1200~\pm~200$ cm-3, $I_{\rm UV} = 11~\pm~2$. The temperature ( $T_{{\rm gas}}$) and electron density ($x_{\rm e}$) gradients are shown in Fig. 2. These values will be adopted in our charge state distribution calculations for HD 147889 (Sect. 4 and onwards) and we believe them to be accurate enough to draw conclusions on the PAH charge state balance along the line of sight toward HD 147889.

  
4 The PAH charge state distribution in HD 147889

In our model we calculate the charge state distribution (fractions of anions, neutrals and cations) of a particular PAH throughout a diffuse or translucent cloud with a given extinction $A_{\rm V}$, electron fraction structure and geometry (slab or sphere). We select the appropriate formalisms (e.g., for reaction rates) from Le Page et al. (2001,2003), Allain et al. (1996b), Salama et al. (1996), Bakes & Tielens (1994), Allamandola et al. (1989) and Omont (1986).

The charge state distribution is calculated via the different reaction rates, which in turn are a function of PAH charge Z, extinction in the cloud $A_{\rm V}$, total extinction of the cloud $A_{\rm V{\rm total}}$, total hydrogen density $n_{\rm H}$ (i.e., related to electron density $n_{\rm e}$ through $x_{\rm e}$), gas temperature $T_{{\rm gas}}$ and specific PAH. The choice of the PAH determines the ionization potential IPZ, electron affinity EA, UV absorption cross section $\sigma_{\rm UV}$ and electron sticking coefficient s(e). We will use s(e) = 1 in the remainder of this paper, unless explicitly stated otherwise. The electron fraction is given by the electron density distribution of the modeled cloud HD 147889 (Sect. 3), and is shown in Fig. 2. We run our model for the predicted total hydrogen density and interstellar radiation field range derived for HD 147889 in Sect. 3.


  \begin{figure}
\par\includegraphics[angle=-90,width=8.8cm,clip]{0331fig2.eps}\end{figure} Figure 2: Temperature T (solid line) and electron abundance $x_{\rm e}$ (dotted line) as a function of optical depth $A_{\rm V}$ (to the center of the cloud).
Open with DEXTER

In effect, choosing a specific PAH (i.e., adopting values for IPZ, $N_{\rm C}$, s(e)) directly limits the number of free parameters. Thus, for a specific PAH, the charge state distribution depends only on the photon flux N(E), the electron density $n_{\rm e}$ and the gas temperature $T_{{\rm gas}}$. For the diffuse interstellar radiation field we use an extension of the Draine field (Draine 1978) by Allain et al. (1996a), who provide an analytical (polynomial) expression.

This approximation of the interstellar radiation field slightly overestimates the average interstellar radiation field as measured by Gondhalekar et al. (1980). It is evident that this uncertainty in the exact radiation field affects the charge state distribution calculations. However, since the UV radiation field  $I_{\rm UV}$ derived from the model for HD 147889 is expressed in units of the Draine field, with an error of 10%; this effect is negligible compared to the uncertainties in the reaction rates as discussed in Appendix A. For a general diffuse cloud we would have no a priori knowledge of the interstellar radiation field. Fortunately, for the computation of accurate absolute values for f(Z), detailed chemical and physical constraints were obtained for HD 147889 (Sect. 3). The derived diffuse UV field is then multiplied with the radiative transfer grid for this specific cloud (see Sect. 3). The radiative transfer grid varies with $A_{\rm V}$, and so does the diffuse UV field passing through the cloud, and thus also the charge state distribution.

For these calculations we will assume that the PAH distribution in the cloud does not perturb the parameters that characterize the cloud. This means that the reaction rates remain unchanged in the calculations. We also assume that the average time for destruction of the PAH is at least an order of magnitude larger than the time needed to reach a charge state equilibrium. In Appendix A we describe how we calculate the necessary reaction rates and how these are used to derive the charge state distribution. Furthermore, we will show in Appendix B the derivation of ionization potentials as a function of point group and number of carbon atoms. The point group of a molecule is defined in group theory as the characteristic set of symmetry operations for that molecule.

It is important to realize that all rates and ionization potentials are uniquely determined for each type of PAH (i.e., point group), and thus require separate calculations. The results of these calculations (Appendix A) with the different values for IP (Appendix B) are presented next.


  \begin{figure}
\par\includegraphics[angle=-90,width=8.3cm,clip]{0331fig3.eps} \end{figure} Figure 3: Total integrated charge fraction of the cation, neutral and anion charge fractions for a given PAH molecule as a function of its number of $\pi $ electrons, i.e., number of carbon atoms ($N_{\rm C}$), for the translucent cloud toward HD 147889 with $n_{\rm H}=1200$ cm-3. For clarity only the results for three point groups are displayed; solid lines for PG 4 (D$_{2\rm h}$ acene series), dashed lines for PG 3 (D2 ovalenes), and thin solid lines for PG 15 (C$_{\rm s}$ benzo-acenes) (see also Table C.1). Other point groups have intermediate values. Temperature and electron density gradients are as in Fig. 2. Figure 4 shows the effects on the charge state fractions of changes in $n_{\rm H}$ and $I_{\rm UV}$ and Fig. 5 the effects of a different cloud geometry.
Open with DEXTER


  \begin{figure}
\par\mbox{\includegraphics[angle=-90,width=8.3cm,clip]{0331fig4.e...
...{5mm}
\includegraphics[angle=-90,width=8.3cm,clip]{0331fig5.eps} }
\end{figure} Figure 4: Effect of variation in hydrogen density and $I_{\rm UV}$ on the total integrated cation, neutral and anion fractions for the HD 147889 cloud (1/2  $A_{\rm V} = 2.3$ mag) as a function of the number of carbon atoms in point group 4. Temperature and electron density gradients as in Fig. 2. The top panel shows the fractions for three different hydrogen densities ( $n_{\rm H} = 1000$, 1200 and 1400 cm-3), but with $I_{\rm UV} = 11$ kept fixed. In the bottom panel the situation is reversed. We fix $n_{\rm H}=1200$ cm-3 and vary the interstellar radiation field ( $I_{\rm UV} = 9$, 11 and 13 eV). In both panels the solid lines represent the derived model parameters, and the dashed or dotted lines the derived values plus or minus their associated uncertainties. It is evident that the errors in $I_{\rm UV}$ and $n_{\rm H}$ have a significant impact on the charge state distribution, making the distinction with other point groups less clear.
Open with DEXTER

In Figs. 3 to 5 we look at variations in the charge state distribution f(Z) due to changes in the choice of point group, the hydrogen density, the UV radiation field and the cloud geometry, respectively. Figure 3 shows the relationship between the the integrated cation, neutral and anion charge fraction and the number of carbon atoms in the PAH for three different point groups (3, 4 and 15, see Table B.1). All other PAH point groups have charge state fraction curves that lie between those shown in Fig. 3. For all point groups, the cation fraction increases with size while the neutral and anion fractions decrease. Within each point group, as the size of the PAH increases, the value of the first ionization potential decreases. As a result, larger PAHs are more easily positively charged. The differences between the different point groups are most pronounced for the cations and neutrals; the total fraction of cations can change by a factor 0.2 between different point groups. The propagations of the uncertainty of the values for the cloud parameters $n_{\rm H}$ and $I_{\rm UV}$ is shown in Fig. 4. The uncertainties (about 20%) introduced by $n_{\rm H}$ and $I_{\rm UV}$ (Sect. 3) give variations in the charge state distribution of 15% or less. An increase in hydrogen density (e.g., increased recombination rate) has a similar effect as a decrease of the UV radiation field (decreased ionization rate), and vice versa. For HD 147889 the choice of geometry (slab or sphere) affects the charge state distribution very little, i.e., less than 10% (Fig. 5). The density, radiation field and point groups have similar quantitative effects on the charge state distribution.


  \begin{figure}
\par\includegraphics[angle=-90,width=8.3cm,clip]{0331fig6.eps} \end{figure} Figure 5: Cation, neutral and anion charge state distribution of a general PAH ( $N_{\rm c} = 50$ and IP = 7 eV) as a function of extinction toward HD 147889 for both slab (dashed) and sphere (solid) geometry. The traces shift vertically almost linearly with increasing/decreasing number of carbon atoms. A change in $N_{\rm c}$ of 6 C atoms will result in a vertical shift of up to 0.1. This effects is largest for smaller PAHs since the ionization potential drops faster for small values of $N_{\rm c}$. For larger PAHs, where the decrease in ionization potential is much smaller for increasing $N_{\rm c}$, the effect of further increasing the size is less apparent.
Open with DEXTER

5 Electronic spectroscopy of selected PAHs

The charge state fractions for each point group that were obtained from the modelling of PAHs toward HD 147889 were used to calculate electronic spectra. Electronic (laboratory) spectra of gas phase PAH molecules are still not available for large PAH molecules. Matrix Isolation Spectroscopy results and quantum-mechanical calculations on PAH molecules will provide us here with the information needed for the calculation of synthetic PAH spectra.

Weisman et al. (2003) and Halasinski et al. (2003) have calculated the electronic transition symmetries and oscillator strength for a few specific neutral and monovalent PAHs using elaborate computational routines (time dependent density functional theory, TDDFT). In this study, where the focus is on the size and charge state distribution of a very large set of PAHs rather than on the energy resolution of a few specific species, we have adopted a semi-empirical calculation routine that permits a large set of data without using extensive computing capabilities. This routine, although not specifically suited for large molecular systems, has been used to obtain UV/Visual spectroscopy information for a large number of PAH molecules. The risk of missing some energy states by oversimplification inherent in this routine is not considered to be a problem since only the lowest transition states were taken into account. The choice of PAHs used here is predominantly based on the availability of data on their ionization potential.

  
5.1 Semi-empirical calculations on neutral and ionized PAHs

Input coordinates for each of the PAHs given in Fig. C.1 were geometrically optimized with a semi-empirical PM3 routine using a Polak-Ribiere conjugate gradient method. Geometries of neutral PAHs and monovalent PAH ions were optimized separately. Following geometry optimization, electronic spectra were calculated with the semi-empirical ZINDO (Zerner's Intermediate Neglect of Differential Overlap) routine. This electronic structure routine has been parameterized for spectroscopic properties of molecules. Calculations were performed with the Restricted Hartree-Fock method with orbital overlap weighting factors of 1.267 for $\sigma{-}\sigma$ bonds and 0.510 for $\pi{-}\pi$ bonds and singly excited configuration interaction with 10 occupied and 10 unoccupied molecular orbitals. We have obtained the two lowest-energy allowed transitions for the first 10 benzene moieties in 12 homologue series of PAHs. We define a homologue series as a series of PAH molecules that differ only in respect to the number of benzene moieties but that all belong to the same point group. As an example note the progression in the acene homologue series belonging to the D$_{2\rm h}$ point group: benzene (C6H6), naphthalene (C10H8), anthracene (C14H10), ... Note that more than one homologue series can belong to the same point group. As an example the acene homologue series and the condensed PAH homologue perylene series: perylene (C20H12), tribenzo[de,kl,rst]pentaphene (C30H16), benzo[1,2,3-cd:4,5,6-c'd']diperylene (C40H20),... both belong to the D$_{\rm 2h}$ point group. The PAH sample set consists of 15 different homologue series divided over 7 point groups: D$_{\rm 6h}$(1), D$_{\rm 3h}$(1), D2(1), D$_{\rm 2h}$(4), C$_{\rm 2h}$(2), C$_{\rm 2v}$(5) and C$_{\rm s}$(1), where the number in parentheses is the number of homologue series per point group. Calculated transition energies and oscillator strengths are given in Table C.1.


  \begin{figure}
\par\includegraphics[width=16.2cm,clip]{0331fig7.eps}
\end{figure} Figure 6: ZINDO-calculated wavelength as function of number of carbon atoms in the PAH molecule. For clarity each panel shows an example of the structure of a PAH in each homologue series as well as the spectral range where DIBs are found (shaded regions, $\sim $4000-10 000 Å). Squares indicate absorptions for PAH anions, circles for neutral PAHs and triangles for cation PAHs. The figure illustrates that for larger molecules (>40 C atoms) the lowest energy transitions fall outside the DIB range. Note that only for the picene homologue series (panel f)) do the absorptions of cations and anions of larger PAHs (>50 C atoms) in this series still fall inside the DIB range.
Open with DEXTER

A linear least squares fit of the ZINDO/S and laboratory data in neon Matrix Isolation (MIS) experiments (Joblin et al. 1999; Ruiterkamp et al. 2002; Salama et al. 1994; Salama & Allamandola 1991,1993; Salama et al. 1999) yields $\Delta E_{\rm calc.}/\Delta E_{\rm exp.} =1.0469$. This level of theory is at least sufficient to give accurate ( R2=0.9778) predictions of the neutral and monovalent cation transitions of PAHs that are available in the literature. We now use this routine to predict the transitions of molecules not available in the literature and will use only calculated values for the calculation of synthetic PAH spectra. Trends in oscillator strength (O.S. or f), as seen in the TDDFT calculations of the terrylene homologue series (Halasinski et al. 2003), are qualitatively reproduced by our ZINDO/S calculations(see Table 4 and Fig. 6). Quantitatively the ZINDO/S routine reproduces the TDDFT results with modest success. A linear least squares fit of the ZINDO/S and TDDFT oscillator strength data yields $f_{\rm ZINDO/S} / f_{\rm TDDFT} =0.8918$.

Note that we further expect a shift due to matrix effects when comparing matrix isolated laboratory spectra with gas-phase spectra (Salama et al. 1999). For neutral PAHs we expect a shift of up to 0.25% in fractional energy to the blue of the PAH band peak energy in the gas-phase as compared to the value in the neon matrix. For ionized gaseous PAHs, we expect shifts of up to 0.5% in fractional energy to the blue of the PAH ion band peak energy as measured in the matrix (Joblin et al. 1999).

We reiterate here that although the ZINDO routine quite accurately reproduces the matrix isolation spectroscopy data for small PAHs (see Table C.1), it is unclear whether this level of theory is sufficient to yield accurate predictions for the larger PAHs in our sample set. However, the lack of spectroscopic data of homologue series and proper caution in interpreting the results justifies the approach in this paper. In this paper we pay no attention to the band profile of the DIBs since our approach does not permit the calculation of bandwidths. However, we refer to Malloci et al. (2003) for a quantitative treatment of DIB profiles. Their results show that the predicted rotational band profiles are remarkably insensitive to ambient conditions. The authors present a model that offers a quantitative link between any given PAH and the observed DIB profiles. We have separately treated the electronic transitions obtained with the ZINDO routine and the transitions that were taken from published laboratory spectra.

Table 4: Comparison between experimental, TDDFT and ZINDO/S transition energies in the perylene, terrylene and quaterrylene homologue series.


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{aa0331f8.eps}\end{figure} Figure 7: Transition wavelength as a function of the number of carbon atoms for condensed series of PAHs. For clarity each panel shows an example of the structure of a PAH in each homologue series as well as the spectral range where DIBs are found (boxed regions). Squares indicate absorptions for PAH anions, circles for neutral PAHs and triangles for cation PAHs. For the condensed pyrene, coronene and ovalene series (panel b)) data were taken from Weisman et al. (2003) and show only cation absorptions. The figure shows that for the condensed PAHs all transitions fall inside the spectral region where DIBs are found. We extrapolated the trend that is also apparent in the condensed series to include circumcircumcircumcoronene (C150H30) to show the maximum size of PAHs within this homologue series.
Open with DEXTER

Figure 8 shows the theoretical UV/Visual spectrum of PAH molecules for six different possible size distributions in the simulated interstellar cloud environment toward HD147889. The size distributions define the relative contribution for all PAHs with a given number of carbon atoms, regardless of the point group. The size distribution of PAHs in interstellar environments is largely unknown but here we assess the effects of six different cases. The corresponding multiplication factors for each mass bin in our sample set were determined by dividing the distribution by the number of samples in each mass bin. For all distributions PAHs 1 through 71 in Fig. C.1. were used and the number of PAHs per mass bin can be found by comparing with Table C.1.

The first distribution, labelled (1/Conf) in Fig. 8 is based on the notion of the number of possible configurations for a given number of benzene moieties in the PAH molecule. As an example there are 7 different PAHs that consist of 4 benzene moieties. These numbers were determined by Klarner (1967) and for planar PAHs with a maximum of 10 benzene moieties these numbers are: 1, 1, 3, 7, 22, 82, 333, 1448, 6572 and 30 490. These numbers in a sense reflect the possible different isomers for a given number of fused carbon hexagons. The corresponding size distribution used here is simply one over the number of configurations. The second distribution labelled (1/Mw) is the inverse of the molecular mass for each PAH per distribution bin. All other cases are based on normalized Gaussian distributions labelled (Gauss x/y) with x the center peak position in number of carbon atoms and y the full width at half maximum (FWHM) in number of carbon atoms. Black bars in Fig. 8 are associated with transition wavelengths of the neutral molecules, red for negatively charged molecules (anions) and blue for positively charged PAHs (cations).

Note from Fig. 8 that the neutral absorptions of our PAHs between 3000 and 4000 Å dominate the spectrum. It is clear that for a given column density of PAHs we expect a major contribution of the absorptions between 3000 and 4000 Å to the UV spectrum of HD 147889. Only a hand full of DIBs are known below 4430 Å but we note that the spectral region between 3000 and 4000 Å in our synthetic spectrum shows so many (neutral) PAH transitions in this spectral region that we can sum discrete contributions to a continuum band. By comparing the integrated optical depth of the UV extinction curve for HD 147889 (kindly provided by K. Gordon) in the 3000-4000 Å region and the integrated optical depth of our PAH transitions we could determine an upper limit for the contribution of our PAHs to the extinction curve. We expressed the PAH spectrum in terms of the optical depth ( $\tau_{\lambda}$) through:

 \begin{displaymath}
\tau_{\lambda} = \frac{\pi \ e^2}{m_{\rm e} \ c} \ f \ \lambda \ N
\end{displaymath} (1)

where (f) is the oscillator strength (see Table C.1) of a transition with center wavelength $\lambda $ (Å) and N (cm-2) is the column density of the molecule that shows this transition (Spitzer 1978a). The constant $\frac{\pi \ e^2}{m_{\rm e} \ c} = 2.654\times10^{-2}$ cm2 s-1. We have assumed a peak profile that is simply a line at peak center ($\lambda $ Å) and height f. When we assume that the total extinction for HD 147889 in this spectral region (3000-4000 Å) is completely due to these PAHs we find an upper limit on the column density of the PAH molecules in our data set per distribution.

We can now draw general conclusions (see Sect. 6), on the effects of the size distribution on the overall appearance of the simulated spectrum compared to the measured DIB spectra toward HD 147889. All distributions in Fig. 8 were determined from the maximum column density allowed when compared to the extinction curve in the 3000-4000 Å range. The effects of the size (1/Configurations) distribution show (Fig. 8) that the steepest declining distribution underestimates the contribution of larger PAHs. The DIB range is dominated by electronic absorptions of only the smallest PAHs such as the anions of benzene, naphthalene, phenanthrene and anthracene. These small PAHs are, however, expected to be less stable under the conditions in the diffuse interstellar medium. Hence, we do not consider this a realistic situation for the PAH size distribution toward HD 147889. The results for the 1/mass distribution show a similar result (Fig. 8), but with slightly higher contributions from larger PAHs. However, this distribution still predicts that large PAHs are not contributing and we discard it as a viable candidate distribution too. Of all the Gaussian size distributions the one centered around 50 carbon atoms and FWHM of 25 carbon atoms (Fig. 8) reproduces the observed DIB spectrum, (see the top panel of Fig. 9) best although the agreement is not convincing to the eye.

Table 5: Contributions of the ionization states to the total EW of the synthetic PAH spectrum.

We used the obtained column density to calculate a synthetic PAH spectrum for the Gauss(50/25) size distribution, shown in Fig. 9. For the size distribution with a Gaussian profile centered at 50 C atoms and with a FWHM of 25 C atoms we calculated a column density of $2.4\times10^{14}$ cm-2 (see Eq. (1); f values have been taken from Table C.1). For this column density these PAHs then contribute half of the observed EW ($\geq$3.8 Å per unit reddening) for DIBs toward HD 147889. Thus, the UV absorption limits the column density of catacondensed PAHs to less than $2.4\times10^{14}$ cm-2 and such low column densities result in less than 50% contribution to the visible DIBs. Hence, catacondensed PAHs could not constitute the bulk of the DIBs. On the other hand they could still be a major contributor to the DIBs observed toward HD 147889.

Table 6: Equivalent widths (mÅ) for the strongest DIBs in HD 147889. Where appropriate the spectrum was corrected for telluric lines with a standard star spectrum (see also Sect. 2.1). Crosses (x) indicate that part of the spectrum was not covered or not of sufficient quality (e.g., fringing or telluric lines) to derive values or upper limits for DIBs. Dashes (-) indicate that the DIB was not detected, because it was too weak or the S/N was not high enough. The 5780/5797 ratio is $\sim $2.3.

Table 5 shows the contribution of the different charge states for the chosen size distributions to the total EW of the PAH spectrum. These values show that even though the total cation fraction can be close to 50% (see Fig. 3) the contribution of these cations to the total EW is only 20%. It is important to keep in mind that our PAH sample set does not include any molecules with more than 42 carbon atoms. However, even for the Gaussian size distribution centered at 100 carbon atoms there is a non-zero contribution of small PAHs (<42 C atoms) to the total EW.

5.2 Matrix Isolation Spectroscopy and TDDFT results

To fully exploit the presented model we have entered published electronic spectra of large PAHs into our model for HD 147889 (Joblin et al. 1999; Halasinski et al. 2000; Ehrenfreund et al. 1992; Ruiterkamp et al. 2002; Salama et al. 1994; Salama & Allamandola 1991,1993; Salama et al. 1999). These laboratory spectra were taken at low temperatures ($\sim $4 K) in inert neon matrices and the spectra closely resemble the gas phase values of the electronic transitions. Recently published TDDFT calculations (Weisman et al. 2003; Halasinski et al. 2003) were also entered into our routine. These homologue series comprise the coronene (C24H12 [72]), circumcoronene (C54H18 [74]) and circumcircumcoronene (C96H24 [75]) series; the ovalene (C32H14 [73]), and circumovalene (C66H20 [76]) series; the pyrene (C16H10 [22]), circumpyrene (C42H16 [77]) and circumcircumpyrene (C80H22 [78]) condensed series of PAHs; and the perylene (C20H12 [14]), terrylene (C30H16 [15]) and quaterrylene (C40H20 [16]) condensed series. Individual entries of PAHs are naphthalene (C10H8 [5]), anthracene (C14H10 [6]), phenanthrene (C14H10 [36]), pentacene (C22H14 [8]), decacyclene (C36H18 [79]), Dibenzo[a, jk]phenanthro[8,9,10,12-cdefgh]pyranthrene (C44H20 [80]) and dicoronylene (C48H20 [81]). For each entry the corresponding number in Fig. C.1 is given in brackets. Figure 10 shows the results for these PAHs in our model for HD 147889 with a Gaussian size distribution with center at 50 C atoms and a width of 50 C atoms. The treatment of these data is identical to the treatment used for the ZINDO results except that we scaled the total EW of these PAHs to match the total EW of the DIB range (4000-8000 Å) toward HD 147889. Since we have no information on the oscillator strengths of neutral pericondensed PAH transitions we cannot scale the column density for large PAHs to the extinction curve (as discussed in Sect. 5.1). Figure 10 shows a comparison of PAH transitions at two different scales while the top panel shows the equivalent widths of diffuse interstellar bands measured toward HD 147889. Since we used published matrix isolation spectroscopy data and quantum-mechanical computations we cannot compare the peak positions in a straightforward way.


  \begin{figure}
\par\includegraphics[width=15.8cm,clip]{0331fig9.eps} \end{figure} Figure 8: Comparison of the effects of different size distributions on the calculated UV/Visual spectrum of the 71 PAHs listed in Table C.1, for all transitions between 2000 and 14 000 Å for HD 147889. Black bars indicate electronic transitions of the neutral molecule, red bars indicate electronic transitions of anions and blue bars indicate electronic transitions of cations calculated with the semi-empirical ZINDO/S routine. The panels show different size distributions: 1/Conf is based on the notion of different configurations of PAHs for a given number of carbon atoms; 1/Mw is a size distribution based on molecular mass; Gauss x/y indicates a Gaussian size distribution with center position x and width y in number of carbon atoms. We conclude from this selection of size distributions that the Gauss 50/25 distribution is the most likely one. Details are explained in Sect. 5.1.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=16.3cm,clip]{0331fig10.eps}
\end{figure} Figure 9: Comparison of the ZINDO-calculated electronic transitions of small catacondensed PAHs (numbers [1-71] in Fig. C.1.) and a compilation of DIB data. The top panel shows the equivalent widths of diffuse interstellar bands measured toward HD 147889 with UVES data, FEROS data and data taken from Thorburn et al. (2003). Bottom panels show the resulting electronic transitions of PAHs based on semi-empirical calculations where the lower panel shows the figure on a smaller vertical scale. The data were scaled by a Gaussian size distribution with center position at 50 carbon atoms and a width at half maximum of 25 carbon atoms. Shaded areas indicate the region where the DIB spectrum overlaps with telluric absorption lines; reliable data are mostly unavailable in this region. Note that our PAH spectrum predicts absorptions in the telluric obscured spectral region, which can at this moment not be checked observationally. The 3$\sigma $ detection limits (20 mÅ) are given for narrow DIBs in clear regions, however they are much higher for very broad DIBs (such as those expected below 4500 Å) which are more difficult to detect, and for DIBs in telluric obscured regions, or in ranges where the CCD has lower sensitivity above 9000 Å.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=16.2cm,clip]{0331fig11.eps}
\end{figure} Figure 10: Comparison of electronic transitions obtained from laboratory measurements or TDDFT quantum-mechanical calulations for large pericondensed PAHs (numbers [14-16] and [72-81] in Fig. C.1.) and a compilation of DIB data. The top panel shows the equivalent widths of diffuse interstellar bands measured toward HD 147889 with UVES data, FEROS data and data taken from Thorburn et al. (2003). Bottom panels show a composite spectrum of PAHs that have been measured in the laboratory where the lower panel shows the figure at a smaller scale (Halasinski et al. 2000; Ehrenfreund et al. 1992; Ruiterkamp et al. 2002; Salama et al. 1994; Salama & Allamandola 1991) or calculated by TDDFT quantum-mechanical methods (Weisman et al. 2003; Halasinski et al. 2003). The 3$\sigma $ detection limits (20 mÅ) are given for narrow DIBs in clear regions, however they are much higher for very broad DIBs (such as those expected below 4500 Å) which are more difficult to detect, and for DIBs in telluric obscured regions, or in ranges above 9000 Å where the CCD has lower sensitivity.
Open with DEXTER

  
6 Results: Properties of the synthetic PAH spectra

As presented in Figs. 9 and 10 our model reproduces the general features of the DIB spectrum for HD 147889 only poorly. For the semi-empirical calculations on PAHs depicted in Fig. 8 we can identify three regions. The short-wavelength part (2000-4000 Å) is dominated by neutral PAH transitions that strongly contribute to the total absorption for most of the considered size distributions. The next wavelength range (4000-10 000 Å) shows minor contributions from neutral PAH molecules but predominantly cation transitions with some strong contributions from PAH anions. The last wavelength range longwards of 10 000 Å shows predominantly anion contributions with a large number of cation contributions. Observations in this spectroscopic region are severely hindered by absorption lines due to the Earth's atmosphere. Although currently no space based observatory is available for this spectral region we expect that in the future a large number of new DIBs will be found longwards of 10 000 Å with high resolution spectrographs in the near-IR.

The strongest absorptions in the calculated DIB spectrum (Fig. 9) are from the largest PAHs in the D$_{\rm 2h}$, C$_{\rm 2v}$ and C$_{\rm s}$ point groups. That there is such a poor agreement between the observed DIB spectrum and our PAH spectrum does not mean that small PAHs with up to 40 carbon atoms cannot contribute at all to the DIBs. In fact we probe only a small portion of the possible PAH molecules with up to 40 carbon atoms. We see that specific PAH point groups and PAH sizes could give rise to strong DIBs. It is the identification of these special PAH series or species that guides our efforts.

Although we cannot assign any DIBs to specific PAHs we can compare the relative contributions for strong DIBs ( ${\it EW} > 100$ mÅ), intermediate DIBs ( $50 < {\it EW } < 100$ mÅ) and weak DIBs ( ${\it EW} < 50$ mÅ). For HD 147889 we find a relative ratio for the number of strong/intermediate/weak DIBs of 1/1.5/3.3. Note that the detected DIBs only provide a lower limit. As much as twice the number of small DIBs is predicted from comparing with other sources that have been observed at higher signal to noise ratios. For the semi-empirical calculated DIB spectrum we find a relative ratio of 1/1.5/70. We see that the DIB spectrum has a much smaller contribution of weak features than the predicted PAH spectrum. The main contribution from our PAH sample set could thus be predominantly to the weak DIBs.

The synthetic PAH spectrum for the laboratory data and TDDFT calculations show a similar picture. As noted by Weisman et al. (2003) there is also a clear indication for progression in the wavelength positions of the homologue condensed PAHs (as determined by TDDFT calculations and shown in Fig. 7, panel b). When we extrapolate this trend to the highest wavelengths where DIBs are found ($\sim $10 000 Å), we find that the largest PAH in this homologue series that still could contribute to the DIB spectrum is circumcircumcircumcoronene (C150H30). For the catacondensed homologue series, shown in panel (a) of Fig. 7 and calculated with the ZINDO routine, the larger PAHs will shift the absorption features further toward the infrared part of the spectrum and outside of the DIB range (Fig. 6). Note that only for the picene homologue series (panel f) the absorptions of cations and anions of larger PAHs (>50 C atoms) in this series still fall inside the DIB range. For the condensed homologue series of PAHs this shift is less prominent, such that even for large PAHs the main absorption lies within the known DIB range.

Note, however, that Fig. 4 indicates that the abundance of cations is larger than 50% in HD 147889 if the size of the PAHs is larger than 50 carbon atoms. This indicates that the cations in the condensed homologue series may form a dominant part of the DIB spectrum, provided that our choice of electron recombination rate is correct (see Appendix A). The behavior of the oscillator strength with the progression of the homologue series, in combination with the abundance (and thus the size distribution) now determines the contribution. The large condensed PAHs have only small contributions to the total EW while the perylene homologue series dominates the spectrum. Due to their high f-value the perylene homologue series can still make a large contribution to the total EW even for small abundances. No clear distinct regions can be assigned, and contributions from neutral, anion and cation PAHs are found throughout the DIB range (this may be due to the lack of laboratory data for homologue series).

When we take into consideration that more DIBs may exist at longer wavelengths than currently known, we still see that the observed spectra show only minor contributions in the 8000-10 000 Å range. This spectral range is accessible for Earth-based telescopes (although this region suffers from contamination by atmospheric water lines) and our model should reproduce this observation. Also for the condensed PAH series these Gaussian size distributions reproduce the possible contributions to the observed DIB spectrum best. The cations of the pericondensed homologue series (see Figs. 7 and 10) do not lead to any significant DIB features above 10 000 Å and as such are quite distinct from the catacondensed PAH species listed in Table C.1.

  
7 Discussion

The translucent cloud toward HD 147889 is a special case that is ideal for the study of PAH charge state distributions. Since the line of sight toward this star probes a single cloud it is especially suited for modelling. Additionally with an E(B-V)= 1.09 mag and $R_{\rm V} = 4.22$ a substantial number of DIBs can be detected in this line of sight. HD 147889 has a relatively large reddening, but still shows a single cloud structure. Most lines-of-sight with significant reddening (>0.5) show multiple cloud structure that is much more difficult to model.

We find that the charge state of specific PAHs at a certain extinction in the cloud depends on their ionization potential which in turn depends on the molecular size (number of $\pi $ electrons) and the molecular structure. The charge state due to the UV field at a certain extinction in the cloud has been determined as a function of both the molecular size and structure. By integration over the line of sight toward HD 147889 we determined the total charge fractions and calculated a synthetic PAH spectrum for this object. We find that different wavelength regions in the DIB spectrum are populated preferentially by different PAH charge states. This also strongly depends on the underlying PAH size distribution and molecular symmetry. For the catacondensed PAH sample set (with 10 to 40 carbon atoms) we could constrain the distribution for HD 147889 by the observed DIB spectrum and found a Gaussian distribution with a mean of 50 carbon atoms and a FWHM of 25 carbon atoms. From the strong neutral PAH absorptions between 3000 and 4000 Å we could constrain the total small catacondensed PAH abundance along the line of sight toward HD 147889 and found an upper limit of $2.4\times 10^{-8}$ compared to the column density of hydrogen. Hence the contribution of these small catacondensed PAHs can never exceed half of the total EW observed in the DIB range along the line of sight toward HD 147889.

If the DIBs are traced by the catacondensed and mixed cata/pericondensed PAH homologue series (with 10 < C atoms < 50) then we predict that more DIB features should be present longwards of 10 000 Å. Note that for PAHs smaller than 40 C we expect a large range of partial dehydrogenation (Le Page et al. 2001; Vuong & Foing 2000). Therefore these PAHs will appear in a large number of dehydrogenated states and isomers, and their spectroscopic signature will be spread over the spectrum in a number of small bands, probably below the detection threshold. If we assume that even larger PAHs that belong to these homologue series would dominate the distribution, we see that their contributions would be well longwards of 10 000 Å. Figure 6 clearly shows this observation for the catacondensed PAH series (panels a through j). Therefore, if the DIBs are caused by the catacondensed PAHs we investigated, we would conclude that only PAHs up to 40 carbon atoms contribute. Note from panel f that this is not the case for the benzo(c)picene homologue series. Larger catacondensed PAH cations and anions would have their strongest transitions in the infrared part of the spectrum. For the mixed cata/pericondensed PAHs (panels k and l in Fig. 6) we find that if the pyrene unit is centered in the molecule, the ion absorptions for molecules up to 50 carbon atoms still fall inside the DIB range (panel k). For the pericondensed PAHs (>40 C atoms) we find that the spectrum is dominated by the specific PAHs (perylene homologue series) while the bulk of the PAHs in the sample set has contributions in the correct EW range. A mix of PAHs could reproduce the small DIBs (<50 mÅ) in the spectrum but not the strong features. We summarize that only the weaker DIBs could be reproduced by a mix of small catacondensed and peri/catacondensed PAHs (<50 C atoms) or large pericondensed PAHs (50 < C atoms < 100), but the strong DIBs can only be reproduced by adding very specific PAH molecules or homologue series to the sample set. Though the existence of small PAHs cannot be excluded (Hudgins & Allamandola 1999) their contribution to the DIB spectrum is less than 50% toward HD 147889.

When we compare the band positions with the PAH spectra to the DIB spectrum toward HD 147889 we observe that there is no direct alignment to any of the strong DIBs. Even when we consider that semi-empirical calculation techniques are know to over-estimate the energies for electronic transitions we find that the only strong DIB that could align with a PAH transition is found at 5780 Å. In this case the DIB at 5780 Å could align with the strongest neutral transition of quaterrylene (benzo[1,2,3-cd:4,5,6-c'd']diperylene, number 16 in Table C.1) at 5608 Å. This molecule has been measured in the laboratory in a neon matrix (Ruiterkamp et al. 2002) and shows the strongest neutral, cation and anion transitions at 607.9, 833.7 and 881.8 nm, respectively. Toward HD 147889 we find no DIBs that match these wavelength positions or relative intensities for the neutral, cation and anion bands. This agrees with our previous conclusion (Ruiterkamp et al. 2002) that it is unlikely that quaterrylene contributes to the DIBs. It is clear from the expected shifts between calculated and measured absorptions that no identification of the small DIBs is possible. We note however that from Fig. 6 it is easy to determine which transitions of what PAH fall inside the DIB range.

In principle, the procedure used in this work can be applied to any line of sight, provided that the physical model for the cloud is reliable. This procedure can also be used to generate DIB spectra for alternative molecular carriers such as carbon chains. In order to draw geometry-independent conclusions, a large number of lines of sight should be modeled, provided that their impact parameter distribution can be determined. One would then have to resort to mean values of the density, electron abundance, and temperature in the physical model which would limit the ability to identify individual features. Nevertheless, this method would lead to statistically more significant conclusions as far as the global characteristics of the DIB spectrum are concerned. High resolution spectra for the spectral regions that are currently unaccessible because of telluric contamination are of crucial importance to fully exploit the predictive power of our procedure since many DIBs may be present at wavelengths longwards 10 000 Å.

Laboratory efforts should focus on homologue series to better constrain the quantum-mechanical calculation efforts. As long as gas phase spectroscopic data remain unavailable, Matrix Isolation Spectroscopy experiments are a valuable alternative. Our results also indicate that the main focus for laboratory spectroscopy should be on larger condensed homologue PAH series that represent better candidates for the DIBs. The discrepancy between the increasing oscillator strength with increasing size in the perylene homologue series (Halasinski et al. 2003), which dominates the PAH spectrum in both Figs. 9 and  10, and the decreasing oscillator strength in the pericondensed PAH series (Weisman et al. 2003) show the importance of studying these systems. High values for oscillator strengths stand out in the simulated PAH spectrum, and since strong DIBs are observed, we predict that specific PAHs or homologue series of PAHs with high oscillator strengths may be the underlying carriers for these strong absorptions. Molecules with a small oscillator strength contribute only to the DIB spectrum if their column density is large.

Apart from laboratory reference data there is an equal need for theoretical quantum-mechanics data of PAH homologue series for reliable oscillator strengths. The pipeline nature of our procedure allows us to easily predict the spectroscopic response of a specific PAH or group of PAHs once these data become available in the literature. To fully exploit the pipeline nature of our approach we plan to expand our PAH transition database. Once we can use many more PAH transitions it becomes useful to look for frequency matching between PAH transitions and DIBs. The current procedure does not yet incorporate mediating processes on the incorporated molecules such as dehydrogenation and destruction rates. This will be implemented in future work.

Acknowledgements
This research was conducted under SRON project MG-049, NWO-VI and ESO programs 67.C-0281 and 64.H-0224. The authors thank ESO for using the VLT, as well as K. Gordon for providing the extinction curve of HD 147889. Furthermore, we thank the referee for his critical comments. This research has made use of: the SIMBAD database, operated at CDS, Strasbourg, France; and the NIST Chemistry WebBook, NIST Standard Reference Database Number 69, March 2003, ed. P. J. Linstrom and W. G. Mallard. The authors thank the referee, John Mathis, for his insightfull and constructive comments.

References

 

  Online Material

   Appendix A: The PAH charge state distribution

Adopting a multiple charge model involves solving Eq. (A.8) from Bakes & Tielens (1994) which gives the probability of finding a grain (e.g., PAH) at charge Ze. Since we are only interested in PAH anions, neutrals and cations, i.e., a 3-charge state distribution, we have:

$\displaystyle \sum_{Z = -1}^{Z = +1} f(Z) = 1.$     (A.1)

Following the argument of Bakes & Tielens (1994) the contribution of the accretion rate of ions is neglected, i.e., it is much larger than the photodetachment rate. This then gives
  
$\displaystyle f(1) = f(0) \frac{\ensuremath{{k}_{\rm ion}} }{\ensuremath{{k}_{\...
...nd} \ f(-1) = f(0) \frac{\ensuremath{{k}_{\rm e}} }{\ensuremath{{k}_{\rm ph}} }$     (A.2)

and we directly obtain f(0) from

 \begin{displaymath}f(0) = \left(1 + \frac{\ensuremath{{k}_{\rm ion}} }{\ensurema...
...uremath{{k}_{\rm e}} }{\ensuremath{{k}_{\rm ph}} }\right)^{-1}
\end{displaymath} (A.3)

where the four reaction rate coefficients are (in units of s-1): If we know these reaction rates, it is then possible to calculate, from Eqs. (A.2) and (A.3), the charge state probability distribution. And since the reaction rates depend on $A_{\rm V}$, $A_{\rm V{\rm total}}$ (i.e., via N(E, $A_{\rm V}$)), $n_{\rm e}$, $T_{{\rm gas}}$ and PAH, so does f(Z).

A.1 Rate coefficients of basic processes

In this paper we will adopt from the literature the following equations for the four reaction rates. We also discuss briefly our motivation for each choice.

photoionization rate (PAH + $h\nu \rightarrow $ PAH $^+ + {\rm e}$).
For the photoionization rate ( \ensuremath{{k}_{\rm ion}}) we adopt (Le Page et al. 2001; Allain et al. 1996b):
 
$\displaystyle \ensuremath{{k}_{\rm ion}} $ = $\displaystyle \int_{\rm IP}^{13.6 ~{\rm eV}} Y_{\rm ion}(E) \ \sigma_{\rm UV}(E) \ {N}_{\rm C} \ N(E) \ {\rm d}E \ \ ({\rm s}^{-1})$ (A.4)

with
$\displaystyle Y_{\rm ion}(E) = 0.8 {\rm exp} \left[-0.00128( (14.89 - {{\rm IP}}_{\rm coronene})/ (14.89 - {{\rm IP}}_{\rm PAH}) (E - 14.89))^4 \right]$     (A.5)

the ionization yield (Le Page et al. 2001; parametrization of Verstraete et al. 1990), where the IP and E are in units of eV, and IP $_{\rm coronene} = 7.29$ eV.

It is evident that laboratory measurements of the UV absorption cross section and the ionization yield (and also electron attachment) would be preferred over a general curve. Unfortunately at this point only limited data are available for individual PAHs (e.g., Jochims et al. 1996). Therefore, until more data are available we will adopt, for consistency, $\sigma_{\rm UV}$ and $Y_{\rm ion}$ curves that are based on studies of PAH mixtures (Joblin et al. 1992b).

The ionization potential is a strong function of the number of $\pi $ electrons and the molecular symmetry (Gallegos 1968). For the ionization potentials we use the values derived in Appendix B. The UV absorption cross section $\sigma_{\rm UV}(E)$ is given per number of carbon atoms, and we adopt the general curve from Joblin et al. (1992a).

We agree with the analysis of Le Page et al. (2003) that it is justified to use the same UV absorption cross section for both neutrals and cations. Also, instead of considering only the photons in the energy band 6-13.6 eV, Le Page et al. (2001) incorporate photons with energies $\leq$5 eV (i.e., they use the extended Draine field), which we adopt here.

Electron recombination rate (PAH $^+ + {\rm e} \rightarrow $ products).

For the electron recombination rate ( \ensuremath{{k}_{\rm rec}}) we use the following equation from Le Page et al. (2003) and Salama et al. (1996):

                                $\displaystyle \ensuremath{{k}_{\rm rec}} $ = $\displaystyle 3 \times 10^{-7} \left(\frac{300}{T}\right)^{1/2} \ s(e) \ {n}_{\rm e} \ \ \ ({\rm s}^{-1})$ (A.6)
  = $\displaystyle {\rm C} \cdot {N}_{\rm C}^{1/2} \ T^{-1/2} \ s(e) \ {n}_{\rm e} \ \ \ ({\rm s}^{-1}) ,$ (A.7)

where the constant C is derived from theoretical calculations or laboratory measurements. For example, Le Page et al. (2003) use $3.5 \times 10^{-7}$ for the recombination rate of naphthalene at 300 K ( $3 \times 10^{-7}$ s-1) or, instead, that of benzene ( $1 \times 10^{-6}$, at 300 K), whereas other authors use still different values: $2.7 \times 10^{-5}$ (Omont 1986), $7.1 \times 10^{-6}$ (Allain et al. 1996b) and $2.4 \times
10^{-5}$ (Vuong & Foing 2000).

We also adopt Le Page et al.'s conclusion that for small PAH+ there is no significant evidence of size dependence. It should be noted that Vuong & Foing (2000) present an alternative expression for the electron recombination rate, based on Spitzer (1978b) and Verstraete et al. (1990), that is applicable to larger, $N_{\rm C} > 30{-}50$, PAHs. They point out that the Le Page expression is more applicable to smaller PAHs. Use of the Vuong & Foing (2000) expression in our model for HD 147889 leads to a PAH charge state distribution with almost no cations at the edge of the model cloud. Since we find this result unlikely, we have chosen to use the Le Page expression for the electron recombination rate.

The sticking coefficient s(e) of a PAH greatly depends on its electron affinity EA (Allamandola et al. 1989). For larger EA (i.e., pericondensed PAHs with more than 20 C atoms) the sticking coefficient s(e) approaches unity (Salama et al. 1996; Allamandola et al. 1989). Allamandola et al. (1989) show s(e) as a function of $N_{\rm C}$. For smaller compact PAHs, s(e) can be smaller than ${\sim} 10^{-3} {-} 10^{-6}$. For coronene and pyrene, which have low electron affinities $\approx$ 0.6 eV, Salama et al. (1996) adopt a low sticking coefficient $s(e) \approx 10^{-3}$ (Allamandola et al. 1989). Pentacene, on the other hand, has a higher electron affinity (1.1 eV) and thus a high sticking coefficient (s(e) = 1). We use the known electron affinity EA to estimate the sticking coefficient s(e). For EA $\geq$ 1 eV we will use s(e) = 1 and for EA $\leq 0.1$we will use $s(e) \sim 10^{-3}$. When applying the fit results (see Appendix B) we will adopt EA = 2.0 eV and s(e) = 1.

Electron attachment rate (PAH + ${\rm e}^- \rightarrow $ PAH-).

The electron attachment rate ( \ensuremath{{k}_{\rm e}}) produces the PAH anions. Following Allamandola et al. (1989) and Omont (1986) we also assume that the electron attaches to the PAH at the Langevin rate.

\begin{displaymath}k_{L} = 2 \pi e (\alpha / {m}_{\rm e})^{1/2} s(e) \ (\rm {cm}^3 \ \rm {s}^{-1}),
\end{displaymath} (A.8)

where $\alpha$ is the polarizability of the neutral PAH, and eand $m_{\rm e}$ the electron charge and mass, respectively. There are different methods to derive the polarizability $\alpha$. Omont (1986) gives $\alpha = 0.9 ({N}_{\rm C})^{3/2}$ Å3 for pericondensed (circular) PAHs and $\alpha =
0.1 ({N}_{\rm C})^{3}$ Å3 for catacondensed (n-cenes) PAHs. According to Allamandola et al. (1989) a better estimate would be $\alpha = 1.5 \times 10^{-24} \ {N}_{\rm C}$ which gives $k_{L} = 1.2 \times 10^{-7} {N}_{\rm C}^{1/2}$ s(e) ( ${\rm cm}^3~ {\rm s}^{-1}$).

We adopt the Omont approximation of the electron attachment rate \ensuremath{{k}_{\rm e}}:

                            $\displaystyle \ensuremath{{k}_{\rm e}} $ = $\displaystyle k_{L} \ {n}_{\rm e} \ ({\rm s}^{-1})$ (A.9)
  = $\displaystyle 9.8 \times 10^{-8} \ {N}_{\rm C}^{3/4} \ s(e) \ {n}_{\rm e} \ ({\rm s}^{-1})$ (A.10)

Electron photodetachment rate (PAH $^- + h\nu \rightarrow $ PAH + e).
The electron photodetachment rate ( \ensuremath{{k}_{\rm ph}}) is due to the absorption of visual or UV photons. Experimental values for $\sigma_{\rm ion}$ are not accurately known, however; it is common practice to assume that $\sigma_{\rm neutral} = \sigma_{\rm ion} \sim 2
\times 10^{-18}$ cm-2 per C atom. \ensuremath{{k}_{\rm ph}} is given as

\begin{displaymath}\ensuremath{{k}_{\rm ph}} = \sigma_{\rm ion} \ {N}_{\rm C} \ {N}_{\rm ph} \ ({\rm s}^{-1}),
\end{displaymath} (A.11)

with $N_{\rm ph}$ the photon flux in the Draine field, between the photodetachment threshold (i.e., the electron affinity EA of the neutral) and 13.6 eV. For EA = 2.3 eV, $N_{\rm ph}
\approx 1{-}10 \times 10^8$ photons eV-1 cm-2 s-1.

  
Appendix B: Ionization potential of PAHs

The electronic spectrum of PAH radical ions shows distinct features compared to the neutral spectrum. These features can be attributed to one-electron excitations involving the occupancy of the singly occupied molecular orbital (Halasinski et al. 2003). It has been proposed that these transitions are responsible for the diffuse interstellar absorption bands (DIBs) (Bakes & Tielens 1994; Salama et al. 1999; Allamandola et al. 1989).

  \begin{figure}
\par\includegraphics[width=12.5cm,clip]{0331fig12.eps}
\end{figure} Figure B.1: Fit of Eqs. (B.2) or (B.3) to literature data for the ionization potential of specific PAH point groups. Values for fit parameters and $\chi ^2$ are given in Table B.1.

The ionization state of interstellar PAHs determines their contribution to the DIB spectrum.

Table B.1: Results for two-parameter fits of the literature data on PAH ionization potentials to Eqs. (B.2) and (B.3).

The classical expression for the ionization potential of a PAH is

$\displaystyle {\rm IP}_Z = w + \left(Z + \frac{1}{2}\right) \frac{e^2}{C} ~{\rm eV},$     (B.1)

with w the work function in eV, Z the charge number, e in esu and C the capacitance in cm. In terms of the number of carbon atoms $N_{\rm C}$ for flat ( $\alpha = 1/2$) and spherical ( $\alpha = 1/3$) PAHs Bakes & Tielens (1994) give:
   
$\displaystyle {\rm IP}_Z = w + \left(Z + \frac{1}{2}\right) \frac{A}{{N}_{\rm C}^\alpha}~ {\rm eV}$     (B.2)
$\displaystyle {\rm IP}_Z = w + \left(Z + \frac{1}{2}\right) \frac{A}{{N}_{\rm C}^\frac{1}{3}} ~{\rm eV}$     (B.3)

were Bakes & Tielens (1994) adopt for the work function wthe value of graphite (4.4 eV Smith 1961). For the constant A they derive 25.1 and 11.1 for the planar and spherical case, respectively. This is appropriate for an infinite graphite sheet but not necessarily for PAH molecules. Additionally, for a material with finite dielectric constant $\epsilon$ the second term has to be multiplied by $\epsilon-1$/$\epsilon$ (Brus 1983). Table C.2 lists available data for the IP and EA of PAHs arranged according to the molecular symmetry. The value of the ionization energy depends strongly on the molecular symmetry and the number of $\pi $-electrons in the PAH molecule (Gallegos 1968). For a large number of PAHs these values are known; for larger PAHs with $40 \leq N_{\rm
C}\leq 120$ no data are currently available.

We therefore fitted the available data (Table C.2) for the IP of PAHs for different point groups with Eq. (B.2), where we assume that the work function w of graphite is not necessarily the proper value for an infinitely large molecule within its point group see Fig. B.1. Since the second term depends on geometry we will fit this value without reference to the underlying geometry of the capacitor, or the dielectric constant. We thus assume that all the PAHs that we consider can be seen as either flat or spherical capacitors with undetermined geometry. This does not mean, however, that the PAH molecules themselves are spherical - in fact all are flat un-curved molecules - we just indicate what geometry fits best the available data for the ionization potential. In Table B.1 we summarize the values of parameters w and A per point group. By assessing the $\chi ^2$ value for the fit of experimental IP's to the predicted IP's from Eqs. (B.2) and (B.3) we determine the best fit per point group, and subsequently determine the most appropriate geometrical representation (smallest $\chi ^2$).

The fitted values as shown in Table B.1 are used in the ionization rate (Eq. (A.4)). They can also be used to obtain the theoretical IP for PAHs with 50-100 carbon atoms. These larger PAHs are more likely to survive the harsh conditions of interstellar space. The strong absorption features in the cation and anion electronic spectra shift towards longer wavelength for larger molecules (Clar & Schmidt 1975). For each homologue series of PAHs (such as the homologue series perylene, terrylene, quaterrylene, ...) the main ionic absorption feature progressively shifts towards the red. This effect causes the feature to fall outside the spectral range where DIB features are detected. There exists an upper limit for the size of a PAH per homologue series (and thus point group) that contributes to the DIBs. However, insufficient laboratory data are available to predict this limit for each series. We also note that for larger PAHs the oscillator strength increases (Halasinski et al. 2003). The fitted IP curves now allow us to predict the IP of large PAHs in specific homologue series or point groups. Note that Salama et al. (1996) from a survey of experimental data find an upper limit of 250 carbon atoms to the molecular size of condensed PAH ions that may contribute to the DIBs.

Appendix C: Tables and PAH geometrical structures

Table C.1: ZINDO-calculated transition wavelengths $\lambda $(Å) and oscillator strength f. The numbers in Col. 1 refer to respective geometrical structures in Fig. C.1. IP values for these PAHs are listed in Table C.2.

Table C.2: First ionization potential (I.P.) for selected PAHs.


  \begin{figure}
\par\includegraphics[width=12.5cm,clip]{0331fig13.eps}
\end{figure} Figure C.1: Geometrical structures of sample PAHs. Numbers correspond to Table entries. Note that entries 72 and 73 are the structures of coronene and ovalene. Explicit hydrogens are not shown.



Copyright ESO 2005