METAL ACCRETION ONTO WHITE DWARFS CAUSED BY POYNTING–ROBERTSON DRAG ON THEIR DEBRIS DISKS

Recent infrared observations with Spitzer and ground-based facilities (Zuckerman & Becklin 1987; Graham et al. 1990; Farihi et al. 2010) revealed near-infrared excesses around more than a dozen metal-rich white dwarfs (WDs). This emission was interpreted as reprocessing of the WD radiation by refractory material residing in an extended, optically thick and geometrically thin disk (Jura 2003; Jura et al. 2007), similar to the rings of Saturn (Cuzzi et al. 2010).

It was hypothesized by Jura (2003) that such compact disks of high-Z material may be naturally produced by tidal disruption of asteroid-like bodies entering the Roche radius of the WD. This idea naturally explains the well-defined outer radii R_out ∼ 1 R_☉ of disks since the Roche radius R_R for tidal disruption of a self-gravitating object of normal density (ρ ∼ 1 g cm⁻³) by M_⋆ ∼ M_☉ central mass is R_R ∼ (M_⋆/ρ)^1/3 ≈ 1 R_☉ (Jura 2003).

Moreover, Jura (2003) went on to suggest that the high-Z material contained in compact disks is responsible for the observed metal enrichment of a significant fraction of WDs. Accretion of this circumstellar material at high enough rate $\dot{M}_Z$ onto WD could maintain a non-zero abundance of metals in the WD atmosphere against rapid gravitational settling. This scenario thus provides a promising alternative to the previously widely discussed interstellar accretion model of WD metal pollution (Dupuis et al. 1993), which is known to have serious problems.

Observed abundances of heavy elements in WD atmospheres and theoretical calculations of their gravitational settling imply typical metal accretion rates $\dot{M}_Z\sim 10^6\hbox{--}10^{10}$ g s⁻¹ (Farihi et al. 2009, 2010). An evolving disk of debris must be able to supply such high $\dot{M}_Z$ to the WD. However, the question of how high-Z elements get transported to the WD atmosphere from the ring of solid particles, which does not extend all the way to the WD surface, has not yet been answered.

In principle, a dense ring of particles should evolve simply because of the angular momentum transfer due to inter-particle collisions, in full analogy with the rings of Saturn. However, the evolution timescale of Saturn's rings due to this process is too long, ∼10⁹ yr (Salmon et al. 2010), and resultant values of $\dot{M}_Z$ are negligible.

Another natural mechanism driving debris toward the WD is due to stellar radiation interacting with the disk and giving rise to the Poynting–Robertson (PR) drag (Burns et al. 1979). Previously, Farihi et al. (2010) claimed that PR drag cannot provide $\dot{M}_Z$ higher than 10³–10⁴ g s⁻¹ and dismissed this process as irrelevant. The goal of this paper is to critically re-examine the effect of radiative forces on the debris disk evolution and to show in particular that the PR drag can give rise to $\dot{M}_Z$ inferred from observations.

We envisage the following conceptual picture of the circum-WD environment. A dense disk (or ring) of particles lies inside the Roche radius R_R of the WD and evolves under the action of external agents, e.g., radiation forces. Particles migrate through the disk toward the sublimation radius R_s, where their equilibrium temperature equals the sublimation temperature T_s:

$$\begin{eqnarray} &&R_s=\frac{R_\star }{2} \left(\frac{T_\star }{T_s}\right)^{2}\approx 22\, R_\star \, T_{\star,4}^2 \left(\frac{1500\,\mbox{K}}{T_s}\right)^2, \end{eqnarray} \tag{ 1 }$$

where R_⋆ is the WD radius, T_s ≈ 1500 K for silicate grains, and T_{⋆, 4} ≡ T_⋆/10⁴ K is the normalized stellar temperature T_⋆. Taking R_⋆ ≈ 0.01 R_☉ typical for massive (M_⋆ ≳ 0.6 M_☉) WDs (Ehrenreich et al. 2011), one finds R_s ≈ 0.2 R_☉, in agreement with observationally inferred inner radii of compact debris disks (Jura et al. 2007, 2009a). Thus, for cool WDs R_out/R_in is several.

Particles sublimate at R_s feeding a disk of metallic gas, which we assume to be transparent to stellar radiation (gaseous component has been detected in several WDs with debris disks, see Gänsicke et al. 2006). Gaseous disk viscously evolves, extending all the way to the WD surface and providing means of metal delivery from R_s to the star. Assuming conventional α-parameterization of viscosity (Shakura & Sunyaev 1973) ν = αc²_s/Ω_K, the characteristic viscous time in the disk of metallic gas is

$$\begin{eqnarray} &&t_\nu \sim \frac{R_s^2}{\nu }\approx 10\, \mbox{yr}\, \frac{\mu _{28}}{\alpha }\frac{1500\, \mbox{K}}{T} \left(M_{\star,1} \frac{R_s}{0.2\, R_\odot }\right)^{1/2}, \end{eqnarray} \tag{ 2 }$$

where M_{⋆, 1} ≡ M_⋆/M_☉, and μ₂₈ is the mean molecular weight of the metallic gas normalized by 28 m_p (value of μ for pure Si). This timescale is short enough for $\dot{M}_Z$ to be determined predominantly by evolution of the disk of solids.

We take the disk of particles to be geometrically thin and characterize it at each point by optical depth τ defined as

$$\begin{eqnarray} &&\tau \equiv \frac{3}{4}\frac{\Sigma _d}{\rho a}, \end{eqnarray} \tag{ 3 }$$

where Σ_d is the surface mass density of debris disk, and ρ and a are the bulk density and characteristic size of the constituent particles (we leave a unspecified for now), respectively. Observations suggest (Jura 2003; Jura et al. 2007) that τ ≳ 1 but we leave τ arbitrary.

Particles in a narrow annulus just outside the sublimation radius are exposed to direct starlight at normal incidence and are heated to sublimation temperature. Farihi et al. (2010) studied the effect of PR drag on particles in this region. Here, we first explore the effect of radiative forces on the rest of the disk, which has only its surface illuminated by the central star at grazing incidence. To account for the finite size of the WD we adopt a "lamp-post" model of stellar illumination, in which radiation is emitted by two point sources of luminosity L_⋆/2 each located at height Z = (4/3π) R_⋆ above and below the disk plane (see Figure 1). Then,

$$\begin{eqnarray} &&\alpha \approx \tan \alpha =\frac{4}{3\pi }\frac{R_\star }{r} \end{eqnarray} \tag{ 4 }$$

is the incidence angle of radiation on the disk surface at distance r from the star, and α ≪ 1 since r ⩾ R_s ≫ R_⋆. The one-sided incident energy flux per unit surface area of the disk is F_r = (L_⋆/6π²)(R_⋆/r³), in agreement with Friedjung (1985).

Zoom In Zoom Out Reset image size

Radiation falling on the disk gives rise to two forces that can change angular momentum of the disk: PR drag and the Yarkovsky force. The latter owes its existence to thermal inertia of spinning objects that causes re-emission of absorbed energy in the direction different from that of incoming radiation (Bottke et al. 2006). Azimuthal components of these forces per unit area of the disk can be represented as

$$\begin{eqnarray} &&f_\varphi = \alpha \frac{L_\star \phi _r}{4\pi r^2 c}\psi, \end{eqnarray} \tag{ 5 }$$

where c is the speed of light and ψ is the factor that characterizes the phase lag of radiative momentum deposition in disk particles and is specified later in Sections 2.1 and 2.2 for PR drag and Yarkovsky force, respectively. Factor ϕ_r characterizes the efficiency of radiative momentum absorption by the disk surface. In the geometrical optics limits (particle size a much larger than characteristic wave length of stellar radiation)

$$\begin{eqnarray} &&\phi _r\approx 1-e^{-\tau _\parallel },\quad \tau _\parallel \equiv \alpha ^{-1}\tau, \end{eqnarray} \tag{ 6 }$$

where τ_∥ is the optical depth encountered by incident photons as they traverse the full disk thickness. Optically thick debris disks have ϕ_r = 1.

Azimuthal force causes radial drift of disk material at speed v_r = 2f_φ/(Ω_KΣ_d), where Ω_K is the Keplerian angular frequency. This gives rise to mass transport at the rate

$$\begin{eqnarray} &&\dot{M}_Z(r)=2\pi r v_r\Sigma _d=\alpha \psi \frac{L_\star \phi _r}{\Omega _{K} r c}. \end{eqnarray} \tag{ 7 }$$

Note that for radiative forces both f_φ and $\dot{M}_Z$ are independent of Σ_d.

It is reasonable to neglect the effect of potentially non-zero WD magnetic field on the evolution of particulate disk. However, simple estimate (Elsner & Lamb 1977) shows that for $\dot{M}_Z\sim 10^8$ g s⁻¹ magnetic field stronger than 1 kG at the WD surface would suffice to disrupt the gaseous disk at R_s. For simplicity, we assume here that the WD has weaker field and disregard complications related to magnetic effects.

2.1. Poynting–Robertson Drag

In the case of PR drag factor ψ is (Burns et al. 1979)

$$\begin{eqnarray} &&\psi _{\rm PR}=\frac{\Omega _K r}{c}\approx 3.3\times 10^{-3} \left(M_{\star,1}\frac{0.2\,R_\odot }{r} \right)^{1/2} \end{eqnarray} \tag{ 8 }$$

leading to

$$\begin{eqnarray} \dot{M}_{\rm PR}(r)&=&\frac{4 \phi _r}{3\pi } \frac{R_\star }{r} \frac{L_\star }{c^2} \nonumber\\ &\approx & 10^8\,\mbox{g s}^{-1}\phi _r \frac{L_\star }{10^{-3}L_\odot }\frac{20}{r/R_\star }, \end{eqnarray} \tag{ 9 }$$

where we used Equations (4) and (7).

The numerical estimate in Equation (9) does not agree with the calculation of Farihi et al. (2010) and we demonstrate the origin of this discrepancy by following these authors and considering the inner edge of the disk where particles sublimate (at r = R_s). At the edge, there is an annulus of the disk directly exposed to starlight (i.e., τ_{∥, a} ≲ 1; see Figure 1) with radial width

$$\begin{eqnarray} &&l\sim \tau _a^{-1}\,h, \end{eqnarray} \tag{ 10 }$$

which follows from the relationship τ_a/h ∼ τ_{∥, a}/l (here h ≪ R_⋆ is the disk thickness, τ_a and τ_{∥, a} ≲ 1 are the vertical and horizontal optical depths of the annulus). Particles inside this annulus experience PR drag as individual objects, with little shadowing from nearby particles. The characteristic decay time of their semimajor axes is (Burns et al. 1979)

$$\begin{eqnarray} &&t_{\rm PR}=\frac{4\pi }{3}\frac{\rho a R_s^2 c^2}{L}. \end{eqnarray} \tag{ 11 }$$

However, it only takes the fully exposed particles time, t_cross ∼ t_PR(l/R_s), to cross the annulus, reach the sublimation point, and evaporate, at which point they stop shadowing objects behind them. Dividing the mass 2πR_slΣ_d contained within the annulus by t_cross, we derive the mass accretion rate inside the annulus:

$$\begin{eqnarray} &&\dot{M}_{a}\sim \frac{2\pi R_s l\Sigma _d}{t_{\rm PR}}\frac{R_s}{l} \sim \tau _a\frac{L}{c^2}, \end{eqnarray} \tag{ 12 }$$

where we used definition (3). Mass conservation requires that $\dot{M}_{a}=\dot{M}_{\rm PR,s}\equiv \dot{M}_{\rm PR}(R_s)$ so that the optical depth inside the annulus is τ_a ∼ α_sϕ_r ≪ 1 (α_s ≡ α(R_s)) and from Equation (10) the width of the annulus is h/(α_sϕ_r). Thus, as long as the disk outside the annulus is not transparent to incident stellar radiation (i.e., τ_∥ ≳ 1 or τ ≳ α) τ drops from its value τ_s outside of R_s + l to τ_a inside. This is because fully exposed particles in the annulus experience stronger PR drag and migrate toward the WD faster than particles outside the annulus. The transition occurs at the distance l outside the disk edge, where τ_{∥, a} ∼ 1 (see Figure 1).

From observational point of view, the most interesting characteristic of the disk is the value of $\dot{M}_{\rm PR}$ at sublimation point R_s, since this is the rate at which disk of solids feeds gaseous disk, from which mass is transported onto the WD atmosphere by viscous torques. Evaluating expression (9) at r = R_s using Equation (1), we find

$$\begin{eqnarray} \dot{M}_{\rm PR, s}&=&\frac{32\phi _r}{3}\sigma \left(\frac{R_\star \, T_\star \, T_s}{c}\right)^2 \nonumber\\ &\approx & 7\times 10^7\,\mbox{g\, s}^{-1}\phi _r\left(R_{\star,-2} T_{\star,4}\frac{T_s}{1500\,\mbox{K}}\right)^2, \end{eqnarray} \tag{ 13 }$$

where R_{⋆, −2} ≡ R_⋆/10⁻² R_☉; then $\dot{M}_{\rm PR}(r)=\dot{M}_{\rm PR, s}(R_s/r)$. Note that $\dot{M}_{\rm PR,s}$ depends on T_⋆ only quadratically.

Previously, Graham et al. (1990) obtained an empirical estimate of $\dot{M}_{\rm PR}$ by assuming accretion to be driven by PR drag and dividing the observationally inferred disk mass (surface area multiplied by ρa, taking disk to be a monolayer of particles) by t_PR, see Equation (11). Even though their numerical estimate of $\dot{M}_{\rm PR}\approx 1.8\times 10^8$ g s⁻¹ is not very different from reality, this method of deriving $\dot{M}_{\rm PR}$ is not well motivated. Later, Farihi et al. (2010) found a value of $\dot{M}_{\rm PR}\lesssim 10^4$ g s⁻¹ much lower than $\dot{M}_{\rm PR,s}$ predicted by Equation (13) because they used t_PR instead of t_cross ≪ t_PR in computing $\dot{M}_{a}$.

2.2. Yarkovsky Force

If circum-WD debris disks behave similarly to dense planetary rings, one expects collisions in the disk to align particle spins predominantly parallel or anti-parallel to the disk normal, with non-zero and positive average spin (prograde mean rotation) and significant random spin component (Salo 1987; Ohtsuki & Toyama 2005). The characteristic spin frequency is ω ∼ Ω_K. Thermal inertia of spinning particles causes re-emission of absorbed stellar energy in the direction different from the radial. This gives rise to azimuthal Yarkovsky force, which causes outward orbital migration of disk particles spinning in prograde sense (Bottke et al. 2006). We now check whether this force can affect the PR-driven inward accretion of solids by estimating its typical magnitude.

Yarkovsky phase lag factor ψ_Y depends on two dimensionless parameters. One of them R' = a/l_ν is the ratio of particle size a to the penetration depth of the thermal wave $l_\nu =\sqrt{K/\rho C_p\omega }$, where K and C_p are the thermal conductivity and specific heat of particle material. Another is the thermal parameter (roughly the ratio of thermal time at depth l_ν to object's spin period) $\Theta =\sqrt{K\rho C_p\omega }/\sigma T^3$, where T is the surface temperature of the body. The actual dependence of ψ_Y on R' and Θ has been calculated in Vokrouhlický (1998, 1999) and these references should be consulted for details. For our current purposes it suffices to note that phase lag factor takes on a maximum value ψ_Y, max ∼ 0.1 when R' ≳ 1 and Θ ∼ 1, and is smaller for other values of R' and Θ.

To evaluate l_ν and Θ, we assume that disk particles have composition typical for terrestrial silicate minerals such as olivine, in agreement with elemental abundances of metal-rich WD atmospheres (Zuckerman et al. 2007; Klein et al. 2010). Then K ≈ 1.5 × 10⁵ erg s⁻¹ cm⁻¹ K⁻¹ and C_p ≈ 10⁷ erg g⁻¹ K⁻¹ at T ≳ 10³ K (Roy et al. 1981). Taking ω = Ω_K, we can write

$$\begin{eqnarray} &&l_\nu =l_{\nu,s}\left(\frac{r}{R_s}\right)^{3/4},\quad \Theta =\Theta _s\left(\frac{r}{R_s}\right)^{3/2}, \end{eqnarray} \tag{ 14 }$$

where l_{ν, s} and Θ_s are the values of l_ν and Θ at sublimation radius r = R_s:

$$\begin{eqnarray} l_{\nu,s}&=& \left(\frac{K}{\rho C_p}\right)^{1/2} \left(\frac{R_\star ^3}{8\,GM_\star }\right)^{1/4} \left(\frac{T_\star }{T_s}\right)^{3/2}\nonumber \\ &\approx & 1\, \mbox{cm} \frac{R_{\star,-2}^{3/4}T_{\star,4}^{3/2}}{M_{\star,1}^{1/4}} \left(\frac{1500\,\mbox{K}}{T_s}\right)^{3/2}, \end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray} \Theta _s&=& \left(\frac{3\pi }{8}\right)^{3/4} \frac{\left(K\rho C_p\right)^{1/2}}{\sigma T_s^3} \left(\frac{GM_\star }{R_\star ^3}\right)^{1/4} \nonumber\\ &\approx & 10 \frac{M_{\star,1}^{1/4}}{R_{\star,-2}^{3/4}} \left(\frac{1500\,\mbox{K}}{T_s}\right)^{3}. \end{eqnarray} \tag{ 16 }$$

Particle size a needed for calculation of R' is not constrained by observations, although Graham et al. (1990) suggested that debris disk may consist of 10–100 cm particles. If a is several centimeters or larger (so that R' ≳ 1) then ψ_Y, max ∼ 0.1 for fully illuminated particles at R_s (see Figure 3 of Vokrouhlický 1998).

However, in the dense debris disk only the upper (or lower) parts of particles, close to their spin axes, are directly illuminated by anisotropic starlight. This is because the latter illuminates the disk at small angle α and most of the particle surface is shadowed by other particles nearby (thermal emission of these particles is on average isotropic in horizontal direction). From heuristic arguments, one then expects stellar heating to be mainly deposited at small colatitude ∼α^1/2 ≪ 1 from the particle spin axis. It is easier for thermal conduction to isotropize surface temperature distribution around the spin axis over this small polar cap region rather than over the whole particle surface (as is assumed in standard calculation of ψ_Y), and this additionally lowers ψ_Y.

Thus, even if material properties of disk particles are favorable for maximizing ψ_Y (as is the case for our estimates here) one still should expect ψ_Y ≲ 0.01. Equations (4) and (7) result in the following expression for Yarkovsky-induced outward mass accretion rate at R_s:

$$\begin{eqnarray} |\dot{M}_{Y,\,s}|&=&\frac{2^{9/2}\phi _r}{3}\psi _Y \frac{\sigma R_\star ^{5/2} T_\star ^3 T_s}{c(GM_\star)^{1/2}} \nonumber\\ &\approx &2.4\times 10^8\,\mbox{g\, s}^{-1}\phi _r\frac{\psi _Y}{0.01} \frac{R_{\star,-2}^{5/2}T_{\star,4}^{3}}{M_{\star,1}^{1/2}} \frac{T_s}{1500\,\mbox{K}}. \end{eqnarray} \tag{ 17 }$$

This is somewhat higher than $\dot{M}_{\rm PR,s}$, thanks to the high adopted ψ_Y = 10⁻² exceeding ψ_PR, see Equation (8). However, this estimate may be very optimistic because of our poor knowledge of material properties and sizes of disk particles.

Finally, we note that observations do not reveal the existence of metal-poor WDs with debris disks around them, which could plausibly exist if outward migration of particles due to the Yarkovsky force were effective at stopping the PR-driven inward mass accretion. This provides additional evidence that Yarkovsky effect does not play a significant role in the circum-WD debris disk evolution.

We now compare our theoretical predictions with data on $\dot{M}_Z$ inferred from observations² of metal-rich WDs. In Figure 2, we plot the values of $\dot{M}_Z$ from Farihi et al. (2009, 2010) versus stellar effective temperature T_⋆ separately for WDs with and without debris disks detected via IR excesses. We also plot our analytical prediction for $\dot{M}_{\rm PR,s}(T_\star)$ for different values of R_⋆ and T_s.

Zoom In Zoom Out Reset image size

It is easy to see from Figure 2 that all systems with detected debris disks exhibit $\dot{M}_Z$ higher than $\dot{M}_{\rm PR,s}$. This is very encouraging since whenever WD is orbited by particulate disk one expects PR drag to set a lower limit of the mass flux at the level of $\dot{M}_{\rm PR,s}$. This lower bound on $\dot{M}_Z$ is hard (if at all possible) to avoid. This constraint naturally explains a strong positive correlation between $\dot{M}_Z$ and the rate of disk occurrence found by Farihi et al. (2010): the fraction of stars with compact debris disks is significantly higher for WDs with high $\dot{M}_Z$ simply because systems with disks cannot have low $\dot{M}_Z<\dot{M}_{\rm PR,s}$. All this strongly suggests that the PR drag indeed plays an important and visible role in circum-WD debris disk evolution.

It is worth emphasizing that the lack of disk-hosting systems at low $\dot{M}_Z$ is not caused by some observational bias: disk-bearing WDs with $\dot{M}_Z\sim \dot{M}_{\rm PR,s}$ such as GD 16 (Farihi et al. 2009), GD 56 (Jura et al. 2009a), and G166-58 (Farihi et al. 2008) have (sometimes highly) significant detections of IR excesses by Spitzer.

It is also clear from Figure 2 that there are metal-rich systems both with and without disks, in which $\dot{M}_Z$ significantly exceeds $\dot{M}_{\rm PR,s}$ for a given T_⋆, sometimes by several orders of magnitude. The most dramatic case is that of DAZB WD GD 362 which has $\dot{M}_Z\approx 2.5\times 10^{10}$ g s⁻¹ (Farihi et al. 2009), exceeding $\dot{M}_{\rm PR,s}$ corresponding to its T_⋆ = 10, 500 K by a factor of ≳ 300. Existence of such objects clearly implies that at least in some cases an additional accretion mechanism must operate on top of the PR drag, giving rise to very high $\dot{M}_Z$. We suggest such a mechanism that naturally operates in the presence of a debris disk in Rafikov (2011).

Systems without detected debris disks tend to occupy a broad range of $\dot{M}_Z$, both above and below $\dot{M}_{\rm PR,s}$. They may be interpreted as WDs that were orbited by compact debris disks in the recent past but have now completely (or largely) lost their disks to accretion. Without an active external source, metals sediment out from atmospheres of these WDs resulting in lower inferred $\dot{M}_Z$ in systems that had longer time since their debris disk disappearance.

Lifetime of a debris disk of mass M_d accreting via the PR drag alone can be estimated as $M_d/\dot{M}_{\rm PR} \approx 2\hbox{--}4$ Myr for M_d = 10²² g (roughly corresponding to the mass of 200 km asteroid) and $\dot{M}_{\rm PR}$ typical for a WD with T_⋆ = 10⁴ K. This estimate only marginally agrees with observational constraints suggesting a several 10⁵ yr disk lifetime (Farihi et al. 2009). However, one should keep in mind that this is an upper limit on the lifetime that is realized if no other processes apart from the PR drag drive accretion of solids through the disk.

In conclusion, we would like to emphasize the robust nature of the mass accretion due to PR drag. Our estimate (13) of $\dot{M}_{\rm PR,s}$ depends neither on the disk properties such as Σ_d (as long as τ_∥ ≳ 1) or M_d nor on the material properties of its constituent particles. Also, $\dot{M}_{\rm PR,s}$ does not vary much as T_s or R_⋆ change within reasonable limits.

Walker & Mészáros (1989) have previously suggested that rapid rotation of the central object can suppress the PR drag on the surrounding disk. However, it is easy to demonstrate that even if the WD rotates at breakup speed, the amount of angular momentum carried by photons absorbed by the disk is still small compared with the angular momentum that the disk loses via the PR drag.

One may worry that our calculation essentially disregards the details of the particle sublimation process. However, $\dot{M}_{\rm PR,s}$ is set by stellar illumination outside the innermost annulus of directly exposed particles at the inner edge of the disk (see Section 2.1). Their sublimation is irrelevant and our calculations are robust. Interior of that point, $\dot{M}_Z$ does not change by continuity and, as a result, the mass accretion rate ends up being insensitive to exactly how directly exposed particles sublimate.

The author thanks Konstantin Bochkarev for useful discussions. The financial support of this work is provided by the Sloan Foundation and NASA via grant NNX08AH87G.