Free Access
Issue
A&A
Volume 538, February 2012
Article Number A20
Number of page(s) 20
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201116919
Published online 27 January 2012

© ESO, 2012

1. Introduction

Protoplanetary discs set the initial conditions for planet formation, and the processing of chemical species and dust grains as discs evolve may ultimately provide the seed for life. Indeed, conditions within circumstellar discs are conducive to a rich chemistry, with a large variety of conditions ranging from icy mantles on grain surfaces in the cold disc midplane to ionised gas bathed in stellar ultraviolet radiation at the disc surface. In order to attempt to probe the gas in discs it is necessary to observe gas emission lines from a variety of chemical species (molecules, atomic species, ions), and couple these observations with detailed thermochemical disc modelling.

One important disc property for planet formation models is the amount of gas present as it evolves. The dust characteristics are generally better-understood than the gas, and the ISM gas/dust mass ratio of 100 is often adopted to give an estimate of the gas mass. However, gas masses derived from millimetre and sub-mm CO emission observations are typically lower than from dust observations with this assumption, sometimes by a factor  ~100 (Thi et al. 2001). This is variously attributed to gas depletion due to photoevaporation or planet formation, and CO freeze-out on grains (Zuckerman et al. 1995; Kamp & Bertoldi 2000). These low level rotational CO transitions are also strongly dependent on the disc size, and often optically thick under certain disc conditions, limiting their ability to probe the deeper disc layers (Panić & Hogerheijde 2009; Hughes et al. 2008). In order to better constrain the disc gas mass it is necessary to observe additional gas species, such as those observable in the far-infrared. The FIR range is a crucial observing window which can help to resolve some of these ambiguities. The fine structure lines of atomic oxygen and ionised carbon probe the gas in the warm disc surface layers, and as products of CO photodissociation provide a good test for this proposed explanation of CO under-abundance. In addition, the high level rotational transitions of CO arise from warm gas in the inner disc, probing different disc regions to the low level mm CO lines.

thumbnail Fig. 1

Observed [Oi]63 μm emission line with Herschel/PACS.

In addition to the discrepancy between gas mass estimates from millimetre continuum and CO observations, there is uncertainty regarding the derived disc outer radii from CO and continuum observations. The outer radius suggested by model fits to dust emission is often found to be smaller than that of the gas disc. In the case of the Herbig Ae star AB Aurigae, Piétu et al. (2005) suggested that a change in dust grain properties, leading to reduced opacity in the outer disc, could explain this discrepancy. For HD 163296, Isella et al. (2007) attributed the apparent difference in disc radii to a sharp drop in surface density, opacity or dust-to-gas ratio beyond 200 AU, while Hughes et al. (2008) proposed that the apparently conflicting dust and gas observations could be reconciled using a different treatment of the disc surface density at the outer edge, motivated by similarity solutions for the time evolution of accreting discs. This study attributed smaller derived disc radii from dust observations to detector sensitivity thresholds, arguing in favour of larger dust discs, albeit with an exponentially-tapered density at the outer edge. A steepening of the density profile in the outer disc has also been observed in the case of DM Tau, LkCa 15 and MWC 480 by Piétu et al. (2007).

The Herschel open-time key program GASPS (GAS in Protoplanetary Systems) aims to characterise the gas in discs at each stage of their evolution, by observing the far-IR lines of a wide range of discs, from young gas-rich accretion discs to gas-poor debris discs (Mathews et al. 2010). This paper represents the second in-depth study of a Herbig disc as part of this program, following after a study focusing on the young Herbig Ae star HD 169142 (Meeus et al. 2010).

Table 1

Observed line data with Herschel/PACS.

HD 163296 is an isolated Herbig Ae star with spectral type A1Ve, mass  ~2.5 M and stellar luminosity  ~38 L (this work), situated at a distance of  pc (van Leeuwen 2007). Scattered light and millimetre continuum observations indicate the presence of an inclined circumstellar disc extending out to a radius of 540 AU (Mannings & Sargent 1997; Grady et al. 2000; Wisniewski et al. 2008). In addition, there is evidence of an asymmetric outflow perpendicular to the disc, with a chain of six Herbig-Haro knots tracing its mass-loss history (Devine et al. 2000). The derived disc dust mass is in the range (5 − 17) × 10-4   M (Natta et al. 2004; Tannirkulam et al. 2008a; Mannings & Sargent 1997; Isella et al. 2007).

Recent near-infrared studies of the inner disc of HD 163296 indicate an inner dust rim feature enclosing a bright emission region extending inwards towards the star. This emission inwards of the dust rim has been attributed to an optically thin inner disc, although there is uncertainty regarding the size and composition of such a disc. Derived radii for the dust rim lie in the range (0.2–0.55) AU (Renard et al. 2010; Eisner et al. 2009; Tannirkulam et al. 2008b; Benisty et al. 2010; Monnier et al. 2006). An increase in opacity at the inner rim due to sublimation of grains could cause the rim to puff up, and it has been suggested that this leads to time-variable self-shadowing of the outer disc (Sitko et al. 2008; Wisniewski et al. 2008). Mid-infrared imaging of warm dust in the disc surface layers seems to indicate little flaring (Doucet et al. 2006), and scattered light brightness profiles are consistent with a non-flared outer disc (Wisniewski et al. 2008). Fits to the ISO-SWS spectra of HAeBe discs led Meeus et al. (2001) to classify HD 163296 as one of their group II objects, in which an optically thick inner disc shadows the outer disc from stellar radiation, resulting in a non-flared geometry.

In this study we aim to fit a wide variety of observed emission lines and the dust spectral energy distribution (SED) simultaneously with a disc with this observed geometry.

2. Herschel/PACS observations

We have obtained observations of HD 163296 in the far-infrared using Herschel/PACS, consisting of spectroscopic line observations and photometry. These include a detection of the [Oi] 63 μm fine structure line, and upper limits for a number of other atomic, ionic and molecular far-IR lines. We observed HD 163296 with the PACS spectrometer in line spectroscopy mode (obsid 1342192161, 1252 s, 3 repetitions) and range scan mode (obsid 1342192160, 5141 s, 3 repetitions), while chopping and nodding.

The spectra were reduced using HIPE version 7.0, with the standard pipeline scripts for a chopped line scan of a point source. The reduction steps include division by the spectral response function and MedianOffsetCorrection, as the new flat fielding task is not yet robust for short range scans such as the present observations (PACS instrument team, priv. comm.). For this paper, we only make use of the data contained in the central spaxel, to which we apply a wavelength-dependent correction factor for the loss of flux due to diffraction losses. With these additional steps, the absolute flux calibration accuracy is 30% (see PACS Spectroscopy Performance and Calibration Document, http://herschel.esac.esa.int/AOTsReleaseStatus.shtml).

The line fluxes and upper limits are given with continuum data in Table 1. The PACS photometry data are given in Table 2. Comparison of the azimuthally averaged radial profiles with the observed point spread function (PSF), using Vesta, show that the disc is unresolved at both 70 μm and 160 μm. The spatial extent of the observed [Oi] 63 μm emission is illustrated in Fig. 2, where it can be seen that the emission is dominated by the central 9.4 × 9.4 arcsec spaxel. This corresponds to a maximum radial extent of 560 AU for the emission.

The pipeline reduction of the [Cii] 158 μm range shows this line in absorption, which is unexpected. Therefore, we also looked at the unchopped spectra in the different chop positions. It is clear that the region around HD 163296 is full of [Cii] 158 μm emission, and that this causes the emission feature at the position of the star to cancel out (and even become an absorption feature). This is likely due to cloud material along the line-of-sight, since HD 163296 is located close to the Galactic plane, and is consistent with the presence of CO emission at radial velocities offset from this object, but unresolvable with PACS (Thi et al. 2001). We also inspected all the other lines for their unchopped appearance, but did not find contamination through the offset positions in those cases.

thumbnail Fig. 2

[Oi]63 μm spectra observed by Herschel/PACS towards HD 163296. Each spectrum corresponds to a 9.4 × 9.4 arcsec pixel centered at the given coordinates. Emission is dominated by the central pixel, corresponding to a maximum outer radius of  ~ 560 AU, assuming a distance of 118.6 parsecs to this object.

Table 2

Photometric fluxes observed with Herschel/PACS.

3. Properties of HD 163296

The study of the disc of HD 163296 is directly linked to the knowledge of the absolute parameters of the star itself and the behaviour of the SED in the ultraviolet (UV), optical and near infrared (near-IR), which is mostly dominated by stellar radiation; in particular, the energy emitted at wavelengths shorter than 2500 Å is a fundamental piece of information used to model the disc, due to its impact on the photo-chemistry and gas heating (Woitke et al. 2010). This section describes the stellar properties and studies the SED at the range mentioned above.

3.1. Stellar parameters

The determination of the stellar parameters of HD 163296 was done in an extensive work devoted to the study of a sample of Herbig Ae/Be stars (Montesinos et al. 2009). Details of the methodology followed and the observations used can be found in that paper. We outline in this section the relevant steps of the process.

An iterative distance-independent algorithm based on the analysis of the SED and optical high-resolution spectra (for the estimation of the effective temperature and metallicity), and mid-resolution spectra of the Balmer lines Hβ, Hγ and Hδ (for the estimation of the gravity), was applied.

The knowledge of the effective temperature, gravity and metal abundance allows us to place the star in a log g − log Teff diagram and superimpose the appropriate set of evolutionary tracks and isochrones for that specific metallicity (Yi et al. 2001). This gives directly – or after a simple interpolation between the tracks and isochrones enclosing the position of the star – the stellar mass and the age. Since there is a one-to-one correspondence between a pair (Teff,   g) on a given track and a pair (Teff,   L/L), the stellar luminosity can also be estimated.

Table 3 gives the stellar parameters of HD 163296 after a slight refinement of the results presented in Montesinos et al. (2009).

3.2. The SED: UV, optical, and near-IR observations

Simultaneous optical UBVRI and near-infrared JHK photometry obtained as part of the EXPORT collaboration (Eiroa et al. 2001; Oudmaijer et al. 2001), were used to build the optical and near-IR photospheric SED of HD 163296. Table 4 shows the wavelengths, magnitudes and corresponding fluxes with their uncertainties. The calibration from magnitudes to fluxes has been done using the zero points given by Bessell (1979) for the optical and Cox (1999) for the near-IR magnitudes.

Ultraviolet IUE spectra of HD 163296 were extracted from the INES archive1. From the collection of 68 SW and LW (for “short” and “long wavelength”) spectra in low resolution, obtained through the large aperture of the spectrograph, all those that looked clean (few or none bad pixels) and with the correct exposure classification codes were selected. The spectra (SW + LW) cover the range 1250–3000 Å.

A high-resolution, far-UV spectrum of HD 163296 was constructed by combining data taken with HST STIS and FUSE. The STIS spectra, covering 1150 Å to 3000 Å, were obtained from the HST STIS Echelle Spectral Catalog of Stars (StarCat2). The FUSE spectrum of HD 163296 is extremely noisy below 970 Å; therefore, we used only the portion between that wavelength and 1150 Å in our high-resolution composite UV spectrum for HD 163296. FUSE spectra contain strong terrestrial airglow emission lines, due to the large size of the observing aperture. We removed these lines from the HD 163296 spectrum before stitching it to the STIS spectrum.

A comparison of the IUE and FUSE+STIS spectra shows that they differ in the level of their continuum. The IUE spectra with the largest flux levels are a factor  ~1.6 higher than the FUSE+STIS spectrum, therefore, regarding the modelling, the ultraviolet contribution of the star will be represented by two different UV input spectra, namely a “low UV” state, corresponding to the FUSE+STIS spectrum, and a “high UV” state, corresponding to the average of the highest IUE spectra3 (see Sect. 4.1).

Table 3

Parameters for HD 163296.

Table 4

Photometry for HD 163296.

thumbnail Fig. 3

Spectral energy distribution for HD 163296. The green dots represent the fluxes corresponding to the UBVRIJHK photometry. The size of the dots is of the order or larger than the uncertainties. The IUE average spectrum and the FUSE+STIS spectrum are plotted in blue and red, respectively; the latter has been slightly smoothed to reduce the noise. The SpeX/IRTF spectrum is plotted in magenta. The black solid line is a photospheric model computed for the specific stellar parameters given in Table 3, reddened with E(B − V) =  + 0.15 and normalized at the flux in V. See text for details.

In addition, a near-infrared spectrum of HD 163296 was obtained at the NASA InfraRed Telescope Facility (IRTF) on Mauna Kea on 2010 July 8 with SpeX (Rayner et al. 2003). Ten individual exposures of the target, each lasting 32 s, were taken using the short-wavelength cross-dispersed (SXD) mode of SpeX. This mode yields spectra spanning the wavelength range 0.8–2.4 μm divided into six spectral orders. The slit width was set to 03, which yields a nominal resolving power of 2000 for the SXD spectra. Observations of HD 156717, an A0 V star, used as a “telluric standard” to correct for absorption due to the Earth’s atmosphere and to flux calibrate the target spectra, were obtained immediately preceding the observations of HD 163296. The airmass difference between the observations of the object and the standard was 0.06. The seeing was estimated to be  ~ at 2.2 μm during the observations and conditions were clear. The statistical S/N varies across the spectral range but is of the order of several hundred over the entire SXD spectrum. The SpeX data were reduced with Spextool (Cushing et al. 2004) and the telluric correction was performed using the procedures and software described by Vacca et al. (2003).

Figure 3 shows the SED of HD 163296. The solid black line is a photospheric model computed for Teff = 9250 K, log g = 4.07 and [Fe/H] = +0.20 from the GAIA grid of models created with the PHOENIX code (Brott & Hauschildt 2005). The model has been reddened with E(B − V) =  +0.15 (RV = 3.1); this will be discussed in Sect. 3.4.

3.3. Spectral type from the ultraviolet and near-infrared spectra

Despite the fact that the precise determination of the stellar parameters has been done through a careful analysis of a variety of observations, as we have mentioned in Sect. 3.1, an attempt to determine the spectral type of HD 163296 in a more qualitative way has been done using the ultraviolet and near-infrared spectra.

A comparison of the ultraviolet data with IUE spectra of stars classified as “spectral-type standards” was done. The spectra were taken from the “IUE Ultraviolet Spectral Atlas of Standard Stars”4 (Wu et al. 1983, 1992). It is interesting to note that below 2000 Å, for a given spectral type, the ultraviolet spectra of the “spectral-type standards” themselves are very different, in particular for A0 V and A1 V. If we consider only the region above 2000 Å, the comparison suggests that HD 163296 is closer to A3 V (and possibly slightly later), however, if a mild extinction is present, the UV spectra would be also close to that of A1 V stars. A larger extinction would push the spectral type determination towards earlier spectral types.

A similar exercise of spectral typing was done by a comparison of the SpeX spectrum in the YJ bands (0.85–1.35 μm) with Kurucz synthetic models (see Montesinos et al. (2009) for references regarding the synthesis codes and models used). The spectrum is complex, the range of Teff’s explored with the Kurucz models (8750–9250 K, corresponding to A3V and A1V) does not produce significant Paschen line variations. The Paschen line wings of HD 163296 are narrower than in any of the models, which would suggest temperatures of 9250 K or higher. The metallic lines are faint, their typical intensity being lower than 5% that of the continuum. Some lines are deeper than in any photospheric spectrum, which suggests the presence of some circumstellar absorption. Some lines are much fainter, which suggests an spectral type A1V or earlier, or dust emission veiling. Good agreement is rarely found. Although we fail to find a unique combination of Teff and log g that provides an adequate fit to both the UV and NIR spectra of HD 163296, all data are consistent with an A1-A3 star, seen through mild foreground extinction.

3.4. The extinction

The determination of the extinction is an intricate problem. The classical method of comparing the observed colours with those of standard stars of similar properties does not lead to any conclusive result, given the slight variability of the star and the photometric (and spectral type) uncertainties. For normal stars, the UV part of the SED can be used to estimate this parameter, but in the case of HD 163296 the ultraviolet spectrum is variable and has contributions both from the photosphere and the accretion shock that are difficult to quantify. The extinction seems to be mostly circumstellar, therefore the exact value of the constant RV is not known. The observations available do not allow us to disentangle completely this problem.

The strategy followed to assign a value of E(B − V) to this star has four steps: 1) introduce extinctions from E(B − V) = 0.0 to 0.2 to the photospheric model; 2) normalize the reddened model photosphere to the flux at V; 3) deredden the normalized model and estimate the integrated stellar flux, F; note that F is the flux that would be observed in the complete absence of extinction. Obviously, different extinctions imply different values of F. A value of RV = 3.1 has been used.

The fourth step is based on the following argument: since the stellar luminosity L/L is known through the distance-independent method outlined in Sect. 3.1, an estimation of the distance, d, can be obtained simply from the expression L = 4πd2F, where it is assumed that the whole star is visible, i.e. and no part of the stellar flux is completely blocked by the disc of other circumstellar material (see a complete discussion of this issue in Montesinos et al. 2009). The value E(B − V) =  +0.15 allows us to recover a distance of 128 pc, which matches the Hipparcos distance to within the uncertainties, therefore it has been adopted as the optimum extinction for consistency with the whole set of parameters. We would like to point out that that value is a simple parameterization of a very complex problem and that a deeper study, outside the scope of this paper, would be needed to study all its peculiarities.

In Fig. 3 the photospheric model reddened with E(B − V) =  + 0.15 has been plotted. It goes through the two ultraviolet spectra but falls below both of them at wavelengths shorter than  ~1260 Å, and it also underestimates the flux at U. Note that the model is purely photospheric, therefore the excess fluxes inferred from the observed data when compared with the model might be due to emission from the accretion shock, i.e. with a non-photospheric origin. The observed SED seems to depart from the photospheric SED at around  ~1 μm, where the contribution from the disc starts to become noticeable.

3.5. Herschel observational constraints

In addition to the new Herschel data we have further constrained our disc models using SMA interferometric observations of the J = 3–2 transition in 12CO and PBI interferometry of the 12CO 2–1 and 13CO 1–0 transitions (Isella et al. 2007), an upper limit for the S(1) transition in H2 obtained by VLT/VISIR (Martin-Zaïdi et al. 2010), the ISO-SWS spectrum in the infrared, and a wealth of photometric data, including SCUBA photometry in the sub-mm (Sandell et al. 2011). See Fig. 4 for a full list of the photometry sources.

thumbnail Fig. 4

Spectral energy distribution for the best-fit fully mixed disc model (“preferred” model), obtained from a simultaneous fit to the observed SED, ISO-SWS spectrum and line emission from HD 163296 (solid black line). Black dotted line indicates the SED for the model with power-law density profile, which is unable to fit the mm continuum image for HD 163296. Blue circles indicate (with increasing wavelength) simultaneous UBVRIJHK photometry (Eiroa et al. 2001; Oudmaijer et al. 2001), LM photometry (de Winter et al. 2001), IRAS photometry (12–100 mic), sub-mm photometry (Mannings 1994), and millimetre photometry (Isella et al. 2007). Also marked are PACS photometric observations (red circles), PACS continua derived from the spectroscopic observations (black circles), SCUBA photometry (green circles) (Sandell et al. 2011), scaled VLT/UVES spectrum (blue line, Martin-Zaidi, in prep.) and the ISO-SWS spectrum (green line). The red line shows the stellar+UV input spectrum. Downwards arrow denotes upper limit. All fluxes were corrected for interstellar reddening using the Fitzpatrick parameterisation (Fitzpatrick 1999) with RV = 3.1 and E(B − V) = 0.15.

Table 5

Fixed disc model parameters.

There is a wealth of observational data available for HD 163296 which constrains to varying extents our choice of model parameters. As well as fitting to the data mentioned in the previous paragraph, we have taken care not to contradict current understanding of the disc geometry. There is at present some uncertainty regarding the radial extent of the disc. Scattered light imaging indicates an outer edge for the dusty disc of  ~540 AU, consistent with that derived for the gas disc from millimetre CO emission by Isella et al. (2007). However, the resolved mm continuum emission seems to indicate a smaller outer disc radius  ~200 AU (Isella et al. 2007). Hughes et al. (2008) propose that the millimetre continuum and CO 3–2 emission can be fit simultaneously by a disc with an exponentially-tapered outer edge. For the purposes of our modelling, we adopt the surface density profile derived by Hughes et al. (2008) for this object, and attempt to fit simultaneously the remaining observational data (SED, line fluxes, CO line profiles) by varying the remaining disc parameters. As an additional test, we compute a model with a power-law density profile, varying in addition the column density power index to fit the same data (SED, line fluxes, CO line profiles). For the purposes of this power-law model the disc outer edge is fixed at the scattered light radius of 540 AU. These adopted disc outer radii are consistent with our far-infrared photometry measurements, in which the disc is unresolved, corresponding to the emission at these wavelengths being dominated by the inner  ~200–300 AU of the disc.

Table 6

Model parameters for four well-fitting disc models.

The disc is also resolved at the inner edge, and while there is some uncertainty regarding the precise structure in this region (Renard et al. 2010; Benisty et al. 2010; Eisner et al. 2009; Tannirkulam et al. 2008b), we adopt the median value of 0.45 AU for the inner disc radius. The disc is assumed to have constant flaring with radius, i.e. the assumed inner rim structure does not cause shadowing. The reference scale heights referred to in our results (see Table 6) are defined at the inner rim, i.e. R = 0.45 AU. The disc inclination is fairly well-constrained by imaging at various wavelengths (Benisty et al. 2010; Wassell et al. 2006; Isella et al. 2007; Tannirkulam et al. 2008a; Grady et al. 2000), and is fixed at 50° throughout the modelling in this paper. See Table 5 for a full list of the parameters which are fixed during the model-fitting process.

While we have not used the spatial CO emission maps as a constraint for our modelling efforts, we do use the extracted line profiles for the three CO lines (Isella et al. 2007) as a constraint for the radial origin of the emission. The interferometry allows us to isolate the emission from the disc, eliminating contamination from cloud material. The bulk of the modelling uses a surface density profile derived from fitting to the spatial CO 3 − 2 emission (Hughes et al. 2008).

4. Modelling

The disc modelling was carried out using the radiation thermo-chemical disc code ProDiMo (Woitke et al. 2009; Kamp et al. 2010). The code solves in turn the radiative transfer problem, chemical network, and gas heating and cooling balance. Finally, the level populations are calculated, followed by line transfer calculations to give the predicted line emission.

The frequency-dependent 2D radiative transfer solver calculates the dust temperature structure and internal continuous radiation field for a given disc and stellar spectrum. For the purposes of this paper it was used to determine an SED-fitting dust model for HD 163296, by varying the total dust mass, grain size parameters and the spatial distribution of dust in the disc, including dust settling as a function of grain size (see Sect. 5.1).

Once the dust temperature structure Tdust(r,z) and internal radiation field Jν(r,z) have been computed, ProDiMo solves the chemical network assuming kinetic equilibrium (such that the net concentration of each species does not change over time), and the gas energy balance. This gives the chemical composition of the disc, the gas temperature Tgas(r,z) and the gas emission lines, and depends on the gas-to-dust ratio and PAH abundance. The chemical network includes 973 reactions between 76 gas phase and solid ice species composed of 10 elements (Woodall et al. 2007; Schöier et al. 2005). There is a detailed treatment of UV-photoreactions (see Kamp et al. 2010), as well as H2 formation on grain surfaces, vibrationally excited H chemistry, and ice formation (adsorption, thermal desorption, photo-desorption, and cosmic-ray desorption) for a limited number of ice species (see Woitke et al. 2009, for details).

The level populations of the various atoms, molecules and ions are calculated using an escape probability method, and used to solve the line transfer for selected spectral lines, in this case those of CO (Flower 2001; Jankowski & Szalewicz 2005; Wernli et al. 2006; Yang et al. 2010), H2O (Barber et al. 2006; Dubernet & Grosjean 2002; Faure et al. 2007; Green et al. 1993; Tennyson et al. 2001), Oi (Abrahamsson et al. 2007; Bell et al. 1998; Chambaud et al. 1980; Jaquet et al. 1992; Launay & Roueff 1977) and Cii (Flower & Launay 1977; Launay & Roueff 1977; Wilson & Bell 2002). See Woitke et al. (2011) for a description of the line transfer calculations, as well as recent improvements to the chemistry and gas heating/cooling balance in ProDiMo.

4.1. Input spectrum

The stellar input spectrum used for the modelling is a composite of the available observed UV data with a model photosphere spectrum in the wavelength range for which detailed spectral data are not available. We have accounted for the observed ultraviolet variability of this object (see Sect. 3) by running two sets of disc models with different UV input spectra, namely, one representing the “low UV” state and the other a “high UV” state, as described in Sect. 3.2. These were provided by the FUSE+STIS and the IUE average spectra, respectively. At the upper wavelength boundary we switch in each case to the GAIA PHOENIX model photosphere computed with Teff = 9250 K, log g = 4.07 and [Fe/H] = +0.20. The observations were de-reddened using the Fitzpatrick parameterization (Fitzpatrick 1999) with E(B − V) =  + 0.15 and RV = 3.1.

The UV spectrum is of central importance to the gas modelling, both in the chemical network and in the gas heating/cooling balance. This requires UV input for wavelengths of 912 Å and above (the wavelength at which atomic hydrogen is ionised). The FUSE spectrum in this case only extends down to 970 Å, and the IUE spectrum down to 1150 Å. In order to give a full spectrum for modelling purposes, the “high UV” input is extended down to 970 Å by scaling the STIS+FUSE (low UV) spectrum accordingly. In both cases a power-law fit to the data is employed for the remaining shortfall. The total photon particle fluxes are 60% greater for the high UV input spectrum than for the low UV state.

4.2. Fitting procedure

An evolutionary χ2-minimisation strategy was employed to simultaneously fit the observed SED, ISO-SWS spectrum and gas line observations (see Sect. 3.5). In the case of the line observations the models were fitted to the observed fluxes, and the line widths in cases where profiles were available. The upper limits were not used as a fitting constraint, but were checked against the predicted best-fit model fluxes for any contradictions.

We have utilised a Monte-Carlo evolutionary method, varying the 11 parameters listed in Table 6 until the minimum χ2 value was reached (see Woitke et al. 2011 for more information regarding the evolutionary strategy and method of χ2 calculation). We note that the best-fit parameters depend in general on the weighting scheme employed for the three groups of observations, and in practice it was found that weighting slightly in favour of the line emission and ISO-SWS observations at the expense of the photometry gave the best fit overall. We also note the possibility that the converged parameter values correspond to a local χ2 minimum, and not a single global minimum, and we do not claim a unique fit to the observations. Indeed, it is clear from our modelling efforts that there exist a lot of degeneracies between the various parameters. However, each of the models discussed in the following sections are the result of several hundred generations of χ2-minimisation.

4.3. Dust properties

The rich solid-state spectrum of HD 163296 was observed by van den Ancker et al. (2000), using the ISO-SWS and ISO-LWS spectrometers. This study found the dust in this object to consist of amorphous silicates, iron oxide, water ice and a small fraction of crystalline silicates, with the presence of large millimetre-sized grains indicated by the continuum temperature.

Table 7

Dust grain composition.

We assume homogeneous and spherical dust grains (Mie theory), with a power-law size-distribution defined by a minimum radius amin, maximum radius amax, and power-law index p. The grain composition is assumed to be constant throughout the disc, and we adopt the dust species mixture determined by Bouwman et al. (2000) for this object, averaging their fractional species abundances over the disc (see Table 7). All grains have the same species composition regardless of their size, and the grain composition is not treated consistently with the physics and chemistry of grains. For instance, the assumed grain water ice fraction is constant throughout the disc, independent of the water ice abundances calculated in the solution of the chemical network. For the dust material mass density, we take the average value for this mix of 3.36 g cm-3. Optical constants for the various dust species were taken from measurements made by the studies listed in Table 7.

Bouwman et al. (2000) carried out a detailed study of this object, fitting to spectral data over a large wavelength range. We note however that the dust opacity law represents a potentially large source of uncertainty when deriving the disc parameters. For instance, the two olivine species FeMgSiO4 and Mg2SiO4, whilst having broadly similar spectra in the wavelength range observable by ISO-SWS, have a factor  >30 difference in absorption opacity at  ~1 micron. This will inevitably affect the derived disc parameters, although the difference in opacities at millimetre wavelengths is negligible. The models referred to in this study have extinction opacities in the range (10.4–12.8) cm2 g-1 (dust) at 1 mm.

4.4. Gas/dust ratio

The models considered have constant dust properties with radius. However, in general the models do allow for differential dust grain scale heights as a function of grain size, to represent the major effects of dust settling. This leads to a local gas/dust ratio which varies with height in the disc, yet at any given radius the ratio of the vertical gas and dust total column densities is constant, and equal to the global gas/dust ratio. The degree of dust settling is determined by a simple parameterisation (see Sect. 5.1), and models in which no dust settling is present are referred to as “fully mixed”.

5. Results

Table 6 gives the derived parameter values for four runs of the evolutionary scheme outlined in the previous section. The first three columns refer to models with a fixed exponentially-tapered density profile. The fourth column refers to a model with power-law density profile, Σ ∝ R − ϵ, in which the outer radius is fixed at 540 AU but the power index, ϵ, is allowed to vary.

thumbnail Fig. 5

1.3 mm continuum emission maps for HD 163296.From left: observed map; emission map for preferred model (Col. 3 in Table 6); residuals for preferred model; emission map for power-law model (Col. 4 in Table 6); residuals for power-law model. Residuals computed by subtracting the model intensity from the observed intensity. Contours spaced at 12 mJy intervals, corresponding to [3,6,9,12,15,18,21,24,27,30,33,36,39,42,45] ×    σ.

The surface density profile in the models with an exponential outer edge (Cols. 1–3) takes the form (1)where γ is analagous to the power law index, and R0 is the scale length over which the disc surface density tapers exponentially. We adopt the previously derived values of γ = 0.9 and R0 = 125 AU (Hughes et al. 2008). For the purposes of our modelling, the disc chemistry etc. was computed out to 850 AU, at which point the column density is negligible.

The first, “low UV” column uses the FUSE+STIS UV spectrum as input for the modelling, and the second “high UV” column uses the average IUE spectrum, as described in Sect. 4.1. Both runs produce settled discs, with variable scale heights for dust grains of different sizes. The third column in Table 6 gives the results for a run with well-mixed dust grains, i.e. no settling is introduced. This run uses the FUSE+STIS (low UV) spectrum, as does the run with power-law surface density referred to in Col. 4.

None of the models described in Table 6 represent a perfect fit to all the available data for HD 163296, but in all cases the models are able to fit almost all the observations. Considering the large number of data used for the modelling, and the range of parameters covered by the well-fitting models, is is clear that there is some degree of parameter degeneracy present.

5.1. Continuum emission and spectral energy distribution

The SED for the best-fit fully mixed disc model (Col. 3 in Table 6) is shown in Fig. 4. With the constraints outlined in previous sections, our aim was to fit the observed SED with a simple dust model, i.e. a continuous disc with constant flaring, and dust species composition constant throughout the disc. All of the models described in Table 6 require the presence of mm-sized grains in order to fit the observed SED. It was not possible to fit both the 10 micron silicate feature and the millimetre tail with a well-mixed disc. This can be seen in Table 6, where the fully mixed model gives a worse fit to the ISO-SWS spectrum than for the models in which dust settling is present. Figure 4 shows that the best-fit fully mixed model gives a smaller than observed silicate emission feature. Sitko et al. (2008) observe a spectral variability of  ~10% in the wavelength range covered by the ISO-SWS data, and while we are unable to fit the observations to within this range with even a settled model (simultaneously with the overall SED and line data), the fit is considerably better than is possible with a well-mixed disc. This result is similar to that found for IM Lupi by Pinte et al. (2008), in which a settled disc was needed to fit the 10 μm silicate feature and SED millimetre tail simultaneously.

On the other hand, the constraint of fitting to the ISO-SWS spectrum leads to models which produce too much emission in the near-infrared (J,H and K bands), giving an overall worse fit to the observed photometry. This effect is less pronounced in the best-fit fully mixed model (Col. 3 cf. Cols. 1 and 2 in Table 6). The flux overprediction in the settled models is likely mainly due to the smaller average grain size in the strongly illuminated disc surface layer. The condition of radiative equilibrium between grains means that this gives a higher average grain temperature at the disc surface, and the re-emitted light peaks at shorter wavelength than for the fully mixed models (in the near-IR as opposed to the mid-IR). This flux overprediction is still present to a smaller extent in the fully mixed model (see Fig. 4). This could be a consequence of our fixing the disc inner radius for modelling purposes, and a better fit in the near-IR would be possible if we had allowed the inner radius to vary. It is also probable that the structure of the inner rim is more complex than allowed by our parameterised model, with evidence for puffing-up of material, which would certainly affect the emission at these wavelengths. There is also significant variability observed for this object in this wavelength region (de Winter et al. 2001; Sitko et al. 2008).

Dust settling enhances the silicate emission feature for an optically thick disc since it removes the large grains from the surface layers, and places them at smaller heights in the disc where they cannot be observed at mid-infrared wavelengths. The flat blackbody opacity of these larger grains is overwhelmed by the characteristic spectrally-varying opacity of smaller micron-size grains which remain at lower optical depth, dominating the observed emission.

In cases such as this where a simple parameterised disc structure is assumed, ProDiMo implements a simple recipe to account for the major effects of vertical dust settling. We assume that the dust grains are distributed vertically with a scale height which decreases for large dust particles: (2)where H(r) is the gas scale height, and s and as are two free parameters (see Woitke et al. 2010 for details).

Introducing two additional parameters to the modelling inevitably improves the fit to the observed data, and indeed the combined fit to the photometry, ISO-SWS spectrum and line emission does improve overall as dust settling is introduced. The fact remains, however, that it is possible to obtain a good fit to almost all the available data with a well-mixed disc, with a gas/dust ratio almost equal to the canonical value of 100, and so we adopt this model, described in Col. 3 of Table 6, as our “preferred” model. However, the inability of a well-mixed disc such as this to fit the observed emission in the mid-infrared means that there remains compelling evidence for a disc exhibiting dust-settling, able to fit the line emission with a depleted gas/dust ratio (see Col. 1 in Table 6). The effect of dust settling on the line emission is further explored in Sect. 5.2.2.

Table 8

Integrated line fluxes for the transitions observed by PACS, and additional CO and H2 lines, as predicted by the models listed in Table 6.

The model with power-law density profile (Col. 4 in Table 6) gives the best overall fit to the SED, ISO-SWS spectrum and line emission (see Fig. 4), with quite different disc parameters to the three models with exponential outer edges (Cols. 1–3). The disc requires a very flat density profile, ϵ    ~    0.085, in order to fit these observations (not including spatial data). The predicted millimetre continuum emission from this model is plotted with the observed maps in Fig. 5. It is clear that the flat density profile leads to an emission deficit in the inner disc in comparison to observations, and so this cannot be considered to be a good model for the disc of HD 163296. A better power-law model could be found by using the maps as a further constraint for the modelling, but that is beyond the scope of this paper. The fully mixed “preferred” model gives a better match to the observed millimetre emission, as should be expected since its density profile is derived from fitting to this data (Hughes et al. 2008). The factors driving the power-law model to such a flat density profile will be discussed in Sect. 5.2.

The power-law density profile model is, however, able to fit the non-spatially-resolved data with a smaller flaring index (~1.02) than the exponentially-tapered discs (~1.07), more in keeping with the observational evidence for a disc that is not strongly flaring (Meeus et al. 2001; Doucet et al. 2006; Wisniewski et al. 2008). Meeus et al. (2001) classify HD 163296 as a group II object, in which the inner disc shadows the outer disc, as opposed to group I objects which have a flared geometry. Meijer et al. (2008) find HD 163296 to lie on the transition between flared and non-flared geometry, i.e. a flaring index close to 1. It is clear that our “preferred” model is unable to account for every observed property of the disc of HD 163296, and it is probable that a more complex model – with non-constant flaring and variable dust properties with radius – would be required to simultaneously fit the full set of observational data available for this object.

The models described in Table 6 have disc dust masses in the range (7−12) × 10-4   M, which is within the range of masses (5−17) × 10-4   M found in the literature for HD 163296. This small spread in dust mass from our fitting efforts is unsurprising given the fixed disc size and constant grain composition.

5.2. Gas properties

thumbnail Fig. 6

Plotted in black are the line profiles observed by Isella et al. (2007) for the 12CO J = 3–2 (upper panel), 12CO J = 2–1 (middle panel) and 13CO J = 1–0 (lower panel) transitions, with the corresponding profiles from the preferred disc model in red. This refers to a single simultaneous fit to the observed continuum and line data (Col. 3 in Table 6). Observed profiles are obtained by integrating over the whole disc.

The predicted line fluxes for the various models are given in Table 8. We have ignored the [Cii] 158 μm result during the model-fitting process, due to systematic uncertainties arising from strong emission at the offset positions.

The observed line profiles for the 12CO J = 3–2, 12CO J = 2–1 and 13CO J = 1–0 transitions (Isella et al. 2007) are plotted with the preferred model profiles in Fig. 6. The observed profiles are double-peaked, consistent with a disc in Keplerian rotation. The CO 3–2 line is observed to be brighter than predicted by the model. This behaviour is repeated across the three exponentially-tapered models, but is less pronounced in the power-law model, whose flat surface density profile gives a larger CO 3–2 flux. Discs with flatter density profiles allow stellar radiation to penetrate deeper into the inner disc, resulting in a slight increase in temperature throughout the disc, including its outer regions. This leads to brighter CO 3–2 lines since this transition is optically thick in all of our models, and is largely dependent on the disc outer radius (which is fixed) and the outer disc temperature at intermediate height (see Fig. 8). The CO 2–1 line also follows this behaviour, with the line flux ratio in our models staying roughly constant at a value of CO 3–2/1–0  ~  3, as expected for lines formed under optically thick LTE conditions (Kamp et al. 2010). This is slightly lower than the observed ratio of CO 3–2/2–1 = 4.3. We note that recent observations of this object give a peak intensity of 25 Jy for the 3–2 line (Hughes et al. 2011)5, compared with 43 Jy in the Isella et al. (2007) data and the 30 Jy predicted by our preferred model. Our predicted flux is within the calibration error margin for both observations, but fitting to the lower value would in general tend to lead to power-law models with steeper density profiles, in contrast to the flat power-law profile we obtain from fitting to the brighter observation. Another factor driving the power-law models to flat density profiles is the overprediction in near-IR continuum emission by our models. The inner radius is fixed in accordance with high resolution imaging, and so the power-law models tend to reduce their near-IR emission by removing material from the inner disc, giving a better fit to the photometry.

PAHs are one of the main sources of gas heating in the disc, as UV photons cause them to emit excited electrons via the photoelectric effect, which thermalise in the gas. For the purposes of our modelling we consider a typical size of PAH molecule (circumcoronene, NC = 54 carbon atoms and NH = 18 hydrogen atoms) and include PAH, PAH, PAH+, PAH2+ and PAH3+ as additional species in the chemical network (see Woitke et al. 2010 for further details). The fractional PAH abundance, fPAH, is defined relative to the standard ISM particle abundance with respect to hydrogen nuclei, (Tielens 2008), so that the PAH abundance, . Small fPAH values (0.007–0.04 relative to ISM abundances) are required in our models in order to fit the observed line data. This is consistent with the findings of Geers et al. (2006), who require PAH abundances fPAH < 0.1 in order to fit the observed PAH emission from Herbig discs. PAH emission was not observed in HD 163296 by Acke et al. (2010), and we have used the Monte Carlo radiative transfer code Mcfost (Pinte et al. 2006) to check that our results are consistent with this non-detection. A disc model with parameters and PAH abundance identical to our preferred model gives an infrared spectrum in which no PAH emission is seen. This is due in part to the low flaring in the disc, β = 1.066, and the large fraction of small dust grains in the disc, with significant opacity in the optical-UV range over which PAH molecules absorb. This results in the PAHs being “hidden” from direct stellar illumination, absorbing a fraction  ~10-4 of the energy absorbed by the disc. In our preferred model 0.3% of the total carbon mass is tied up in PAH molecules.

The total hydrogen number density and gas temperature within the disc are plotted in Fig. 7, and the spatial origin of the various emission lines is visualised in Fig. 8. These plots all refer to the “preferred” disc model, i.e. one in which no dust settling is present. The CO J = 3–2 and [Oi] 63 μm lines are optically thick throughout the disc, and cannot be used alone to trace the gas mass. Even the 13CO line is optically thick throughout much of the disc, only becoming optically thin outside of  ~400 AU. Also, there is evidence that this line (and most others) can be affected by the degree of dust settling in the disc (see Sect. 5.2.2). The [Cii] 157.74 μm line is optically thin, but traces only the warm ionised gas in the disc surface. Clearly care must be taken when attempting to use individual emission lines as a tracer of gas mass. Note the contrast in disc radii at which line and continuum become optically thin in the case of CO 3–2, reflecting the conflicting derived radii for the gas and dust discs in HD 163296. The line optical depths in Fig. 8 were computed using an escape probability formalism, and do not take into account the effects of Keplerian shear in the disc.

As well as optical depth effects, disc gas mass estimates from CO emission alone can be affected by the extent of CO freeze-out on grain surfaces. ProDiMo calculates the grain adsorption and desorption as part of the solution of the chemical network, and the spatial abundances of gas-phase and ice CO in our best-fit model are plotted in Fig. 9. The total mass of gas-phase CO in the disc is 1.51  ×  10-5   M, and the CO ice mass is 4.32    ×    10-5   M. We would expect a smaller fraction of CO to be frozen out on grain surfaces in discs with flatter density profiles, due to weaker dust shielding in the inner disc leading to warmer temperatures throughout the disc, and less CO freeze-out.

The derived chemical abundances and gas properties in our preferred model have been checked by re-computing the chemistry using a time-dependant solver, and comparing the results to those obtained with the assumption of kinetic equilibrium. We assume molecular cloud initial abundances, and after running the solver for 4 Myr there is no major departure from the equilibrium chemistry. The biggest change in line flux is a reduction of 20% in the 180.4 μm water line, while the [Oi] 63 μm line decreases by 4% and the CO lines by  <1%. The assumption of a constant 13CO/12CO ratio is also valid since dust shielding dominates over CO self-shielding in the models, with negligible change to the results when CO self-shielding is swtiched off. We would therefore expect there to be no fractionation effect present in this disc.

All of the disc models considered are passive discs, i.e. the viscous “α” heating parameter is set to zero. This is not strictly consistent with the presence of accreting gas in this object, since this implies some form of viscosity to transport angular momentum through the disc. However, viscous heating is likely to be unimportant in the case of Herbig Ae discs such as this, which are thought to be dominated by radiative heating by the central star (D’Alessio et al. 1998).

thumbnail Fig. 7

The total hydrogen number density (left panel) and gas temperature (right panel) are plotted for a vertical cross-section through the preferred disc model (Col. 3 in Table 6). Dashed lines indicate contours of gas temperature and visual extinction AV, as marked.

thumbnail Fig. 8

Spatial origin of the various gas emission lines in the preferred model (Col. 3 in Table 6). From top-left clockwise: CO J = 3−2, CO J = 18−17, [Cii] 157.74 μm, o-H2 S(1) (with red Tg contours), [Oi] 63.18 μm, 13CO J = 1–0. In each panel, the upper plot shows the line optical depth as a function of radius (blue line) and the continuum optical depth at the corresponding wavelength (black line). The middle plot shows the cumulative contribution to the total line flux with increasing radius. The lower plot shows the gas species number density, and the black line marks the cells that contribute the most to the line flux in their vertical column. The two vertical dashed lines indicate 15% and 85% of the radially cumulative face-on line flux respectively, i.e. 70% of the line flux originates from within the two dashed lines.

thumbnail Fig. 9

CO abundances in the preferred model (Col. 3 in Table 6). Left panel: gas-phase CO abundances. White and blue dashed lines plot the gas and dust temperature contours respectively. Right panel: CO ice abundance.

5.2.1. Effect of UV variability

We have investigated the effects of UV variability on the gas in the disc in two ways. Firstly, by computing a test model with parameters identical to our “preferred” model, but with the high UV state spectrum used as input. This allows us to isolate the effect of increasing the UV intensity on the line emission, since all other model parameters remain the same. The results of this test are summarised in Table 9. The second approach consists of a separate evolutionary run of models using the high UV state as input, and has been covered in Sect. 5. This method allows us to assess how the disc model parameters change in order to fit the same line data with an increased UV intensity (see Table 6). By examining firstly the physical effect of the UV on the gas in the disc, and secondly the effect of the UV on the model parameter fit, we can estimate the degree of uncertainty introduced by UV variability in this object.

An increased UV strength affects the gas chemistry by promoting photochemical reactions, and the photodissociation of H2 and CO. It also leads to increased desorption of ice species from grain surfaces in regions where the UV is able to penetrate. It also strongly affects the gas heating, via the photoelectric effect in dust grains and PAHs.

The “high UV” input spectrum represents a 60% increase in UV intensity over the low UV input. Table 9 summarises the effects of this increase on the gas in the preferred model. There is a slight overall increase in the disc gas temperature, and an increase in the atomic hydrogen and C+ masses, due to photoelectric effects and PAH heating, and increased photodissociation of gas molecules. The total gas-phase CO mass actually increases slightly, since the increased rate of photodissociation is balanced by an increase in the rate of ice desorption. The J = 36–35 line in particular sees the biggest increase in flux,  ~50%. By considering only the gas in the warm layers of the inner disc, where the high J CO lines form, it can be seen that this increase in line flux is caused by a substantial increase in the gas temperature in this region. The [Oi] 63 μm flux increases by  ~30%, roughly equal to the calibration uncertainty margin, and the low J CO lines increase by  ~5%. All non-detected lines stay within the observed upper limits. In summary, it would not be possible to distinguish between the high and low UV states from these disc observations, since all fluxes are within the observational error margins. The general influence of the UV intensity on the disc gas is discussed by Woitke et al. (2010), in the context of a large grid of models.

The derived parameters for the evolutionary run using the “high UV” spectrum as input are listed in Table 6. They are largely similar to those derived from the low UV run (Cols. 2 cf. 1), as would be expected from the small fractional change in UV intensity. The slight decrease in gas mass and flaring will tend to counteract any increase in line emission. The PAH fractional abundance decreases slightly as might be expected. The main change lies in the grain size distribution, with a reduction in the minimum grain size. This makes it harder for the UV to penetrate the disc, balancing the increase in UV intensity. This also produces even more continuum emission in the near-IR, worsening the SED fit still further in comparison to the “low UV” model.

5.2.2. Effects of dust settling

Table 9

Effects of UV variability. Comparison between the preferred model, and an identical model illuminated by an input spectrum typical of a “high UV” state for this object.

thumbnail Fig. 10

The fractional enhancement in line flux is plotted against the minimum grain size affected by dust settling, relative to the preferred (fully mixed) model (as = amax). All models have a settling parameter of 0.5. Top panel: CO lines; middle panel: water lines; lower panel: [Oi]63 μm, [Cii]158 μm, OH 79.11 μm and H2 S(1) lines.

The effect of dust settling on various line fluxes can be seen in Fig. 10. This shows the increase in line flux across the various transitions as smaller and smaller grains are allowed to settle towards the midplane. The models are otherwise identical to the preferred (fully mixed) model, but with a settling parameter of 0.5. The models run from strongly settled discs through to entirely well-mixed, i.e. identical to the preferred model. For the purposes of this exercise we have computed some extra line fluxes in addition to those observed in HD 163296, namely the 12CO 6–5 and 1–0 transitions, the 13CO 3–2 transition, and the o-H2O 538.3 μm transition.

Kamp et al. (2011) computed a large grid of disc models using ProDiMo, in tandem with Mcfost (Pinte et al. 2006). This study noted a trend of increasing disc dust temperature with settling, due to a combination of reduced emissivity at long wavelengths in the disc surface, and increased illumination of the dust as the stellar radiation is able to penetrate further into the disc. Our models of HD 163296 follow this same behaviour, albeit with an accompanying increase in gas temperature. The effect of dust settling on the gas and dust temperature contours in our models is demonstrated in Fig. 11, where the chemical abundances of various chemical species are plotted. The chemical structure is seen to “follow” the settled dust grains towards the midplane, since the UV penetrates further into the disc, causing the characteristic layering of the various gas species to move closer to the midplane. This is despite the gas scale heights being identical in both models. In the settled models the warm gas in the disc resides at lower heights where the gas density is higher, leading to a general brightening of the emission lines. This seems to be a firm result for this object, with serious implications for the derivation of a precise disc gas mass from the emission lines.

The brightening of the line fluxes seems to be a general result for all the transitions considered. The exact behaviour of the various line fluxes with variable dust settling, and the extent to which they are affected, depends on the height and radial position in the disc from which they originate. This is well-illustrated by the CO lines in Fig. 10. As the largest grains begin to settle, the 13CO 1–0 line is the first to show an increase in flux. This is due to an increase in thermal desorption of CO ice from grain surfaces in the outer disc midplane, which this line traces (see Fig. 8). This behaviour is echoed in the other 13CO line, the J = 3–2 transition. In contrast, the 12CO 3–2 line, which is also formed in the outer disc, shows a smaller increase in flux. This is because this line is optically thick, and formed higher in the disc where it remains unaffected by the thermal desorption of ice in the midplane. The 12CO 2–1 and 1–0 lines follow this same behaviour (not plotted). As smaller grains begin to settle, the 13CO lines increase less rapidly than for the lines formed in the warm gas higher in the disc. This is most pronounced in the high J CO lines, which form in the inner disc where the temperature and density gradients are more extreme, and the lines more sensitive to the downwards shift in chemical structure. The 12CO 6–5 line forms at intermediate distance in the disc, and this is reflected in the level of flux enhancement (Fig. 10; top panel).

The behaviour of the water lines (middle panel in Fig. 10) is less straightforward. The o-H2O 179.5 μm line shows a stronger increase in flux than both the 180.5 μm line, which is at higher excitation, and the 538.3 μm line, which is at lower excitation. The p-H2O 89.99 μm and o-H2O 78.7 μm lines are at higher excitation, forming in the inner disc, and these also show a large increase in flux as the dust grains settle. It is likely that the relative behaviour of the various water lines is a complex function of the changing temperature structure and H2O chemical abundance structure in the disc, affecting the excitation of the various levels.

The largest increase in line flux occurs in the H2 S(1) transition (Fig. 10; lower panel). The emission region for this line is seen to follow the 300 K gas temperature contour (see Fig. 8), which intersects a greater molecular hydrogen mass in the settled models. This leads to a dramatic enhancement in flux in the strongly settled models, a factor of  ~30 increase on the fully mixed model. The [Oi] 63 μm line also starts to increase as the smallest grains are allowed to settle, coinciding with an increase in atomic gas, but the enhancement is less extreme than since this line is formed further out in the disc (see Fig. 8). The OH 79.11 μm line is formed at intermediate distance between H2 S(1) and [Oi] 63 μm, and its behaviour with increasing dust settling reflects this. The [Cii] 158 μm line is less strongly affected by settling. This line is formed in a thin layer at the disc surface, with emission dominated by the outer disc, and so it is less sensitive to changes in the internal disc temperature structure.

The settling of dust grains might be expected to drive line formation conditions closer to LTE, with the various species residing lower in the disc, in regions of higher density, with more frequent collisions between particles. However, there is no evidence for this in our models. The majority of lines show only minor departures from LTE, in both settled and well-mixed models. The largest departure from the line fluxes assuming LTE level populations occurs in the water lines, but the LTE line fluxes increase at roughly the same rate as those computed using escape probability, so the ratio FNLTE/FLTE remains roughly constant. The critical density for the water lines, ncrit ~ 108−1010 cm-3. These densities are reached only in the inner disc in our models, with each of the lines arising in regions with n < ncrit, even in the settled models. The water line which is closest to LTE in our models is the 538.3 μm line, which has the lowest level of excitation and the smallest critical density. In general, the critical densities decrease rapidly at large optical depths, driving most of the transitions towards LTE.

These results are in contrast to the findings of Jonkheid et al. (2007), who find a trend of decreasing gas temperature and line emission as the dust grains are allowed to settle towards the midplane. However, a direct comparison is not appropriate, since their models assumed a simultaneous reduction in the dust/gas ratio and PAH abundances with dust settling. The reduction in gas temperature can largely be attributed then to the drop in photoelectric heating from dust grains and PAH molecules. Jonkheid et al. (2007) also note that the dust temperatures in their models are high enough to prevent CO freeze-out, whereas considerable freeze-out is present here. Meijerink et al. (2009) also predict an increase in line emission with dust settling, from models in which the effect of settling is simulated by simply increasing the global gas/dust ratio, as opposed to considering the vertical distribution of grains in the disc.

thumbnail Fig. 11

Effect of dust settling on the disc chemical abundances. Left column indicates abundances in a fully mixed model, and right column represents a strongly-settled model with otherwise identical parameters. Top row: CO abundance with dust (blue dashed lines) and gas (white dashed line) temperature contours. Middle row: H2 abundance with 300 K gas temperature contour (black dashed line) and visual extinction contour (white dashed line). Bottom row: H2O abundance. White dashed lines indicate Tgas contours enclosing “hot water” region, red dashed lines indicate contours of UV field strength per hydrogen nucleus, log (χ/nH) as defined in Draine & Bertoldi (1996), enclosing the cool water belt.

Our models indicate that the 13CO 1–0 line flux can change by a factor  ~2 in strongly settled models. This is in a high mass disc for which the line is optically thick throughout most of the disc, and so this settling flux enhancement might be expected to be even stronger in a lower mass disc. In any case, variable grain freeze-out would seem to limit the ability of this line to trace the total disc gas mass. In general, the effect of dust settling on the vertical thermal structure of the gas in discs is seen to introduce new degeneracies when attempting to fit the disc parameters to the observed line emission. While our preferred model fits the observed data with a well-mixed disc with the canonical gas/dust ratio, we cannot rule out the possibility of a gas-depleted disc in which dust settling gives an enhancement in the various line fluxes. Indeed, such a disc with gas/dust  ~20 allows a better fit to the line data and the observed 10 micron silicate emission (Col. 1 in Table 6).

5.3. Effect of X-rays

It is currently unclear how important a role X-rays play in determining the gas chemistry and temperature structure in discs. Aresu et al. (2011) studied the effects of stellar X-ray emission on models of T Tauri discs computed with ProDiMo. This study found that while Coulomb heating by X-rays introduced an extended hot gas surface layer to the discs, the [Oi] and [Cii] fine structure lines were only affected for X-ray luminosities LX > 1030 erg s-1. One would expect this to hold for the warmer gas in Herbig discs, and Günther & Schmitt (2009) derive an X-ray luminosity of LX = 1029.6 erg s-1 for HD 163296, lower than but close to the threshold value noted by Aresu et al. (2011) for T Tauri discs. HD 163296 represents an interesting test case for the influence of X-rays in Herbig discs.

We have computed a disc model with parameters identical to our preferred model, but with an additional X-ray component in the input spectrum. The X-ray luminosity was set equal to the value of LX = 1029.6 erg s-1 observed by Günther & Schmitt (2009) for HD 163296. The effect of this on the gas temperature is illustrated by Fig. 12 (cf. Fig. 7). The X-ray heating processes produce an extended hot surface gas layer, as noted by Aresu et al. (2011), with gas temperatures  >5000 K extending out to the outer disc.

There is little discernible effect on the disc chemistry, with none of the chemical species masses or predicted line fluxes changing by more than a few percent. In all cases the effect of the additional X-rays is less than produced by switching to the “high UV” input spectrum . We would therefore expect the UV to dominate the gas chemistry in a Herbig disc of the sort considered here.

6. Summary and conclusions

This paper presents new observations of the far-IR lines of the Herbig Ae star HD 163296, obtained using Herschel/PACS as part of the GASPS open time key program. These consist of a detection of the [Oi] 63.18 μm line, as well as upper limits for eight additional lines. We have supplemented these observations with additional line and continuum data, as well as high resolution observations of the variable optical spectra obtained as part of the EXPORT collaboration.

We have computed radiation thermo-chemical disc models using the disc code ProDiMo, and employed an evolutionary χ2-minimisation strategy to find the best simultaneous model fit to the continuum and line data. The stellar parameters and UV input spectrum for the modelling were determined through detailed analysis in the UV, optical and near-infrared.

We obtain reasonable fits to the observed photometry, ISO-SWS spectrum, far-IR line fluxes and millimetre CO profiles for a variety of discs, and note that parameter degeneracies preclude the precise derivation of the disc properties. In particular, the effects of dust settling on the vertical thermal structure of the gas in discs strongly influence the line emission, placing further limits on the derivation of a disc gas mass from the line data. We are able to fit the photometry and line data with a well-mixed disc with gas/dust  ~100, but an equally good fit is possible with a disc with gas/dust  ~22, in which dust settling is present. The main advantage of the settled model is that the fit to the observed ISO-SWS spectrum is greatly improved, and it seems probable that some degree of dust settling and gas-depletion is present in the disc. These discs have a radial surface density profile as derived by Hughes et al. (2008), with an exponentially-tapered outer edge.

A model fit to the same data with power-law radial density profile gives a gas/dust ratio of  ~9. This power-law model gives the best fit to the observed SED data, although it does very badly when it comes to fitting the radial intensity distribution of the resolved millimetre continuum. This is due to the flat density power index required to fit the SED and the observed CO 3–2 emission, and presents a challenge for future modelling efforts.

thumbnail Fig. 12

Gas temperature structure for the preferred model (Col. 3 in Table 6), irradiated by an additional X-ray spectral component with LX = 1029.6 erg s-1, as observed by Günther & Schmitt (2009) in HD 163296. This is in comparison to the right hand panel of Fig. 7, in which no X-ray spectral component is present.

The emitted line fluxes are in general sensitive to the degree of dust settling in the disc. This is a firm result in our models, and has serious implications for attempts to derive the disc gas mass and other properties from line observations. This settling flux enhancement arises from changes to the vertical temperature and chemical structure, as settled dust grains allow stellar UV to penetrate deeper into the disc. The effect is strongest in lines which are formed in the warm gas in the inner disc (e.g. a factor  ~30 increase in the H2 S(1) line), but the low excitation molecular lines are also affected, e.g. a factor  ~2 increase in the 13CO 1–0 line.

We explore the effects of the observed UV variability in this object on the gas chemistry in our models, and conclude that this effect is not large enough to affect the observable line fluxes beyond the current range of instrumental uncertainty.

We also examine the effect of X-rays on the gas chemistry of our models, and find that while X-rays present a significant source of gas-heating in the disc surface layers, the observed X-ray luminosity of LX = 1029.6 erg s-1 does not significantly alter the gas chemistry or line emission. Any effects are smaller than those expected as a result of the observed UV variability in this object.

It is difficult to reach any firm conclusions regarding the evolutionary state of the disc of HD 163296. There is some evidence to suggest that the disc is gas-depleted. This is in contrast to the result found for the Herbig Ae star HD 169142 by Meeus et al. (2010), where despite indications that the disc is transitional, it was found to be gas-rich. We note that there are uncertainties associated with our gas/dust values arising from uncertainties in the disc dust composition, and the associated opacity law. This is reflected by the factor of  ~3 spread in the range of derived dust masses found in the literature for this object (Natta et al. 2004; Tannirkulam et al. 2008a; Mannings & Sargent 1997; Isella et al. 2007). All of our derived dust masses are within the range (5−17) × 10-4   M from the literature. We also note the restrictions placed on our conclusions by the assumption of constant dust grain properties throughout the disc, and it has been suggested that the dust properties in discs should in general be variable with radius (Birnstiel et al. 2010; Guilloteau et al. 2011).

It would appear that the modelling of HD 163296 is a far from straightforward task, and it is difficult to fit it into the standard evolutionary picture developed over the past decade. We are unable to find a single disc model which fits perfectly the entire wealth of observational data for this object. As well as the gas/dust ratio, there is uncertainty regarding the disc flaring, which itself is strongly tied in to the disc gas heating and line emission, as well as being indicative of the disc’s evolutionary state. Evidence of in-falling material and the presence of a bipolar outflow in HD 163296 seem to indicate that the star is actively accreting material, typical of a young object with a gas-rich disc. However, it has been suggested that the observed disc geometry and flaring could indicate a star at a later stage of its evolution. This would be consistent with the evidence for grain growth from our modelling, with all the models requiring large grains to fit the observations. There is also possible further evidence for particle growth from high resolution optical spectra of this object (see Appendix). This is also consistent with possible evidence for dust settling and gas-depletion from this study, but is hard to reconcile with evidence for substantial ongoing activity in the inner disc. It is clear that this object is at a fascinating stage in its evolution.

Appendix: Échelle spectra analysis

thumbnail Fig. 13

Transient absorption components in 16 May 1998. The high resolution observed spectrum of several lines (solid lines) is plotted along with the synthetic photospheric spectrum (dashed lines; Kurucz 1993). The short vertical dashes mark the approximate location of the metallic-only transient absorptions at +50 km s-1. A transient red-shifted component of velocity v ~ 50 km s-1 is clearly seen in the metallic but not the hydrogen lines. In addition a blue-shifted absorption is identified in all lines for v   ~   −100 km s-1.

thumbnail Fig. 14

Transient absorption components in 16 July 1998. The high resolution observed spectrum of several lines (solid lines) is plotted along with the photospheric spectrum (dashed lines). No clear transient absorption is detected in the metallic lines. On the other hand a significant blue-shifted absorption is detected in all Balmer lines at v   ~  −100 km s-1.

Section 5 provides possible evidence of a low gas-to-dust ratio in the HD 163296 circumstellar disc as compared to interstellar abundances. This result is also supported in this Appendix from the analysis of visible high resolution échelle spectra.

The EXPORT collaboration (Eiroa et al. 2001) obtained high resolution échelle spectra (R = λλ = 48   000) in the range 3800–5900 Å during five nights in May and July 1998 using the Utrecht Echelle Spectrograph (UES) on the William Herschel Telescope, La Palma.

Significant spectral variability was found in the optical spectrum. The EXPORT data monitored those changes on time scales of hours, days and months (Mendigutía et al. 2011). Almost any strong photospheric line is variable to some extent. In particular, transient absorption components (TACs) are observed in many hydrogen (Balmer series) and metallic lines, such as Ca ii K & H, Na i D2 and D1 and the Fe ii 4924, 5018 and 5169 Å triplet. The Appendix summarizes the main results on the dynamics of the circumstellar gas obtained from the TACS.

thumbnail Fig. 15

Kinematics of the circumstellar gas. Features with similar radial velocity observed in different lines and different dates are assumed to share a common origin and grouped into events. Each event is represented by a box enclosing a number of points. Each point represents the weighted average of the radial velocity simultaneously measured for several lines. The error bar is the standard deviation. The weighted number of lines is given next to each point. Red boxes represent in-falling material, blue boxes out-flowing material and green boxes in-falling material only detected in metallic lines.

Many TACs observed are similar to those studied by Natta et al. (2000) and Mora et al. (2002, 2004) in similar age Herbig Ae stars. These authors propose magnetospheric accretion of primordial gas as the most likely origin. That is, the TACs in Herbig Ae stars typically trace the kinematics of material with gas to dust ratios similar to that of the interstellar medium.

However, a distinctive feature of HD 163296 is the presence of low velocity (v  ~  50  km s-1) red-shifted TACs only observed in refractory elements: Fe i, Fe ii, Ti ii and Sc ii. These “metallic” TACs have no clear counterpart neither in H nor Ca ii and Na i. They were observed during all nights except 28 July 1998. These transient absorptions could be better explained in terms of the model developed by Lagrange-Henri et al. (1988) for the evaporation of solid bodies in the β Pic debris disc and subsequently postulated by Grinin et al. (1991, 1994) for accretion discs around Herbig Ae/Be stars.

Figure 13 shows selected line profiles observed in 16 May 1998. Wavelengths have been converted to radial velocities in the stellar reference system. The approximate location of the metallic-only transient absorptions at +50 km s-1 is displayed with a short vertical dash. The underlying photosphere is represented by a synthetic Kurucz model computed with Teff = 9250 K, log g = 4.07 and [Fe/H] = +0.20, using the ATLAS9 and SYNTHE codes by Kurucz (1993). In this work we have used the GNU Linux version of the codes available online6 (Sbordone et al. 2004). It is clear that no Balmer line absorption counterpart can be identified, even for the higher Hδ and Hζ lines in the series for which no significant core line emission is apparent. On the other hand, blue-shifted absorption with v ~  −100 km s-1 is clearly seen both in metallic and hydrogen lines. The metallic only transient events are not an artifact generated by poor photosphere characterisation, as revealed by Fig. 14, which shows the spectra obtained in 28 July 1998, the only night where no metallic only events where identified.

Similar results were found for WW Vul by Mora et al. (2004), but not for the other four stars studied. In addition, both hydrogen-rich and hydrogen-less TACs are simultaneously detected in HD 163296 and WW Vul. If particle growth indeed explains the metallic-only events, that would suggest the circumstellar disc around HD 163296 (and incidentally WW Vul) has an evolutionary status between a pure primordial gas rich accretion disc and a gas-depleted second generation debris disc. That finding is fully consistent with the results obtained in Sect. 5 from the analysis of Herschel/PACS far-IR spectroscopy and ground-based mm CO observations.

Kinematics of the circumstellar gas

A detailed study of the dynamics of the circumstellar gas around HD 163296 using the method presented in Mora et al. (2004) is out the scope of this paper. The main results are, however, summarized here.

The following lines were studied: the Balmer series (Hβ, Hγ, Hδ, Hϵ and Hζ), The Na i D and Ca ii HK doublets, the Fe ii 4924, 5018 and 5169 Å triplet, Fe i 4046 Å, Ti ii 4444 and 4572 Å and Sc ii 4247 Å. Three types of variability were always simultaneously observed: core line emission (most prominent in the Balmer, Na i and Ca ii lines), transient red-shifted absorptions (associated to in-falling material, observed in all lines) and transient blue-shifted absorptions (associated to outflows, observed in all lines). HD 163296 was more active than the stars studied by Mora et al. (2002, 2004). In particular, the larger line core emission made the analysis more difficult.

The in-falling material always span a large range of velocities, from near zero to +300 km s-1 for the red end of the higher velocity components. The stellar winds maximum velocities can also be high. A maximum of –400 km s-1 was observed in 17 May 1998.

Figure 15 shows a summary of all the transient absorptions identified. Features with similar radial velocity observed in different lines are assumed to share a common origin. The average radial velocity and the standard deviation are displayed as points and error bars in the figure. The number of lines used to compute the average is given next to each point. The fainter lines are given half weight in the average (see Mora et al. 2004 for additional details). Points with similar velocity on nearby nights are also assumed to be produced by the same clump of circumstellar material and are grouped in “events”. Events are displayed in the figure with boxes enclosing a certain number of average points. Red boxes represent in-falling material, blue boxes outflows and green boxes in-falling material only detected in metallic lines. Metallic-only events have smaller radial velocities than hydrogen-rich in-falling material.

One of the most interesting events observed is #4, a strong wind observed in plenty of Balmer and metallic lines during July 1998. The most relevant feature is an almost constant acceleration rate of  ~–0.4 m/s2. The average radial velocity starts at –120 km s-1 on 28 July 1998 and grows up to –235 km s-1 in 31 July 1998. The width of the absorption features has an almost constant value of 120 km s-1, supporting the hypothesis of a common origin of all the transient absorptions in the event.


3

The average IUE spectrum was build with the following individual spectra: SWP 28805, SWP 28811, SWP 28813 and LWR 02065.

5

This CO 3–2 data also seems to indicate the presence of turbulence of order 0.3 km s in the disc of HD 163296 (Hughes et al. 2011). We have not yet explored the effect of this possible stronger turbulence on our results.

References

  1. Abrahamsson, E., Krems, R. V., & Dalgarno, A. 2007, ApJ, 654, 1171 [NASA ADS] [CrossRef] [Google Scholar]
  2. Acke, B., Bouwman, J., Juhász, A., et al. 2010, ApJ, 718, 558 [NASA ADS] [CrossRef] [Google Scholar]
  3. Aresu, G., Kamp, I., Meijerink, R., et al. 2011, A&A, 526, A163 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Barber, R. J., Tennyson, J., Harris, G. J., & Tolchenov, R. N. 2006, MNRAS, 368, 1087 [NASA ADS] [CrossRef] [Google Scholar]
  5. Bell, K. L., Berrington, K. A., & Thomas, M. R. J. 1998, MNRAS, 293, L83 [NASA ADS] [CrossRef] [Google Scholar]
  6. Benisty, M., Natta, A., Isella, A., et al. 2010, A&A, 511, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Bessell, M. S. 1979, PASP, 91, 589 [NASA ADS] [CrossRef] [Google Scholar]
  8. Birnstiel, T., Dullemond, C. P., & Brauer, F. 2010, A&A, 513, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Bouwman, J., de Koter, A., van den Ancker, M. E., & Waters, L. B. F. M. 2000, A&A, 360, 213 [NASA ADS] [Google Scholar]
  10. Brott, I., & Hauschildt, P. H. 2005, in The Three-Dimensional Universe with Gaia, ed. C. Turon, K. S. O’Flaherty, & M. A. C. Perryman, ESA SP, 576, 565 [Google Scholar]
  11. Chambaud, G., Levy, B., Millie, P., Roueff, E., & Tran Minh, F. 1980, J. Phys. B Atomic Mol. Phys., 13, 4205 [NASA ADS] [CrossRef] [Google Scholar]
  12. Cox, A. N., ed. 1999, Allen’s Astrophysical Quantities, 4th ed. (New York: Springer-Verlag) [Google Scholar]
  13. Cushing, M. C., Vacca, W. D., & Rayner, J. T. 2004, PASP, 116, 362 [NASA ADS] [CrossRef] [Google Scholar]
  14. D’Alessio, P., Canto, J., Calvet, N., & Lizano, S. 1998, ApJ, 500, 411 [NASA ADS] [CrossRef] [Google Scholar]
  15. de Winter, D., van den Ancker, M. E., Maira, A., et al. 2001, A&A, 380, 609 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. Devine, D., Grady, C. A., Kimble, R. A., et al. 2000, ApJ, 542, L115 [NASA ADS] [CrossRef] [Google Scholar]
  17. Dorschner, J., Begemann, B., Henning, T., Jaeger, C., & Mutschke, H. 1995, A&A, 300, 503 [NASA ADS] [Google Scholar]
  18. Doucet, C., Pantin, E., Lagage, P. O., & Dullemond, C. P. 2006, A&A, 460, 117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Draine, B. T., & Bertoldi, F. 1996, ApJ, 468, 269 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Dubernet, M.-L., & Grosjean, A. 2002, A&A, 390, 793 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Eiroa, C., Garzón, F., Alberdi, A., et al. 2001, A&A, 365, 110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Eisner, J. A., Graham, J. R., Akeson, R. L., & Najita, J. 2009, ApJ, 692, 309 [NASA ADS] [CrossRef] [Google Scholar]
  23. Faure, A., Crimier, N., Ceccarelli, C., et al. 2007, A&A, 472, 1029 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Fitzpatrick, E. L. 1999, PASP, 111, 63 [NASA ADS] [CrossRef] [Google Scholar]
  25. Flower, D. R. 2001, J. Phys. B Atomic Mol. Phys., 34, 2731 [Google Scholar]
  26. Flower, D. R., & Launay, J. M. 1977, J. Phys. B Atomic Mol. Phys., 10, 3673 [Google Scholar]
  27. Geers, V. C., Augereau, J., Pontoppidan, K. M., et al. 2006, A&A, 459, 545 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Grady, C. A., Devine, D., Woodgate, B., et al. 2000, ApJ, 544, 895 [Google Scholar]
  29. Green, S., Maluendes, S., & McLean, A. D. 1993, ApJS, 85, 181 [NASA ADS] [CrossRef] [Google Scholar]
  30. Grinin, V. P., Kiselev, N. N., Chernova, G. P., Minikulov, N. K., & Voshchinnikov, N. V. 1991, Ap&SS, 186, 283 [NASA ADS] [CrossRef] [Google Scholar]
  31. Grinin, V. P., The, P. S., de Winter, D., et al. 1994, A&A, 292, 165 [NASA ADS] [Google Scholar]
  32. Guilloteau, S., Dutrey, A., Piétu, V., & Boehler, Y. 2011, A&A, 529, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Günther, H. M., & Schmitt, J. H. M. M. 2009, A&A, 494, 1041 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Henning, T., Begemann, B., Mutschke, H., & Dorschner, J. 1995, A&AS, 112, 143 [NASA ADS] [Google Scholar]
  35. Hughes, A. M., Wilner, D. J., Qi, C., & Hogerheijde, M. R. 2008, ApJ, 678, 1119 [NASA ADS] [CrossRef] [Google Scholar]
  36. Hughes, A. M., Wilner, D. J., Andrews, S. M., Qi, C., & Hogerheijde, M. R. 2011, ApJ, 727, 85 [NASA ADS] [CrossRef] [Google Scholar]
  37. Isella, A., Testi, L., Natta, A., et al. 2007, A&A, 469, 213 [CrossRef] [EDP Sciences] [Google Scholar]
  38. Jager, C., Mutschke, H., & Henning, T. 1998, A&A, 332, 291 [NASA ADS] [Google Scholar]
  39. Jankowski, P., & Szalewicz, K. 2005, J. Chem. Phys., 123, 104301 [NASA ADS] [CrossRef] [Google Scholar]
  40. Jaquet, R., Staemmler, V., Smith, M. D., & Flower, D. R. 1992, J. Phys. B Atomic Mol. Phys., 25, 285 [Google Scholar]
  41. Jonkheid, B., Dullemond, C. P., Hogerheijde, M. R., & van Dishoeck, E. F. 2007, A&A, 463, 203 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Kamp, I., & Bertoldi, F. 2000, A&A, 353, 276 [NASA ADS] [Google Scholar]
  43. Kamp, I., Tilling, I., Woitke, P., Thi, W., & Hogerheijde, M. 2010, A&A, 510, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Kamp, I., Woitke, P., Pinte, C., et al. 2011, A&A, 532, A85 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  45. Kurucz, R. L. 1993, VizieR Online Data Catalog, 6039, [Google Scholar]
  46. Lagrange-Henri, A. M., Vidal-Madjar, A., & Ferlet, R. 1988, A&A, 190, 275 [NASA ADS] [Google Scholar]
  47. Launay, J.-M., & Roueff, E. 1977, J. Phys. B Atomic Mol. Phys., 10, 879 [Google Scholar]
  48. Mannings, V. 1994, MNRAS, 271, 587 [NASA ADS] [Google Scholar]
  49. Mannings, V., & Sargent, A. I. 1997, ApJ, 490, 792 [NASA ADS] [CrossRef] [Google Scholar]
  50. Martin-Zaïdi, C., Augereau, J., Ménard, F., et al. 2010, A&A, 516, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Mathews, G. S., Dent, W. R. F., Williams, J. P., et al. 2010, A&A, 518, L127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Meeus, G., Waters, L. B. F. M., Bouwman, J., et al. 2001, A&A, 365, 476 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Meeus, G., Pinte, C., Woitke, P., et al. 2010, A&A, 518, L124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Meijer, J., Dominik, C., de Koter, A., et al. 2008, A&A, 492, 451 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Meijerink, R., Pontoppidan, K. M., Blake, G. A., Poelman, D. R., & Dullemond, C. P. 2009, ApJ, 704, 1471 [NASA ADS] [CrossRef] [Google Scholar]
  56. Mendigutía, I., Eiroa, C., Montesinos, B., et al. 2011, A&A, 529, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Monnier, J. D., Berger, J., Millan-Gabet, R., et al. 2006, ApJ, 647, 444 [NASA ADS] [CrossRef] [Google Scholar]
  58. Montesinos, B., Eiroa, C., Mora, A., & Merín, B. 2009, A&A, 495, 901 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Mora, A., Natta, A., Eiroa, C., et al. 2002, A&A, 393, 259 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Mora, A., Eiroa, C., Natta, A., et al. 2004, A&A, 419, 225 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  61. Natta, A., Grinin, V. P., & Tambovtseva, L. V. 2000, ApJ, 542, 421 [NASA ADS] [CrossRef] [Google Scholar]
  62. Natta, A., Testi, L., Neri, R., Shepherd, D. S., & Wilner, D. J. 2004, A&A, 416, 179 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  63. Oudmaijer, R. D., Palacios, J., Eiroa, C., et al. 2001, A&A, 379, 564 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  64. Panić, O., & Hogerheijde, M. R. 2009, A&A, 508, 707 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Piétu, V., Guilloteau, S., & Dutrey, A. 2005, A&A, 443, 945 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Piétu, V., Dutrey, A., & Guilloteau, S. 2007, A&A, 467, 163 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Pinte, C., Ménard, F., Duchêne, G., & Bastien, P. 2006, A&A, 459, 797 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Pinte, C., Padgett, D. L., Ménard, F., et al. 2008, A&A, 489, 633 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  69. Posch, T., Kerschbaum, F., Fabian, D., et al. 2003, ApJS, 149, 437 [NASA ADS] [CrossRef] [Google Scholar]
  70. Rayner, J. T., Toomey, D. W., Onaka, P. M., et al. 2003, PASP, 115, 362 [NASA ADS] [CrossRef] [Google Scholar]
  71. Renard, S., Malbet, F., Benisty, M., Thiébaut, E., & Berger, J. 2010, A&A, 519, A26 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Sandell, G., Weintraub, D. A., & Hamidouche, M. 2011, ApJ, 727, 26 [NASA ADS] [CrossRef] [Google Scholar]
  73. Sbordone, L., Bonifacio, P., Castelli, F., & Kurucz, R. L. 2004, Mem. Soc. Astron. Ital. Suppl., 5, 93 [Google Scholar]
  74. Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F., & Black, J. H. 2005, A&A, 432, 369 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Servoin, J. L., & Piriou, B. 1973, Physica Status Solidi (b), 55, 677 [Google Scholar]
  76. Sitko, M. L., Carpenter, W. J., Kimes, R. L., et al. 2008, ApJ, 678, 1070 [NASA ADS] [CrossRef] [Google Scholar]
  77. Tannirkulam, A., Monnier, J. D., Harries, T. J., et al. 2008a, ApJ, 689, 513 [NASA ADS] [CrossRef] [Google Scholar]
  78. Tannirkulam, A., Monnier, J. D., Millan-Gabet, R., et al. 2008b, ApJ, 677, L51 [NASA ADS] [CrossRef] [Google Scholar]
  79. Tennyson, J., Zobov, N. F., Williamson, R., Polyansky, O. L., & Bernath, P. F. 2001, J. Phys. Chem. Ref. Data, 30, 735 [NASA ADS] [CrossRef] [Google Scholar]
  80. Thi, W. F., van Dishoeck, E. F., Blake, G. A., et al. 2001, ApJ, 561, 1074 [NASA ADS] [CrossRef] [Google Scholar]
  81. Tielens, A. G. G. M. 2008, ARA&A, 46, 289 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  82. Vacca, W. D., Cushing, M. C., & Rayner, J. T. 2003, PASP, 115, 389 [NASA ADS] [CrossRef] [Google Scholar]
  83. van den Ancker, M. E., Bouwman, J., Wesselius, P. R., et al. 2000, A&A, 357, 325 [NASA ADS] [Google Scholar]
  84. van Leeuwen, F., ed. 2007, Hipparcos, the New Reduction of the Raw Data, Astrophys. Space Sci. Lib., 350 [Google Scholar]
  85. Warren, S. G., & Brandt, R. E. 2008, J. Geophys. Res., 113, D14220 [Google Scholar]
  86. Wassell, E. J., Grady, C. A., Woodgate, B., Kimble, R. A., & Bruhweiler, F. C. 2006, ApJ, 650, 985 [NASA ADS] [CrossRef] [Google Scholar]
  87. Wernli, M., Valiron, P., Faure, A., et al. 2006, A&A, 446, 367 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  88. Wilson, N. J., & Bell, K. L. 2002, MNRAS, 337, 1027 [NASA ADS] [CrossRef] [Google Scholar]
  89. Wisniewski, J. P., Clampin, M., Grady, C. A., et al. 2008, ApJ, 682, 548 [NASA ADS] [CrossRef] [Google Scholar]
  90. Woitke, P., Kamp, I., & Thi, W. 2009, A&A, 501, 383 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. Woitke, P., Pinte, C., Tilling, I., et al. 2010, MNRAS, 405, L26 [NASA ADS] [CrossRef] [Google Scholar]
  92. Woitke, P., Riaz, B., Duchêne, G., et al. 2011, A&A, 534, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  93. Woodall, J., Agúndez, M., Markwick-Kemper, A. J., & Millar, T. J. 2007, A&A, 466, 1197 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  94. Wu, C., Ake, T. B., Boggess, A., et al. 1983, NASA IUE Newsl., No. 22, 2+324 pp., 22 [Google Scholar]
  95. Wu, C., Reichert, G. A., Ake, T. B., et al. 1992, NASA Reference Publication, 1285 [Google Scholar]
  96. Yang, B., Stancil, P. C., Balakrishnan, N., & Forrey, R. C. 2010, ApJ, 718, 1062 [NASA ADS] [CrossRef] [Google Scholar]
  97. Yi, S., Demarque, P., Kim, Y.-C., et al. 2001, ApJS, 136, 417 [NASA ADS] [CrossRef] [Google Scholar]
  98. Zuckerman, B., Forveille, T., & Kastner, J. H. 1995, Nature, 373, 494 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]

All Tables

Table 1

Observed line data with Herschel/PACS.

Table 2

Photometric fluxes observed with Herschel/PACS.

Table 3

Parameters for HD 163296.

Table 4

Photometry for HD 163296.

Table 5

Fixed disc model parameters.

Table 6

Model parameters for four well-fitting disc models.

Table 7

Dust grain composition.

Table 8

Integrated line fluxes for the transitions observed by PACS, and additional CO and H2 lines, as predicted by the models listed in Table 6.

Table 9

Effects of UV variability. Comparison between the preferred model, and an identical model illuminated by an input spectrum typical of a “high UV” state for this object.

All Figures

thumbnail Fig. 1

Observed [Oi]63 μm emission line with Herschel/PACS.

In the text
thumbnail Fig. 2

[Oi]63 μm spectra observed by Herschel/PACS towards HD 163296. Each spectrum corresponds to a 9.4 × 9.4 arcsec pixel centered at the given coordinates. Emission is dominated by the central pixel, corresponding to a maximum outer radius of  ~ 560 AU, assuming a distance of 118.6 parsecs to this object.

In the text
thumbnail Fig. 3

Spectral energy distribution for HD 163296. The green dots represent the fluxes corresponding to the UBVRIJHK photometry. The size of the dots is of the order or larger than the uncertainties. The IUE average spectrum and the FUSE+STIS spectrum are plotted in blue and red, respectively; the latter has been slightly smoothed to reduce the noise. The SpeX/IRTF spectrum is plotted in magenta. The black solid line is a photospheric model computed for the specific stellar parameters given in Table 3, reddened with E(B − V) =  + 0.15 and normalized at the flux in V. See text for details.

In the text
thumbnail Fig. 4

Spectral energy distribution for the best-fit fully mixed disc model (“preferred” model), obtained from a simultaneous fit to the observed SED, ISO-SWS spectrum and line emission from HD 163296 (solid black line). Black dotted line indicates the SED for the model with power-law density profile, which is unable to fit the mm continuum image for HD 163296. Blue circles indicate (with increasing wavelength) simultaneous UBVRIJHK photometry (Eiroa et al. 2001; Oudmaijer et al. 2001), LM photometry (de Winter et al. 2001), IRAS photometry (12–100 mic), sub-mm photometry (Mannings 1994), and millimetre photometry (Isella et al. 2007). Also marked are PACS photometric observations (red circles), PACS continua derived from the spectroscopic observations (black circles), SCUBA photometry (green circles) (Sandell et al. 2011), scaled VLT/UVES spectrum (blue line, Martin-Zaidi, in prep.) and the ISO-SWS spectrum (green line). The red line shows the stellar+UV input spectrum. Downwards arrow denotes upper limit. All fluxes were corrected for interstellar reddening using the Fitzpatrick parameterisation (Fitzpatrick 1999) with RV = 3.1 and E(B − V) = 0.15.

In the text
thumbnail Fig. 5

1.3 mm continuum emission maps for HD 163296.From left: observed map; emission map for preferred model (Col. 3 in Table 6); residuals for preferred model; emission map for power-law model (Col. 4 in Table 6); residuals for power-law model. Residuals computed by subtracting the model intensity from the observed intensity. Contours spaced at 12 mJy intervals, corresponding to [3,6,9,12,15,18,21,24,27,30,33,36,39,42,45] ×    σ.

In the text
thumbnail Fig. 6

Plotted in black are the line profiles observed by Isella et al. (2007) for the 12CO J = 3–2 (upper panel), 12CO J = 2–1 (middle panel) and 13CO J = 1–0 (lower panel) transitions, with the corresponding profiles from the preferred disc model in red. This refers to a single simultaneous fit to the observed continuum and line data (Col. 3 in Table 6). Observed profiles are obtained by integrating over the whole disc.

In the text
thumbnail Fig. 7

The total hydrogen number density (left panel) and gas temperature (right panel) are plotted for a vertical cross-section through the preferred disc model (Col. 3 in Table 6). Dashed lines indicate contours of gas temperature and visual extinction AV, as marked.

In the text
thumbnail Fig. 8

Spatial origin of the various gas emission lines in the preferred model (Col. 3 in Table 6). From top-left clockwise: CO J = 3−2, CO J = 18−17, [Cii] 157.74 μm, o-H2 S(1) (with red Tg contours), [Oi] 63.18 μm, 13CO J = 1–0. In each panel, the upper plot shows the line optical depth as a function of radius (blue line) and the continuum optical depth at the corresponding wavelength (black line). The middle plot shows the cumulative contribution to the total line flux with increasing radius. The lower plot shows the gas species number density, and the black line marks the cells that contribute the most to the line flux in their vertical column. The two vertical dashed lines indicate 15% and 85% of the radially cumulative face-on line flux respectively, i.e. 70% of the line flux originates from within the two dashed lines.

In the text
thumbnail Fig. 9

CO abundances in the preferred model (Col. 3 in Table 6). Left panel: gas-phase CO abundances. White and blue dashed lines plot the gas and dust temperature contours respectively. Right panel: CO ice abundance.

In the text
thumbnail Fig. 10

The fractional enhancement in line flux is plotted against the minimum grain size affected by dust settling, relative to the preferred (fully mixed) model (as = amax). All models have a settling parameter of 0.5. Top panel: CO lines; middle panel: water lines; lower panel: [Oi]63 μm, [Cii]158 μm, OH 79.11 μm and H2 S(1) lines.

In the text
thumbnail Fig. 11

Effect of dust settling on the disc chemical abundances. Left column indicates abundances in a fully mixed model, and right column represents a strongly-settled model with otherwise identical parameters. Top row: CO abundance with dust (blue dashed lines) and gas (white dashed line) temperature contours. Middle row: H2 abundance with 300 K gas temperature contour (black dashed line) and visual extinction contour (white dashed line). Bottom row: H2O abundance. White dashed lines indicate Tgas contours enclosing “hot water” region, red dashed lines indicate contours of UV field strength per hydrogen nucleus, log (χ/nH) as defined in Draine & Bertoldi (1996), enclosing the cool water belt.

In the text
thumbnail Fig. 12

Gas temperature structure for the preferred model (Col. 3 in Table 6), irradiated by an additional X-ray spectral component with LX = 1029.6 erg s-1, as observed by Günther & Schmitt (2009) in HD 163296. This is in comparison to the right hand panel of Fig. 7, in which no X-ray spectral component is present.

In the text
thumbnail Fig. 13

Transient absorption components in 16 May 1998. The high resolution observed spectrum of several lines (solid lines) is plotted along with the synthetic photospheric spectrum (dashed lines; Kurucz 1993). The short vertical dashes mark the approximate location of the metallic-only transient absorptions at +50 km s-1. A transient red-shifted component of velocity v ~ 50 km s-1 is clearly seen in the metallic but not the hydrogen lines. In addition a blue-shifted absorption is identified in all lines for v   ~   −100 km s-1.

In the text
thumbnail Fig. 14

Transient absorption components in 16 July 1998. The high resolution observed spectrum of several lines (solid lines) is plotted along with the photospheric spectrum (dashed lines). No clear transient absorption is detected in the metallic lines. On the other hand a significant blue-shifted absorption is detected in all Balmer lines at v   ~  −100 km s-1.

In the text
thumbnail Fig. 15

Kinematics of the circumstellar gas. Features with similar radial velocity observed in different lines and different dates are assumed to share a common origin and grouped into events. Each event is represented by a box enclosing a number of points. Each point represents the weighted average of the radial velocity simultaneously measured for several lines. The error bar is the standard deviation. The weighted number of lines is given next to each point. Red boxes represent in-falling material, blue boxes out-flowing material and green boxes in-falling material only detected in metallic lines.

In the text

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

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

Initial download of the metrics may take a while.