© ESO, 2011

1. Introduction

A significant fraction of white dwarfs (WD) that have cooled below T_eff ≈ 25   000 K (20−30%) displays photospheric absorption lines from metals (Zuckerman et al. 2010). These polluted WDs must actively accrete matter at rates from 10^-15 to 10^-17   M_⊙/yr, otherwise the atmospheres would have been purified by gravitational settling of heavy elements (Paquette et al. 1986; Koester & Wilken 2006; Koester 2009). Until the recent past it was assumed that these stars accrete matter from the interstellar medium (Dupuis et al. 1992). This view has changed considerably in the past few years because it was discovered that many of the polluted WDs host dust discs within the stellar tidal radius (Zuckerman & Becklin 1987; Becklin et al. 2005; Kilic et al. 2005, 2006; Farihi et al. 2007, 2009; Farihi et al. 2010). It is now widely accepted that WDs are polluted by matter accretion from these discs. It is thought that they contain material from tidally disrupted asteroids that were scattered towards the central stars as a consequence of dynamical resettling of a planetary system in the post-main sequence phase (Debes & Sigurdsson 2002; Jura 2003).

The photospheric metal abundance pattern in the polluted WDs allows us to indirectly measure the composition of the accreted matter. This opens up the exciting possibility of studying the composition of extrasolar planetary material. The first impressive results have already been obtained, showing that the orbiting debris broadly mimics the terrestrial matter of the inner Solar System, including the possibility of water (Zuckerman et al. 2007; Klein et al. 2010; Dufour et al. 2010; Farihi et al. 2011). These results are based on our knowledge of metal diffusion rates in WD atmospheres and envelopes, which are difficult to obtain (Koester 2009). Several uncertainties can affect the resulting composition of accreted material as concluded from the photospheric abundance pattern. One example is the assumption of stationary accretion on several diffusion timescales of all elements involved. It is therefore highly desirable to exploit alternative possibilities of determining the chemical composition of the accreted material.

Such an alternative method is offered by the recent discovery that three circumstellar dust discs around polluted WDs also host gaseous metal components that are interpreted as the collisional remains of solid material (Gänsicke et al. 2006, 2007, 2008). We are developing accretion disc models to derive the chemical abundances in the gas discs, so that we should be able to directly measure the composition of the parent planetary material. This method is complementary to the measurement of photospheric abundances because the disc spectra might reveal trace elements that are not seen in the WD spectra. It also provides a means to test our understanding of diffusion processes in the WD atmospheres and envelopes. Building on our preliminary work (Werner et al. 2009), we present here new results for our modelling efforts of gas discs.

In Sect. 2, we briefly introduce our method, followed by a description of the object that we study in detail (Sect. 3). In Sect. 4, we first summarise our results concerning simple models for pure-calcium discs, investigating the influence of effective temperature and surface-mass density on the emergent spectrum. We then present vertical structures and spectra of discs composed of an asteroid-like mixture of light metals, comprising C, O, Mg, Si, and Ca. Finally, we present results for non-axisymmetric disc geometries in order to explain the observed asymmetry of the double-peaked line profiles. Based on hydrodynamical simulations, we investigate the time evolution of this asymmetry in comparison with observations.

2. Accretion-disc modelling

For calculating geometrically thin accretion-disc models, we use our code AcDc (Nagel et al. 2004). We assume axial symmetry, so that we can separate the disc into concentric rings of plane-parallel geometry. In that way, the radiative transfer becomes a one-dimensional problem. By integrating the spectra of the individual rings, we obtain a complete disc spectrum for different inclination angles.

The free parameters of one ring with radius R are effective temperature T_eff(R), surface mass density Σ(R), chemical composition, and the WD mass M_WD. For the energy equation we assume that the emitted radiation is viscously generated, so the Reynolds number Re (or α) enters as an additional parameter. In the case of viscous α-discs, the radial run of T_eff(R) and Σ(R) can be expressed in terms of the mass-accretion rate and mass M_WD and radius R_WD of the central star (Shakura & Sunyaev 1973). For comparison with observations, the emergent spectra from ring segments are Doppler-shifted to account for Keplerian motion, hence R_WD and disc inclination i appear as additional parameters.

For each disc ring, the following set of coupled equations were solved simultaneously under the constraints of particle number and charge conservation:

• radiation transfer for the specific intensity I at frequency ν(1)with the absorption coefficient χ, the emission coefficient η, the disc height z above the midplane, and μ = cosθ, with θ the angle between the ray and z;

• hydrostatic equilibrium of gravitation, gas pressure P_gas, and radiation pressure (2)with ρ denoting the mass density and H the Eddington flux. Here, we also introduced the column-mass density m as (3)

• energy balance between the viscously generated energy E_mech and the radiative energy loss E_rad(4)with (5)and (6)with the angular velocity ω, the mean intensity J, and w the kinematic viscosity written following Lynden-Bell & Pringle (1974): (7)For the models presented here we assume Re = 15   000;

• NLTE rate equations for the population numbers n_i of the atomic levels i(8)where P_ij denotes the rate coefficients, consisting of radiative and electron collisional components.

Detailed information about the involved atomic data is provided in the form of a model atom (cf. Rauch & Deetjen 2003). The model atoms we used for our NLTE calculations are summarised in Table 1. They were taken from TMAD, the Tübingen Model Atom Database¹.

The principal problem for any modelling attempt is posed by the question of what heats the Ca ii emission line region. It cannot be gravitational energy released through viscosity because the required mass-accretion rate would have to be of the order of 10^-8   M_⊙/yr, which is by many orders of magnitude larger than the accretion rate invoked for the presence of settling metals in DAZ photospheres (≈10^-15   M_⊙/yr, Koester & Wilken 2006). A speculation by Jura (2008) was additional heating by energy dissipation through disc asymmetries, which are driven by some external unseen planet. Alternatively, Melis et al. (2010) suggested a “Z ii” model in analogy to H ii regions. In the case of the discs, the metal-dominated material is photoionised and heated (hence the name Z ii) by absorbing photons from the WD and cools through optically thick emission lines. Given this lack of knowledge, we need to use T_eff, which is a measure of the vertically integrated dissipated energy (Eq. (5)), hence the amount of energy radiated away from the disc surface per unit time and area, as a free parameter.

3. SDSS J122859.93+104032.9

Our models are tailored to SDSS J122859.93+104032.9 (henceforth SDSS J1228+1040). This metal-polluted WD was the first one discovered to be surrounded by a gaseous metal disc (Gänsicke et al. 2006). It is a DAZ white dwarf with atmospheric parameters T_eff = 22   020  ±  200 K and log g = 8.24  ±  0.04, and the derived stellar mass and radius are M_WD = 0.77  ±  0.02   M_⊙ and R_WD = 0.0111  ±  0.0003   R_⊙.

The Ca ii infrared triplet (λλ 8498, 8542, 8662 Å) with double-peak emission line profiles is the hallmark of the gaseous metal discs (Fig. 1). In the case of SDSS J1228+1040, Gänsicke et al. (2006) measured a peak-to-peak separation of 630 km s^-1, i.e. the Keplerian rotation velocity is v   sini = 315 km s^-1. From a spectral analysis with a kinematical LTE emission model, it was concluded that we see a geometrically thin, optically thick gaseous disc at high inclination (i = 70°). Two other weak emission features of Fe ii λλ 5018, 5169 Å were seen by Gänsicke et al. (2006). Subsequent observations by Melis et al. (2010) failed to detect the Fe ii λ 5018 Å line, possibly because of the lower signal-to-noise ratio of their high-resolution spectra as compared to the low-resolution spectrum of Gänsicke et al. (2006). On the other hand, very weak emissions from the Ca ii H and K lines were discovered by Melis et al. (2010).

There is a clear asymmetry in the emission strengths of the double-peak line profiles of SDSS J1228+1040 (Gänsicke et al. 2006). A similar phenomenon is well known from Be star discs (Carciofi 2010) and is ascribed to one-armed spiral waves. In addition, Melis et al. (2010) observe that the asymmetry in SDSS J1228+1040 had changed such that the stronger of the two emission peaks has switched from the red side of the double-peaked emission complex, as seen in Gänsicke et al. (2006), to the blue side. We describe these characteristics in more detail in Sect. 4.3, where we investigate the temporal variability predicted by our models.

4. Results

4.1. Parameter study for pure calcium discs

In a first exploratory study (Werner et al. 2009), we calculated disc models composed only of calcium, with two values for Σ and three values for T_eff. The inner and outer disc radii were set to 1.0 R_⊙ and 1.2 R_⊙, respectively. We found that the emission strength of the Ca ii triplet decreases with T_eff increasing from 5000 K to 7000 K (Fig. 1). The reason is the shifting Ca ii/Ca iii ionisation balance. A closer comparison of the three line components shows that their emission strengths become equal with increasing T_eff, a behaviour that constrains T_eff. A similar trend is seen when Σ is reduced from 0.6 g/cm² to 0.3 g/cm² at T_eff = 6500 K. We stress that the models have a considerable continuum flux compared to the line-emission peak heights. The relative strength of the profile depression between the double-peaks increases with increasing inclination, the double peak structures become broader. The disc models are optically thin in terms of the Rosseland optical depth.

Fig. 1 Normalised spectra of three pure Ca disc models (top) with different Teff compared to the observed spectrum of SDSS 1228+1040 taken from SDSS (bottom). The model spectra are shifted vertically for clarity.

The observed spectrum does not show H_α emission. This can be used to determine an upper limit for the hydrogen abundance. We varied the H content (H = 1%, 0.1%, 0.01%, by mass) and found that with an abundance of 1%, the H_α peak height is comparable to that of the Ca ii triplet, so would be detectable in the spectrum of SDSS J1228+1040.

The effective temperature of the disc is well constrained by three models, with T_eff  ≈  5800 K (Fig. 1), but the asymmetry of the line profiles is of course not matched by our symmetric models. The cooler model (T_eff = 5000 K) is perhaps more favourable because of the larger line-to-continuum emission ratio, while the hotter model (T_eff = 7000 K) has the advantage that the relative strengths of the three line components are reproduced better.

4.2. O-Si-Mg-C-Ca-H disc models

In the next step, we expanded the set of considered chemical elements in order to achieve a composition comparable to CI chondrites in the Solar System or a bulk-Earth mixture (Klein et al. 2010). New species included are O, Si, Mg, and C. Iron poses special numerical problems and will be introduced in future work. The initially chosen element abundances are representative of the class of CI chondrites. In detail, they are: H = 10^-8, C = 4.6, O = 65.5, Mg = 13.5, Si = 15.1, Ca = 1.3 (% mass fraction).

We investigated the influence of the radial disc extent on the spectrum. At its largest, the disc model consists of 11 rings extending from R_i = 2 R_WD = 0.022 R_⊙ to an outer radius of R_o = 136 R_WD = 1.5 R_⊙. The ring radii and effective temperatures are listed in Table 2. The model for ring 8 did not converge; its spectrum was set equal to that of ring 9. The entire disc is assumed to have a radially constant surface density of 0.3 g/cm².

In Fig. 2 we present the vertical run of temperature, Rosseland optical depth τ_ross, mass density, gravity, and geometrical height of the disc rings. Only the inner rings are optically thick. Figure 3 shows spectra of three different rings, at 2.0, 71.4, and 136.4 R_WD. They all show emission lines, getting stronger in the outer parts of the disc.

Fig. 2 Vertical structure of every other of the 11 disc rings along the column-mass density (increasing from outer layers towards the disc’s midplane).

We combined the 11 rings to several axisymmetric discs with varying radial extent. Their spectra are shown in Fig. 4. The dip between the triplet line components of the Ca ii lines seen in observed spectra can best be reproduced if the disc does not reach too far inwards, resulting in a minimum inner radius of about 60 R_WD. Variation in the outer disc radius has almost no effect on the Ca ii triplet, but reducing the outer radius decreases the line strength of a C ii doublet (λλ 8685, 8699 Å) near the red component of the Ca ii triplet. While the lack of observed C ii emission may suggest a constraint on the outer disc radius, we note that decreasing the carbon abundance by an order of magnitude also removes this feature from the model spectra (see Fig. 5). Given the carbon-poor nature of the observed discs (Reach et al. 2005; Jura et al. 2009) and polluted WD atmospheres (Jura 2006), abundance may play a substantial role.

Our model strongly overpredicts the strength of the Ca ii H & K emission lines. This problem has already been noted in the pure-Ca disc described above (Werner et al. 2009). This problem did not disappear in the present multi-element disc model. In contrast, our model features emission lines from other metals, markedly from Mg ii that are not observed. Another weakness of our model is the relatively strong continuum flux that would detectable by distorting the WD spectrum. This problem is alleviated when the effective temperature is reduced.

4.3. Asymmetric disc models

To investigate the asymmetry of the line profiles in the spectra of SDSS J1228+1040, we modified our method combining the disc rings to receive a complete disc spectrum. In the surface integrating step, we used only parts of the rings in order to construct a spiral-arm like or an eccentric shape of the disc. The ring segments are still assumed to undergo Keplerian rotation for the calculation of the spectral Doppler shift. Depending on the orientation towards the observer, these non-axisymmetric accretion discs result in asymmetric line profiles that can be compared to the observations.

Fig. 3 Spectra of three (non-rotating) circumstellar rings at i = 77° and: a) 2.0; b) 71.4; and c) 136.4 RWD. For clarity, the spectra are vertically shifted by amounts as colour-coded numerically on the plots. With the exception of the infrared Ca ii triplet and the adjacent C ii doublet, as well as Ca ii H & K, no fine structure splitting was taken into account.

We performed hydrodynamic calculations with the FARGO code (Masset 2000) in order to motivate our particular choice of geometries. This code was originally developed to compute the hydrodynamic evolution of a protoplanetary disc on a fixed polar coordinate system, but it is also suitable for other kinds of sheared fluid discs. A possible scenario could be that an asteroid coming within the tidal radius of the WD is disrupted, forming a locally concentrated debris cloud. Matter then spirals inwards. We start the simulations by putting a gas blob at the tidal radius (R = 137   R_WD) onto a circular orbit (P_orb = 3.6  ×  10⁴ s) around the WD. The gas mass is 7  ×  10²¹ g. The initial surface density distribution is Gaussian with a blob radius of 2.8  ×  10⁵ km. We chose an open boundary condition and fixed the simulation rim at a value of R = 194   R_WD.

The Keplerian rotating material gets smeared out into a spiral-arm like structure within a short time of t = 1.21  ×  10⁴ s. After another Δt = 3.1  ×  10⁴ s, the spiral is catching up with its own starting point and an eccentric closed disc forms and retains this shape for the rest of the simulation. For both situations the surface mass density is rather homogeneous with Σ  ≈  0.3 g/cm². Examples for both geometries found by the simulations are shown in Fig. 6. In Fig. 7 we display the corresponding assembly of our ring segment models.

The upper panel of Fig. 8 shows the resulting spectrum in the case of a spiral-arm like shape with a minimum radius of 58 R_WD and a maximum radius of 136 R_WD for an observer’s position according to Fig. 7 (left). The strength of the Ca ii lines and the asymmetry of the blue and red parts of each line component are well reproduced. Only the dip between the Doppler-shifted parts of each line are not as deep as in the observation. The mass of such a spiral arm would be about 3.9  ×  10²¹ g. In the lower panel of Fig. 8, the resulting spectrum of an eccentric disc for an observer’s position according to Fig. 7 (right) is shown. The fit quality to the observed spectrum is similar to the spiral-arm case. Such an eccentric disc would have a mass of about 6.3  ×  10²¹ g.

Fig. 4 Effect of variation in the disc’s radial extent at i  =  77°. Upper panel: the inner radius is fixed at 58 RWD and the outer radius varies: 84, 110, and 136 RWD. The C ii line disappears with decreasing outer radius. Lower panel: the outer radius is fixed at 110 RWD and the inner radius varies: 2, 30, and 58 RWD. Broad line wings develop when the inner radius becomes smaller.

Fig. 5 Variation in the carbon abundance. Shown are two models (a): 4.6% and b): 0.46% by mass) and the observed spectrum (red). The radial disc extent is 58−136 RWD, i  =  77°.

Fig. 6 Two exemplary structures from the hydrodynamical simulations. The initial gas blob first evolves into a spiral-arm structure (left) and after some time into an eccentric disc (right). The surface density is colour-coded as indicated.

Fig. 7 Sketch of the two non-axisymmetric geometries used for our spectral models. Left: a spiral-arm like structure. Right: an eccentric accretion disc. Both are composed of disc ring segments (red). For the unmarked ring segments, the emergent flux is set to zero. The observer’s position is towards the right.

Fig. 8 Emergent spectrum for the a) spiral-arm like and b) eccentric disc structures as shown in Fig. 7. The modelled profiles (thick lines) display an asymmetry as observed (thin line). For the spiral arm the radial extent is 58−136 RWD, whereas in the case of the eccentric disc it is 45−136 RWD. For both, i  =  77°, φ = 285°, and Σ = 0.3 g/cm2.

Fig. 9 Time series of the computed Ca ii triplet emission, covering 100 min. For the spiral-arm geometry (left panel), a change in the line-profile asymmetry evolves very fast, whereas for the eccentric disc (right) it takes much longer.

4.4. Time variability

Gänsicke et al. (2006) present two time series of spectra from SDSS J1228+1040 taken on June 30 and July 1, 2006. They indicate that the asymmetry of the Ca ii line profiles did not change during the two 100 min observations. To investigate this with our models, we assumed a Keplerian rotating spiral arm and an eccentric disc with ω = 0.01°/s, and calculated synthetic spectra in steps of Δφ = 6° concerning the orientation towards the observer, starting at φ = 231° and centred on 285°, for which we found the best fit to the emission lines. The resulting spectral time series are shown in Fig. 9, covering 100 min in order to be comparable with the lower right part of Fig. 1 in Gänsicke et al. (2006). In the case of a spiral-arm like geometry, the asymmetry of the line profile would change significantly within this time interval (Fig. 9, left panel). The separation of the maxima of the double peaks decreases, and at the same time the Ca ii line asymmetry reverses in wavelength. In contrast, for the eccentric disc the relative strength of the asymmetry changes quite slowly (Fig. 9, right) and the double-peak separation remains almost constant, which is in better agreement with the non-variable observations. On the other hand, Melis et al. (2010) found a switch of the asymmetry on a longer timescale, between the Gänsicke et al. observations in 2006 and their own in 2008.

5. Summary and conclusion

We performed non-LTE modelling of gaseous metal discs around WDs in order to study their spectral signatures in comparison to observations of SDSS J1228+1040. The modelling was done in three steps.

At first, pure calcium models were constructed to constrain the disc characteristics by fitting the observed infrared Ca ii emission triplet. Qualitatively good fits can be obtained with a geometrically and optically thin, Keplerian viscous gas disc ring at a distance of 1.2 R_⊙ from the WD, with T_eff  ≈  5800 K and a low surface mass density Σ  ≈  0.3 g/cm². The disc is hydrogen-deficient (H  < 1% by mass), and it is located within the tidal disruption radius (R_tidal = 1.5 R_⊙).

In the second step, we constructed axisymmetric disc models composed of elements in Chondritic abundance, namely C, O, Si, Mg, and Ca. We found that the inner radius of the observed, emitting, Ca ii gaseous component of the disc is well constrained at  ≥0.65 R_⊙ = 58 R_wd. The outer radius can be constrained by the emission strength of a C ii doublet (λλ 8685, 8699 Å) that is not seen in the observations. An alternative explanation could be a reduced carbon abundance, which would be a hint that the disc in SDSS J1228+1040 has a bulk-Earth like composition instead of a CI chondritic one.

In the third step, we investigated asymmetric disc structures by assuming spiral-arm and eccentric disc shapes as suggested by hydrodynamical simulations. Both geometries can qualitatively explain the asymmetry observed in the double-peak line profiles well. An investigation of the time variability of the computed line profiles suggests that the eccentric disc model displays less significant variability than the spiral-arm geometry. Considering the current observational material, the eccentric disc model is more realistic. Mass estimates for the circumstellar gas material using the two geometric models results in 3−6  ×  10²¹ g, which are equivalent to the mass of a 135-km diameter Solar System asteroid.

Fig. 10 Mg ii resonance doublet in the HST/COS spectrum, overplotted with our (rotating) disc model and a pure-H WD model, plus the co-added model. Fine-structure splitting of the line is considered in the disc model.

Following the acceptance of our paper, an HST/COS observation of SDSS J1228+1040 became public (dataset LB5Z040, observation date 2010-02-19). This is the first available UV observation. The spectrum shows a strong emission line of Mg ii that our disc model, as displayed in Fig. 3, predicts to be the strongest UV line. This discovery is essential because Mg is the third element (after Ca and Fe) that is observed in the metal disc, supporting the idea of ground planetary material. Figure 10 shows the observed Mg ii resonance doublet (λλ 2796.35, 2803.53 Å) in emission with central absorption components of photospheric origin. Overplotted is a pure-H WD model scaled to fit the continuum and our disc model (including all rings 5−11 as listed in Table 2, i = 77°) scaled arbitrarily, plus the co-added WD+disc model spectrum. The optical depth in the Mg ii line varies across the disc and has a maximum of τ ~ 10⁶ in the outermost disc region.

Acknowledgments

We thank Tobias Müller from the Computational Physics group at our institute for providing us with the FARGO simulations. S.H. and T.R. are supported by DFG (grant We1312/37-1) and DLR (grant 50 OR 0806), respectively.

References

All Tables

All Figures

	Fig. 1 Normalised spectra of three pure Ca disc models (top) with different Teff compared to the observed spectrum of SDSS 1228+1040 taken from SDSS (bottom). The model spectra are shifted vertically for clarity.
In the text

	Fig. 2 Vertical structure of every other of the 11 disc rings along the column-mass density (increasing from outer layers towards the disc’s midplane).
In the text

	Fig. 3 Spectra of three (non-rotating) circumstellar rings at i = 77° and: a) 2.0; b) 71.4; and c) 136.4 RWD. For clarity, the spectra are vertically shifted by amounts as colour-coded numerically on the plots. With the exception of the infrared Ca ii triplet and the adjacent C ii doublet, as well as Ca ii H & K, no fine structure splitting was taken into account.
In the text

	Fig. 4 Effect of variation in the disc’s radial extent at i  =  77°. Upper panel: the inner radius is fixed at 58 RWD and the outer radius varies: 84, 110, and 136 RWD. The C ii line disappears with decreasing outer radius. Lower panel: the outer radius is fixed at 110 RWD and the inner radius varies: 2, 30, and 58 RWD. Broad line wings develop when the inner radius becomes smaller.
In the text

	Fig. 5 Variation in the carbon abundance. Shown are two models (a): 4.6% and b): 0.46% by mass) and the observed spectrum (red). The radial disc extent is 58−136 RWD, i  =  77°.
In the text

	Fig. 6 Two exemplary structures from the hydrodynamical simulations. The initial gas blob first evolves into a spiral-arm structure (left) and after some time into an eccentric disc (right). The surface density is colour-coded as indicated.
In the text

	Fig. 7 Sketch of the two non-axisymmetric geometries used for our spectral models. Left: a spiral-arm like structure. Right: an eccentric accretion disc. Both are composed of disc ring segments (red). For the unmarked ring segments, the emergent flux is set to zero. The observer’s position is towards the right.
In the text

	Fig. 8 Emergent spectrum for the a) spiral-arm like and b) eccentric disc structures as shown in Fig. 7. The modelled profiles (thick lines) display an asymmetry as observed (thin line). For the spiral arm the radial extent is 58−136 RWD, whereas in the case of the eccentric disc it is 45−136 RWD. For both, i  =  77°, φ = 285°, and Σ = 0.3 g/cm2.
In the text

	Fig. 9 Time series of the computed Ca ii triplet emission, covering 100 min. For the spiral-arm geometry (left panel), a change in the line-profile asymmetry evolves very fast, whereas for the eccentric disc (right) it takes much longer.
In the text

	Fig. 10 Mg ii resonance doublet in the HST/COS spectrum, overplotted with our (rotating) disc model and a pure-H WD model, plus the co-added model. Fine-structure splitting of the line is considered in the disc model.
In the text