Evaporation and accretion of extrasolar comets following white dwarf kicks

Abstract

1 INTRODUCTION

Metal absorption lines have been observed in white dwarf (WD) atmospheric spectra for many decades, but only recently has their importance for the study of exoplanetary systems become apparent. Observations of metals in WD atmospheres were initially puzzling: WDs with either hydrogen (DA) or helium (DB) dominated atmospheres possess relatively short sedimentation time-scales for heavier elements. Given that early efforts to explain these observations through convective instability and dredge-up have been demonstrated unlikely, the atmospheric metal content of WDs must come from an external source. Direct accretion from the interstellar medium fails to reproduce observed abundances by several orders of magnitude, and furthermore would predict metal pollution to be correlated with the Galactocentric orbits of WDs in ways that are not observed (Farihi et al. 2010). The general consensus at present is that these metals are actively accreted from the tidal disruption of small rocky bodies, likely asteroids that survived their host star's post-main-sequence (MS) evolution and then were excited to high eccentricity by dynamical encounters with surviving planets (Debes & Sigurdsson 2002; Jura 2003). This assessment is strengthened by the detection of dusty debris discs around a significant fraction of WDs; these discs are only seen in association with metal-polluted WDs and are estimated to occur around ∼1–3 per cent of all single WDs (Farihi, Jura & Zuckerman 2009).

Metal pollution in WDs has already been used to study the chemical composition and planetary architecture of rocky bodies in post-MS systems (Jura et al. 2007). Less has been said on the volatile-rich bodies which likely form in most stellar systems beyond the ice line. This is largely because measured abundances (Klein et al. 2011; Gänsicke et al. 2012) have disfavoured tidal disruption of hydrogen-rich objects as the primary source of metal pollution (Jura & Xu 2012), although Farihi, Gänsicke & Koester (2013) provides a notable exception. However, it is still plausible that tidal disruption of icy planets or planetesimals may contribute in a subdominant way to the total metal influx. This class of sources could be subdominant at all times, or perhaps only in a time-averaged sense. Past investigations of this possibility have focused on perturbations to surviving OCAs from external sources such as passing stars (Debes & Sigurdsson 2002), or scattering from planets that have survived post-MS evolution (Bonsor & Wyatt 2012).

We further explore the disruption and accretion of volatile-rich bodies in this paper, focusing on the orbital evolution of comets in exoplanetary Oort cloud analogues (hereafter OCAs). In particular, anisotropic mass-loss during the birth of a WD has the potential to radically alter OCA orbits, unbinding most (Parriott & Alcock 1998) but putting a small fraction on nearly radial orbits. This paper is motivated by two puzzling observations concerning accretion on WDs.

• Since the pioneering discovery of an IR excess around the central WD of the Helix nebula (Su et al. 2007), subsequent observations have established that ∼15 per cent of very young WDs possess IR excesses that can be explained by large-scale debris discs (Chu et al. 2011). The young WDs referred to here generally have effective temperatures ≳10⁵ K and ages ≲1 Myr, and often still have associated planetary nebulae (PN). Their dust discs have a large inferred radial extent, |${\sim } 10\hbox{-}100\ {\rm {\rm au}}$|⁠, in stark contrast with the |${<}0.01\ {\rm {\rm au}}$| discs that are observed around older DAZ or DBZ WDs.

• Many WDs are extremely hydrogen-deficient, with less H in their atmospheres (M_H < 10⁻¹⁵ M_⊙) than would be delivered by a single, modestly sized comet (Bergeron et al. 2011). This fact has been used in the past to constrain the density of interstellar comets (Jura 2011), but has not so far been applied to post-MS dynamics of OCAs.

The evaporation of a significant fraction of an OCA's mass following a natal WD kick would have clear implications for the first of these observations, and would be strongly constrained by the second. However, evaporation of comets is not the only explanation for IR excesses seen around young WDs. Debris discs around hot young WDs have been explained by enhanced collisions (Su et al. 2007) among analogues to Kuiper Belt Objects (KBOs), driven by orbital perturbations from surviving gas giant exoplanets (Dong et al. 2010). However, such a mechanism both requires survival of an appropriate exoplanet and of a sizeable population of KBO analogues. The latter of these assumptions may be questionable. An alternative explanation for these excesses is the formation of a dust-rich disc due to binary interactions with the asymptotic giant branch (AGB) wind (Bilíková et al. 2012), in analogy with the hot, dusty discs seen around post-AGB stars. This possibility has gained credibility with the recent discovery of binary companions around a majority of young WDs with IR excesses (Clayton et al. 2014); we discuss both of these alternate hypotheses in Section 5.

We emphasize that the well-studied IR excesses seen around older, metal-polluted WDs are not the focus of this paper. We are instead concerned with the poorly understood IR excesses seen around newborn WDs. The dynamical model we introduce for the response of an OCA to natal WD kicks is meant to produce large amounts of dust on |${\sim } 10\ {\rm {\rm au}}$| scales shortly after the birth of a WD, and has no connection to the |${\lesssim } 0.01\ {\rm {\rm au}}$| size dusty discs seen in a subpopulation of DAZ/DBZ stars.

The following sections explore the orbital dynamics, destruction, debris evolution, and observability of this deeply plunging subpopulation of comets. Specifically, in Section 2 we produce Monte Carlo samples of OCA objects (Section 2.1), evolve their orbital elements in response to both gradual and impulsive mass-loss (Section 2.2), and verify the validity of the impulsive approximation (Section 2.3). We next consider the fate of deeply plunging comets (Section 3), which can be destroyed either through tidal disruption (Section 3.1), or, more likely, through sublimation (Section 3.2). The properties of the cometary debris are considered in Section 3.3. In Section 4, we model the evolution and observable properties of the gaseous and solid-state debris discs formed from the leftovers of comet sublimation. Section 5 provides a general discussion of the observational implications for OCAs around newborn WDs. These include observable (and perhaps already observed) debris discs (Section 5.1), implications for WD natal kicks (Section 5.2), hydrogen abundances in tension with the most hydrogen-depleted WDs observed (Section 5.3), and surviving exoplanetary systems (Section 5.4). We summarize our conclusions in Section 6.

We pay special attention to the two sets of observations mentioned above, focusing on both the observability of the solid-state debris, and the total H accretion on to the WD. As we will show, cometary debris discs around young WDs are promising sources for the James Webb Space Telescope (JWST). The properties of such discs provide a probe of the properties of extrasolar Oort clouds and hence, indirectly, of the architectures of extrasolar planetary systems necessary to produce such a reservoir of icy bodies.

2 RESPONSE OF COMETARY ORBITS TO mass-loss

As an AGB star sheds its envelope, the orbits of planets and planetesimals evolve in response. Purely axisymmetric mass-loss (that is also reflection-symmetric about the equatorial plane) increases the energy of surrounding orbits, while leaving all other orbital elements unchanged. However, deviations from axisymmetry (or reflection symmetry) enable orbital eccentricities to increase. Likewise, the rate of mass-loss defines three qualitatively different regimes: orbits close to the star evolve adiabatically, orbits far from the star see the mass-loss as impulsive, and in between lies a more complicated transition regime.

Non-axisymmetric mass-loss may occur during this stage, for example due to magnetohydrodynamic (MHD) instabilities in AGB winds arising from magnetic cool spots on the stellar surface (Thirumalai & Heyl 2010). Dust-driven AGB superwinds account for the plurality of post-MS mass-loss for WDs with zero-age MS masses M_ZAMS ≳ 1.5 M_⊙, and contribute significantly to total mass-loss for M_ZAMS ≳ 1.2 M_⊙ (Wachter et al. 2002). The extremely short duration of the superwind phase (t_SW ∼ 3 × 10⁴ yr) renders it the only plausible epoch capable of accounting for the (sometimes disputed) dynamical evidence for WD kicks observed in globular and open clusters (see Section 5.2).

In contrast to planets and the Kuiper Belt, the majority of the Solar Oort cloud is concentrated at radii ≳1000 au ≳ a_imp, placing the response of OCAs to AGB mass-loss into the impulsive regime to first order. As a result, OCA objects whose velocity vectors align with |$\boldsymbol {v}_{\rm k}$|⁠, the kick velocity vector of the WD, can find their post-kick angular momentum decreased to very low values. In this section, we use Monte Carlo simulations to explore the parameter space of plausible extrasolar Oort clouds and WD kicks, paying special attention to the subpopulation of comets placed on to deeply plunging orbits.

2.1 Extrasolar Oort clouds

Although observations of debris discs around young stars provide some evidence for extrasolar comets (Ferlet, Vidal-Madjar & Hobbs 1987; Beust et al. 1990; Melnick et al. 2001; Rodigas et al. 2014), it is unclear whether the progenitors of these discs originate in a coplanar, KBO-like configuration, or in a more spherical OCA. Due to our limited information about the planetary dynamics of these systems, little can be said about the masses or radial distribution of their extrasolar comet reservoirs. For these reasons, our models of OCAs are based on our (limited) knowledge of the Oort cloud surrounding our own Solar system.

The Solar Oort cloud is a loosely bound, collisionless collection of icy bodies; based on observations of long-period comets, it is estimated to contain ∼10¹² objects over a kilometre in size (Weissman 1996). The Oort cloud is roughly spherically symmetric, with an inner boundary of |$R_{\rm in} \approx 10^3\ {\rm {\rm au}}$| set during its formation and modified by the Sun's orbit through the galaxy (Kaib, Roškar & Quinn 2011). Often an inner Oort cloud beginning at this radius is distinguished from an outer Oort cloud starting at |${\sim } 10^4\ {\rm {\rm au}}$| (Fernández 1997). Mass estimates for the inner Oort cloud, the region most relevant for this paper, are highly uncertain because of the limited influence of the Galactic tide in this region (Dones et al. 2004) and its consequent underproduction of observed long-period comets. Theoretically, the inner regions of the Oort cloud are less thermalized in both angles and eccentricity than the outer regions, with a bias towards radial orbits somewhat aligned with the plane of the ecliptic.

The total estimated mass of the entire cloud ranges between 2 and 60 M_⊕ (Weissman 1996; Kaib & Quinn 2008; Kaib et al. 2011), although the most recent dynamical work on this subject has favoured lower estimates for the total Oort cloud mass, ∼7 M_⊕ (Kaib & Quinn 2009). The Oort cloud is truncated by both the Galactic tide and encounters with passing stars, setting an outer edge |$R_{\rm out} \approx 5\hbox{--}20\times 10^4 \ {\rm {\rm au}}$| (Duncan, Quinn & Tremaine 1987). Although the detailed formation of the Oort cloud is not fully understood, a schematic picture is sketched in Appendix A. In this simplified analytic model (Tremaine 1993), a coplanar belt of planetesimals whose orbits cross those of large planets are perturbed to more weakly bound energies, and are eventually isotropized by the galactic tide.

More detailed properties of the Solar Oort cloud are predicted by simulations of planetesimal scattering off gas giants during the early history of the Solar system. These simulations typically predict that OCA objects possess a power-law distribution of semimajor axes n(a) ∝ a^−γ with γ ∼ 2–4, with numerical formation simulations generally finding γ ≈ 3.5 (Duncan et al. 1987; Wiegert & Tremaine 1999; Brasser, Duncan & Levison 2006; Kaib & Quinn 2008).

Our fiducial model assumes that cometary orbits possess an a^−3.5 semimajor axis distribution and a spherically symmetric distribution of inclinations, nodal angles, and lines of pericentre. We take |$R_{\rm in}=1000\ {\rm {\rm au}}$| and |$R_{\rm out}=5\times 10^4\ {\rm {\rm au}}$|⁠. For each value of a, we sample from a modified thermal distribution of eccentricities. The standard thermal distribution, P(e) = 2e, is truncated at high e so that no pre-kick pericentres below |$500\ {\rm {\rm au}}$| are sampled. This is the rough pericentre value for which our approximation of impulsivity qualitatively breaks down (Section 2.3), and is only a factor of ≈2 greater than the ice line during the star's AGB phase (Section 3.2). We take a final WD mass of M_WD = 0.6 M_⊙, but assume that it has lost half of its mass prior to the birth kick (we alter the a and e distributions to be consistent with this adiabatic and symmetric mass-loss). We also consider non-fiducial models that vary the assumption of spherical symmetry, γ, and R_in (our results are not sensitive to R_out).

2.2 The impulsive limit

Having generated a large (N = 10⁶) Monte Carlo sample of pre-kick comet orbits using the above prescriptions, these orbits are evolved in response to stellar mass-loss off the MS. This section considers the analytically tractable regime of impulsive mass-loss, leaving the more general case to the next section.

Orbital elements generated through Monte Carlo sampling are converted into six-dimensional Cartesian coordinates in position and velocity space. These coordinates are then translated evenly in velocity space to account for the WD kick, before being converted back to standard orbital elements. The transformations used by this procedure can be found in most celestial mechanics textbooks and are reviewed in Appendix B. This procedure is used to explore the effect of different combinations of kick velocity v_k and Oort cloud parameters {R_in, R_out, γ} on the distribution of post-kick cometary properties.

Fig. 1 shows the fraction of the OCA mass, f_OCA = dN/dlnR_p, deposited per log pericentre radius R_p following the kick, calculated for a range of kick velocities, assuming our fiducial model for the initial OCA properties. For the low pericentre radii of relevance, the distribution of comets is found to be roughly flat per unit pericentre, independent of kick velocity. However, there is an important dichotomy: when v_k ≲ 1 km s⁻¹, almost all of the low-R_p orbits are elliptical, but when v_k ≳ 2 km s⁻¹, the vast majority are hyperbolic. This has potentially important implications for debris evolution, although we argue in Section 3.2 that the large thermal spread in debris energy reduces the importance of this distinction.

Figure 1.

Open in new tabDownload slide

Fraction of the mass of the OCA, f_OCA = dN/dlnR_p, deposited on orbits with pericentre radius R_p following the natal WD kick. Distributions shown are calculated for kick velocities v_k = 0.05 (black), 0.1 (dark blue), 0.25 (blue), 0.5 (cyan), 1 (purple), 2 (red), and 4 (magenta) km s⁻¹, respectively. The initial properties of the OCA adopted are motivated by models of the Solar Oort cloud: inner radius |$R_{\rm in}=1000 \ {\rm {\rm au}}$|⁠, outer radius |$R_{\rm out}=5\times 10^4\ {\rm {\rm au}}$|⁠, and density profile n(a) ∝ a^−γ for γ = 3.5. Solid lines represent all comets, while dotted lines only include the fraction of comets that remain bound to the WD after the kick. The dashed green line, shown for reference, corresponds to an equal distribution of comets per linear pericentre distance.

The ‘sweet spot’ for impulsive velocity perturbations occurs for v_k ≈ 0.5 km s⁻¹. Kick velocities significantly below this value are too weak to maximize the fraction of OCA objects put on to nearly radial orbits. For kick velocities significantly above this, almost all OCA objects are on hyperbolic orbits, and the velocity space ‘loss cone’ for nearly radial orbits shrinks in size as v_k increases further. Intuitively, this sweet spot corresponds to the orbital velocity near the inner edge of the Oort cloud (or at the transition radius between impulsive and semi-adiabatic kicks, whichever is larger).

Fig. 2 shows the post-kick OCA mass distribution calculated for steeper (γ = 4.0) and shallower (γ = 2.5) density profiles than the fiducial case (Fig. 1). By concentrating more comets at the inner edge, the steeper distribution gives an enhanced fraction |$f_{\rm OCA}^{\prime }$| (≈50 per cent greater than the fiducial case) of OCA objects perturbed on to orbits with a given pericentre R_p ≪ R_in. Likewise, the shallower power-law index γ = 2.5 results in a ≈50–80 per cent decrease in |$f_{\rm OCA}^{\prime }$| for v_k ≳ 0.1 km s⁻¹, with the decrease becoming noticeably larger for greater kick velocities. A different value of γ also shifts the sweet spot for v_k; the steeper profile increases the kick velocity needed to maximize |$f_{\rm OCA}^{\prime }$|⁠, while the shallower profile reduces it.

Figure 2.

Open in new tabDownload slide

Fraction of the mass of the OCA, |$f_{\rm OCA}^{\prime }={\rm d}N/{\rm dln}R_{\rm p}$|⁠, deposited on orbits with pericentre radius R_p following the natal WD kick, calculated for two (non-fiducial) assumptions regarding the radial distribution of the OCA mass. The result is expressed as a fraction of the fiducial distribution f_OCA (Fig. 1) and shown with the same colour coding. Solid lines correspond to a more centrally concentrated OCA (γ = 4.0), while the dashed lines represent a shallower OCA profile (γ = 2.5), than the fiducial case, γ = 3.5. The limited resolution of the Monte Carlo sample is apparent as artificial noise at low values of R_p.

An important consideration in these calculations is the net angular momentum of the comets placed on to nearly radial orbits. If (as we will argue in the following sections) a large fraction of these comets sublimate and their orbital energy is dissipated, then their gaseous debris will eventually circularize at a radius determined by both the net post-kick angular momentum, |$\boldsymbol {J}_{\rm tot}$|⁠, and the time-scale for different gaseous debris streams to redistribute their own angular momentum in shocks. A stable, circular disc will form if |$|\boldsymbol {J}_{\rm tot}|>0$|⁠, as is possible only if the initial OCA possesses a net angular momentum that is not aligned with the direction of the WD kick. Even if |$|\boldsymbol {J}_{\rm tot}|=0$|⁠, short-lived gas discs may be able to form with very small radii due to the finite number of evaporating comets, and their resulting Poissonian excess in angular momentum. Multiple eccentric gas discs may be able to form if shocks are inefficient at redistributing angular momentum, but the dissipative interaction of these with each other will likely combine them into a single disc or inflow quickly. If |$\boldsymbol {J}_{\rm tot}=0$| due to the underlying symmetry of the system, gaseous cometary debris will likely accrete very quickly, either falling directly on to the WD surface or forming a very compact accretion disc.

Observations suggest that the net angular momentum of the Solar Oort Cloud is very small, as determined by the fraction 0.501 ± 0.051 of long-period comets that reside on prograde orbits with respect to the ecliptic (Wiegert & Tremaine 1999). This observation is in tension with theoretical models of the formation of the Oort Cloud, which generally predict an excess of comets on retrograde orbits due to the preferential ejection of prograde comets as a result of their longer gravitational encounters with massive planets (Wiegert & Tremaine 1999). The numerical simulations of Brasser et al. (2006) have found some net angular momentum to the inner Oort Cloud, albeit with error bars comparable in size to the angular momentum in question. More recent analysis of observations found no significant net angular momentum in observed Oort comets and suggested that past claims of a prograde bias have been due to selection effects (Wang & Brasser 2014).

We have replicated our fiducial scenario using an OCA with varying degrees of net rotation, by placing a fraction f_rot of comets on orbits with a preferred angular momentum direction.¹ We find that the resulting distribution of comets on nearly radial orbits preserves that specific angular momentum, such that the set of all comets with pericentres near R_p possess an average specific angular momentum |${\approx } f_{\rm rot}\sqrt{2GM_{\rm WD}R_{\rm p}}$|⁠.

2.3 Non-impulsive mass-loss

If WD natal kicks are delivered during an AGB superwind phase of duration ≈3 × 10⁴ yr, then we can compare this time-scale to an object's orbital time to estimate the validity of the impulsive approximation. As in Veras & Wyatt (2012), we define the adiabaticity parameter Ψ = t_orb/t_SW. The impulsive limit described in the prior section formally corresponds to Ψ ≫ 1, but for an object orbiting a central mass of M = 0.6 M_⊙ at the inner edge of an OCA (⁠|$a=2000\ {\rm {\rm au}}$|⁠), Ψ = 3.8 and the star's mass-loss will appear only moderately impulsive. To determine the error this introduces into our estimates of f_OCA, we have conducted simple two-body integrations for test particle orbits around central stars being kicked in a given direction over an adjustable time-scale t_k.

In general, we find that a narrow window of initial true anomalies f can be excited to near-parabolic eccentricities by central kicks. The precise location of this band of f does not asymptote to the impulsive limit until Ψ ≳ 10, but the width of the band does not change dramatically as long as it exists. Indeed, in the trans-adiabatic regime, the bias is towards a widening. However, once a kick becomes sufficiently non-impulsive, the initial orbital phase f does not matter and comets with all possible initial phases reach the same final eccentricity e^′.

We plot the dependence of e^′ on f for a wide range of Ψ in Fig. 3, in the special case where initial eccentricity e = 0 and the kick is within the orbital plane. We find that an impulsively excited band of f will produce e^′ arbitrarily close to 1 so long as Ψ > Ψ_c ≈ 1.5. Test integrations with different initial eccentricities and mutual inclinations (between the kick and the orbit) change the critical adiabaticity Ψ_c by factors of a few. For the innermost edge of our fiducial model after mass-loss (⁠|$a=2000\ {\rm {\rm au}}$|⁠), and orbits with inclinations i = π/4, we find that the minimum pericentre qualitatively consistent with the impulsive approximation (i.e. capable of producing radial orbits) is |${\approx } 500\ {\rm {\rm au}}$|⁠.

Figure 3.

Open in new tabDownload slide

Final eccentricity e^′ as a function of initial orbital phase f, for test particles on |$a=1000\ {\rm {\rm au}}$|⁠, e = 0 orbits around a central body receiving a kick v_k = 0.5 km s⁻¹. The different coloured lines represent different kick time-scales t_k; the red, blue, cyan, green, yellow, orange, and red lines indicate t_k/yr values of 10², 10³, 10⁴, 10^4.25, 10^4.5, 10^4.75, and 10⁵, respectively. The corresponding adiabaticity parameter values are labelled in the plot. The thick, dashed grey line (overlapping with the Ψ = 408 curve) represents the analytic prediction for the fully impulsive limit, as can be calculated using Appendix . All values of Ψ > 1.29 produce eccentricity excitation to near-parabolic values, and the Ψ = 1.29 case barely misses the window.

We note that a much more thorough treatment of anisotropic mass-loss in non-impulsive regimes has been recently given by Veras et al. (2013). The analytic formalism presented in that paper would likely be necessary to investigate OCAs with |$R_{\rm in} \lesssim 1000\ {\rm {\rm au}}$|⁠, but that is beyond the scope of this work.

3 DESTRUCTION OF DEEPLY PLUNGING COMETS

The minimum of the size distribution on the MS will be R_{c, min} ≈ 200 μm, as set by orbital decay due to Poynting–Robertson (PR) drag over the OCA formation time ∼1 Gyr (assuming a solar type star and distance |$50\ {\rm {\rm au}}$|⁠; see Section 4.3). This may not be a firm minimum due to collisional replenishment of small grains during OCA formation (Bonsor & Wyatt 2010); however, post-MS evolution will eliminate many of the smaller OCA particles. In particular, radiation pressure blowout on the AGB branch and drag from dust-driven winds will (at |${\sim } 1000\ {\rm {\rm au}}$| distances) both remove particles ≲ 1 cm in size (Bonsor & Wyatt 2010). We therefore take R_{c, min} = 1 cm, and set the maximum size to R_{c, max} = 1000 km, roughly the upper limit of observed KBOs in the Solar system. More details of OCA formation are discussed in Appendix A.

For the remainder of this paper, all our comets are treated as objects with constant density ρ_c = 0.6 g cm⁻³, to match observations of Solar system comets (Davidsson & Gutiérrez 2004; A'Hearn et al. 2005; Davidsson & Gutiérrez 2006). This is significantly lower than the density of their constituent solid-state silicates and ices, ρ_d = 3 g cm⁻³. This discrepancy arises because comets are generally loosely bound rubble piles with high internal porosity.

3.1 Tidal disruption

Much more precise treatments of tidal disruption in the solid-state limit exist (Davidsson 1999, for example), but these do not change the fact that the tidal disruption radius is extremely small (⁠|$R_{\rm t}^{\rm SG} \sim 3\times 10^{-3}\ {\rm {\rm au}}$|⁠, for the Ceres example). In practice, a much greater number of comets are destroyed through sublimation, which we focus on hereafter.

3.2 Evaporation

Comets with pericentres obeying R_ev < R_p < R_ice can also sublimate, but less efficiently and over many pericentre passages. Appendix C shows that this results in an enhancement of at most a factor of a few to the total mass lost from the largest comets; we include this modest enhancement to R_ev (assuming λ = 1.25, R_ev increases by a factor of ≈1.5).³

Fig. 4 shows the mass of the kicked-in OCA that actually sublimates in a given range of pericentre radii, calculated according to the criterion R_c < R_ev for the two mass functions (equations 2 and 3) which bracket our ignorance of the amount of collision processing experienced by the OCA. The less collisionally processed distribution sees ∼50 per cent of the kicked-in mass at a given pericentre sublimate, but for the top-heavy mass distribution of equation (3), only ∼10 per cent sublimates.

Figure 4.

Open in new tabDownload slide

Distribution of mass sublimated from comets as a function of the pericentre radius R_p in comparison to the total distribution of kicked-in comets f_OCA = dN/dlnR_p (solid lines; same colour labelling as Fig. 1). The fraction of comets that sublimate is calculated according to criteria given in equation (12) taking into account the size distribution of comets. Dotted lines are calculated assuming the optimistic ‘primordial’ mass function equation (2), while the dashed lines show the pessimistic, ‘collisionally processed’ mass function equation 3). Both dotted and dashed lines truncate at the ice line for L_WD = 10² L_⊙; if R_p > R_ice, no sublimation is possible.

The debris from sublimating comets can only produce observable consequences, such as forming a gaseous or dusty disc, if it remains bound to the WD. However, Fig. 1 shows that the vast majorities of comets are placed on to hyperbolic trajectories for kick velocities v_k ≳ 1 km s⁻¹. As we now discuss, the sublimation processes imparts a spread in the specific orbital binding energy Δϵ of the sublimated gaseous/dusty debris which can in some cases exceed the centre of mass specific orbital energy |ϵ₀|. In this case a significant fraction, ∼1/2, of the debris (that imparted with negative energy) will remain bound to the WD, even if the initial centre of mass orbit of the comet was not.

For orbits with R_p < R_hyper, roughly half of the cometary debris will become unbound, while the other half will remain bound to the WD with a wide range of orbital eccentricities and semimajor axes. For gas, this is generally the case when v_k ≲ 3 km s⁻¹. If the rocket effect scatters a large fraction of solid debris, the same criterion holds; otherwise, most sublimating orbits will only satisfy R_p < R_hyper if v_k ≲ 1 km s⁻¹. Conversely, when R_p > R_hyper, all the debris will remain ‘frozen in’ to orbits with roughly the same parameters as the centre of mass trajectory, which are bound (unbound) for v_k ≲ 1 km s⁻¹ (v_k ≳ 1 km s⁻¹). Because R_hyper > R_ice is generally satisfied for v_k ≲ 3 km s⁻¹, we conclude that to first order 50 per cent of all sublimated OCA gas will remain bound to the WD for these kick velocities.⁶ When v_k ≳ 3 km s⁻¹, a large majority of sublimated debris will instead exit the system, producing much weaker observational signatures. For the remainder of this paper, we assume 50 per cent of debris remains bound, but we note here that if the rocket effect is unimportant, most solid debris will exit the system for v_k ≳ 1 km s⁻¹.

The many critical radii discussed in this subsection are plotted as functions of time in Fig. 5. As the WD ages and cools, the sublimation radii for comets of fixed size contract quickly, while the ice lines move inward more slowly.

Figure 5.

Open in new tabDownload slide

Critical radii as a function of time (or, equivalently, the WD age). Radii shown include the tidal radius for a self-gravitating body similar to Ceres, |$R_{\rm t}^{\rm SG}$| (equation 4; black solid), sublimation radius for ice R_ice (equation 9; blue solid), sublimation radius for silicates, R_rock (equation 10; brown solid), and the radius interior to which most debris from comets on initially hyperbolic orbits remains bound to the WD, R_hyper (equation 15; purple). Two different R_hyper curves are shown: one for sublimated volatiles (purple solid; this may also correspond to the dust hyperbolic radius if the rocket effect is strong), and one for ejected dust (purple dashed). Red lines show the sublimation radii R_ev (equation 12) for comets of size R_c = 10 km (solid), 1 km (dashed), and 0.1 km ( dotted), respectively. We note that the sublimation radii have no physical meaning when located above their respective ice lines (for either water ice or silicates); for example, when R_ev > R_ice, portions of parameter space between these curves are where the comet absorbs enough (eventually re-radiated) energy at pericentre to fully sublimate, but cannot because its temperature is too low.

3.3 Cometary debris

Once one or more comets have sublimated, multiple streams of gas and solid debris will orbit the WD on highly eccentric trajectories. We parametrize the relative composition of these two constituents with f_gas, the volatile mass fraction of the comet prior to sublimation. We also account for the hydrogen mass fraction of the sublimated gas, f_H. We hereafter make the simplifying assumptions that f_gas = 0.5 and that water ice is the dominant volatile, implying f_H = f_gas/8 = 0.06.

Electron microscopy analysis of interplanetary dust grains indicates a lower size limit of b_min ∼ 0.1 μm (Bradley 1994), which we take as our fiducial value. The maximum dust particle size, b_max, is hard to measure remotely due to its generally subdominant contribution to reradiated sunlight from sublimated comets. However, NASA's Deep Impact mission to the comet 9P/Tempel 1 excavated a large crater on the comet's surface, enabling a rough calculation of b_max by comparing the ejecta mass estimates from volumetric considerations to those from IR emission (Küppers et al. 2005). Although this calculation is subject to large uncertainties, the most recent analysis of Tempel 1 results indicates that 1 μm ≲ b_max ≲ 100 μm (Gicquel et al. 2012). Subsequent Deep Impact observations of the coma of comet 103P/Hartley indicate the presence of much larger particles, up to ∼1 m in size. However, the inferred power-law slope of these ‘boulders’ is so steep that they would contribute negligibly to both the mass and surface area budget of ejected dust. In addition, density arguments favour an icy composition for the boulders, indicating that they would not be present for fully destructive sublimation events (Kelley et al. 2013). As we will show in Section 4.3, only particles of size b ≳ 50 μm can survive around young luminous WDs, motivating us to take b_max = 200 μm as fiducial. If such large dust grains are not in fact present, then both the hypothesis of this paper (OCA origins) and the KBO model in Su et al. (2007) would have serious difficulty in explaining the Spitzer observations of extended dust haloes around young WDs.

4 WHITE DWARF NATAL discS

4.1 Gaseous discs: evolution

Assuming that a fraction f_ev ∼ f_OCA ∼ 10⁻⁴–10⁻³ of the total mass M_OCA ∼ 10 M_⊕ ∼ 6 × 10²⁸ g of the OCA sublimates (i.e. those comets with pericentre radii interior to the sublimation radius R_ev; Fig. 4), then the resulting gas will circularize into a disc of mass M_gas = f_gasf_evM_OCA ∼ 10²⁵ g and characteristic radius |$R_{\rm circ} \sim f_{\rm rot}^{2}R_{\rm ice} \sim 0.3(f_{\rm rot}/0.1)^{2}L_{2}^{1/2}$| au, where f_rot is the rotation parameter of the original OCA (see Section 2.2).

The disc mid-plane will be photoionized by ultraviolet radiation from the WD, as is justified because the radial optical depth of the disc τ ∼ (Σ_g/H)κ_esr is less than unity at radii of interest, where κ_es ≈ 0.2 cm² g⁻¹ is the electron scattering opacity. Heating from ionization drives a powerful outflow from the disc exterior to the characteristic radius |$R_{\rm g} \simeq GM_{\rm WD}/c_{\rm s}^{2} \sim 30(T_{\rm g}/10^4\,{\rm K})^{-1}$| au (Hollenbach et al. 1994), with a mass-loss rate |$\dot{M}_w$| that is exponential in r/R_g (more detailed calculations find that |$\dot{{\rm M}}_w$| declines slowly with decreasing radius, until cutting off sharply at radii |$r <R_{\rm g}^{\prime }\approx 0.15R_{\rm g} \approx 5(T/10^4\,{\rm K})^{-1}\,{\rm {\rm au}}$|⁠; Adams et al. 2004). The rate of mass-loss due to photoevaporation is sufficiently high that the outer disc (⁠|$r \gtrsim R_{\rm g}^{\prime }$|⁠) will dissipate on a time-scale much shorter than the age of the system. Although |$R_{\rm g}^{\prime } \gg R_{\rm circ}$|⁠, viscous spreading will eventually bring some mass to radii vulnerable to photoevaporation. However, the bulk of the gas that circularizes interior to |$R_{\rm g}^{\prime }$| will accrete on to the WD, guaranteeing that a minimum H fraction be available to pollute the WD atmosphere.

In this subsection, we have treated these gas discs in an approximate manner, with the largest uncertainty being the unknown value of f_rot. However, in the following two sections we shall see that even a very approximate treatment is adequate, because (i) even if a mass M ≪ M_gas accretes on to the WD, there will still be observational consequences, and (ii) the effect of these discs on solid-state debris is negligible.

4.2 Gaseous discs: observational implications

The atmospheres of many WDs possess extremely low hydrogen abundances (e.g. Bergeron et al. 2011); DB atmospheres with less than 10⁻¹⁶ M_⊙ of hydrogen are not uncommon, and the most hydrogen depleted systems have upper limits below 10⁻¹⁸ M_⊙.⁷ The accretion of a single moderately sized comet would overproduce the observed hydrogen in many DB atmospheres, a fact which is puzzling in light of the comet accretion mechanisms described in previous sections. Because hydrogen does not sediment away as is the case with metals, extremely hydrogen-deficient WD atmospheres provide tight constraints on the joint parameter space of WD kick velocities and OCAs.

How much hydrogen will a natal kick deliver to the WD surface? The total hydrogen mass, M_H, accreted on to the WD via the transient gas disc discussed in the previous section is just that which will circularize interior to the photoevaporation radius |$R_{\rm g}^{\prime } \approx 5$| au. Since |$R_{\rm circ} \ll R_{\rm g}^{\prime }$| generally, M_H ≈ M_gasf_H where f_H ≈ 0.06 is the H fraction of typical OCA bodies. For a characteristic value M_gas = 10²⁵ g, we have M_H ∼ 6 × 10²³ g. This expected level of pollution is many orders of magnitude higher than the observed H abundances in DB WD atmospheres, and at least comparable to (often higher than) H abundances in DBA atmospheres (Bergeron et al. 2011).

4.3 Debris clouds: evolution

Fig. 6 shows the time evolution of the minimum grain size set by different processes at several characteristic radii. In general, radiation pressure blowout is the limiting factor at early times and large radii, while PR drag plays a larger role close to the WD and at late times, as the central luminosity declines. The characteristic minimum size is b_min ∼ tens of μm on time-scales t ∼ few 10⁵ yr, which characterize the free fall time for most OCA objects on nearly parabolic orbits.⁸

Figure 6.

Open in new tabDownload slide

Minimum grain size b_min that survives PR drag (solid lines; equation 26) and radiation blow-out (dotted black line; equation 27) as a function of the WD age, calculated using the WD cooling evolution in equation (8) for orbits with different pericentre radii: |$R_{\rm p,0} = 5\,{\rm {\rm au}}$| (red), |$10\,{\rm {\rm au}}$| (orange), and |$50\,{\rm {\rm au}}$| (blue), respectively. All cases are calculated assuming a PR efficiency Q_PR = 1, solid density ρ_d = 3 g cm⁻³, WD luminosity normalization L₀ = 10² L_⊙ (equation 8), and initial semimajor axis |$a_0=500\,{\rm {\rm au}}$|⁠, the latter characteristic of the debris following energy freeze-in during evaporation.

We estimate the time-scale for an eccentric collisional cascade using equation (16) of Wyatt et al. (2010), which determines whether or not the size distribution remains fixed in time. If we take characteristic values b_max = 100 μm, total dust mass ∼10⁻³ M_⊕, and semimajor axis |$a=100\,{\rm {\rm au}}$|⁠, we can estimate the collisional cascade time-scale if we know the catastrophic disruption threshold |$Q_{\rm D}^*$| and the dimensionless function S(e). We take S(e) = 1 (for 1 − e = 0.3; see Wyatt et al. 2010, fig. 3), and extrapolate fig. 3 of Benz & Asphaug (1999) to ∼ 100 μm scales to estimate |$Q_{\rm D}^* \sim 5\times 10^4\,{\rm J\ kg}^{-1}$| for rocky dust. These parameters give a collisional cascade time-scale t_cc ≈ 5 × 10⁶ yr, which is greater than the age of the systems we are concerned with. We assume a constant particle size distribution for the remainder of this paper, but note that collisional grinding of the solid debris over longer time-scales could produce a reservoir for eventual metal accretion on to the WD. Small dust grains produced over time-scales t > t_cc may eventually accrete due to PR drag.

4.4 Debris clouds: observational implications

The above discussion establishes that, after an initial phase of dissipation and circularization, any gaseous accretion disc will drain rapidly before it has time to entrain or circularize significant amounts of solid debris. The solid debris from cometary disruption will therefore remain on highly eccentric orbits, and will re-radiate absorbed WD luminosity at cooler temperatures. This eccentric debris cloud will manifest observationally as an IR excess in the WD spectrum. Such excesses are detected as point sources around ∼10 per cent of hot, young WDs (Su et al. 2007; Chu et al. 2011). It is also possible that such emission could be resolvable by future high-resolution IR imaging (e.g. JWST).

To calculate the detailed properties of IR emission from the cloud, we project its luminosity density to produce a 2D surface brightness profile, I(s). We make the following assumptions in our calculation.

• The cloud is spherically symmetric. This assumption is not entirely accurate, as the overlap between pre-kick and post-kick loss cones will preferentially place comets with orbital angular momenta parallel to the kick direction on to the lowest angular momentum trajectories. Deviations from spherical symmetry grow with v_k. Analysis of our MC samples indicates that for kicks in the sweet spot of v_k ≲ 1 km s⁻¹, the sphericity assumption is correct to leading order, but is incorrect for significantly greater v_k.

• The cloud is optically thin. This assumption is verified post-hoc for all models presented in this section.

• Dust particles travel on closed Keplerian ellipses, such that emission at radius r is calculated as the time-weighted average of every dust particle with a pericentre r_p < r and apocentre r_a > r.

• The dust is in thermal equilibrium with absorbed WD radiation at all times, so that |$T_{\rm d}^4=L_{\rm WD}/(16\pi \sigma r^2)$|⁠, where σ is the Stefan–Boltzmann constant and both emission and absorption cross-sections are geometric. The corresponding Wien peak wavelength is generally smaller than the characteristic dust particle size set by radiation pressure blowout and PR drag (2πb ≳ λ), justifying this assumption for λ ≲ 10² μm.

• All dust particles have a semimajor axis distribution dN/da given by equation (16).

• All dust particles are spheres with a power-law distribution of radii n(b) ∝ b^−L for L = 3.5 as in equation (18).

Fig. 7 shows the bolometric surface brightness profile calculated for WD luminosities of 1, 10, and 100 L_⊙, and kick velocities of 0.5, 1, and 2 km s⁻¹. In all of these profiles, we take as fiducial parameters: b_min = 0.1 μm; b_max = 100 μm; OCA mass M_OCA = 10 M_⊕; and time since the WD kick (important for calibrating dust depletion due to PR drag) of 10⁵ yr. In general, we find projected surface density profiles that peak between 1 and 10 au and are relatively constant interior to that radius. This constant-brightness core is created by a combination of projection effects and the truncation of the dust distribution at small radii by PR drag: the outer edge of the core corresponds to the radius where PR drag has eliminated almost all of the dust from our grain size distribution.

Figure 7.

Open in new tabDownload slide

Surface brightness I(s) versus projected radial distance s, calculated for WD luminosities L_WD = 100 L_⊙ (purple), 10 L_⊙ (blue), and  L_⊙ (cyan). Solid lines represent v_k = 1 km s⁻¹, dashed lines represent v_k = 0.5 km s⁻¹, and dot–dashed lines represent v_k = 2 km s⁻¹. In all sets of curves, we take initial minimum particle size b_min = 0.1 μm, maximum particle size b_max = 100 μm, and M_OCA = 10 M_⊕, although b_min is modified in a pericentre-dependent way to account for radiation forces. The vertical black lines show diffraction-limited JWST angular resolution for the two nearest detected debris discs around young WDs, at both 10 and 20  μm.

Figs 8 and 10 show sample SEDs at distances of 200 pc, calculated for different assumptions regarding the mass M_OCA, the WD luminosity, and the natal kick velocity v_k. Sensitivity curves for Spitzer and JWST are shown for comparison. Since F_ν ∝ M_OCAϒ, detectability is sensitive to both OCA mass and to the distribution of sublimated dust grain sizes n(b). In Fig. 8, we can see that decreasing the WD luminosity moderately dims and reddens the resulting SED. This reduction is less than might be expected, because the decreased temperature of the dust is partially compensated for by the increased value of b_min. Fig. 10 shows the comparably modest dependence of the SEDs on v_k. Reasonable choices for M_OCA and v_k produce SEDs that are observable by Spitzer at a distance 200 pc similar to that of the Helix, while JWST will be sensitive to clouds produced for a much wider range of parameters to significantly larger distances.

Figure 8.

Open in new tabDownload slide

Spectral flux density F_ν versus wavelength λ, calculated for two values of the WD luminosity L_WD = 100 L_⊙ (solid coloured) and 10 L_⊙ (dotted coloured). Blue lines show an ‘optimistic’ scenario corresponding to an initial OCA mass M_OCA = 50 M_⊕. The orange and red curves represent ‘standard’ and ‘pessimistic’ scenarios for which M_OCA = 5 M_⊕ and M_OCA = 0.5 M_⊕, respectively. All three scenarios assume an initial ice line (as determines the outer edge of the sublimated mass distribution dN/dR_p) of |$R_{\rm ice} = 30\,{\rm {\rm au}}$|⁠, a kick velocity v_k = 1 km s⁻¹, and a distance d = 200 pc corresponding to that of the Helix nebula. The limiting sensitivities of Spitzer (black dashed) and JWST (black solid) are shown for comparison.

Fig. 9 shows Spitzer data points for the 10 young WDs⁹ observed to contain significant IR excesses, in combination with predicted spectral flux densities F_ν. The plotted F_ν contains both dust emission from our model and the blackbody emission from the central WD; in both cases, L_WD, R_WD, and d are taken from the observations described in Chu et al. (2011). The primary free parameter in our dust model is M_OCA, which we fit to the 24  μm observations. Kicks of 0.5 km s⁻¹ are assumed. In general, we are able to reproduce the observed long-wavelength IR excesses with reasonable values of OCA mass. Specifically, we find that 24 μm emission from CSPN K 1−22, CSPN NGC 2438, WD 0103+732, WD 0109+111, WD 0127+581, WD 0439+466, WD 0726+133, WD 0950+139, WD 1342+443, and WD 2226−210 can be reproduced with OCA masses (in M_⊕) of roughly 145, 1225, 215, 2.0, 20, 26, 21, 365, 19, and 175, respectively. Half of the observed systems are fitted with Solar-type masses of comets (1.5 M_⊕ ≲ M_OCA ≲ 20 M_⊕). The other half require values of M_OCA an order of magnitude greater (almost two orders of magnitude for the exceptional case of CSPN NGC 2438), although three of these anomalously large M_OCA values may be artefacts of the extreme luminosities of the central WDs (CSPN K1−22, CSPN NGC 2438, and WD 0103), which make our model extremely sensitive to the upper cutoff in the assumed dust grain size distribution.

Figure 9.

Open in new tabDownload slide

Spectral flux density F_ν versus wavelength λ, calculated for all 10 young WD systems observed to contain an IR excess by Chu et al. (2011). The spectral flux density curves include both dust emission from equation (34) and blackbody emission from the central WD (and in the case of CSPN K1-22, blackbody emission from an unresolved companion star). The primary free parameter is the total mass of the Oort Cloud analogue, M_OCA, which is adjusted to fit the 24  μm data point. The 3–8  μm data points are not fit for because of their sensitivity to time-dependent effects not accounted for in our model. For five WD systems, this fitting gives Solar-like values for the total OCA mass, with 1.5  M_⊕ ≲ M_OCA ≲ 20  M_⊕. The other five systems require much larger OCA masses, with M_OCA/M_⊕ = {145, 1225, 215, 365, 175} for CSPN K1-22, CSPN NGC2438, WD 0103, WD 0950, and WD 2226, respectively. As in prior plots, the Spitzer and JWST sensitivity curves are shown as dashed and solid black lines.

Several of the sources in the Spitzer sample have secondary 3-10 μm IR peaks, which we do not fit for because of their sensitivity to time-dependent effects not accounted for in our exploratory model. Specifically, at high frequencies dust emission will be:

• decreased by the higher temperatures of the WD at earlier times; comets sublimated at these times will see a much larger fraction of their debris eliminated by stronger PR drag forces;

• increased by the steady-state PR flow to small radii.

Which of these effects wins out is not obvious, and would require development of a time-dependent PR drag model to estimate. The four systems observed in Chu et al. (2011) which possess a ∼ 5 μm emission bump may be indicative of high-temperature dust; we speculate that this may correspond to a steady-state PR flow at scales of |${\sim } 1\,{\rm {\rm au}}$|⁠.

5 DISCUSSION

The Spitzer observations of 10 young WDs possessing IR excesses were the primary motivation for the theoretical model produced in Section 4. As discussed in Section 4.4, the combination of gentle WD kicks and OCAs is capable of reproducing the 24 μm Spitzer data points seen in all 10 systems. Half of these systems are well fitted with M_OCA values in the Solar range (2 ≲ M_OCA/ M_⊕ ≲ 30), but the other half require OCA masses roughly an order of magnitude larger. This may reflect selection biases, as the largest OCAs would produce the most luminous debris discs around young WDs. However, our model is generally less capable of reproducing the hotter, second IR excesses seen in four of the Spitzer sample's WDs (Chu et al. 2011). It is possible that this could be remedied by incorporating the time-dependent PR flow of dust into our model, or it could indicate that OCAs are unable to explain IR excesses around some fraction of young WDs. We therefore summarize two alternative proposals in this section.

A commonly discussed explanation for debris discs around young WDs is a collisional cascade among KBO analogues (Su et al. 2007). Such collisions may be triggered by the destabilization of planetesimal orbits due to interactions with surviving planets, whose orbits would themselves have been altered by mass-loss from the central star (Dong et al. 2010). Although JWST may help discriminate between the KBO scenario and our model as described above, other arguments disfavour the KBO hypothesis. When the host star crosses the horizontal branch its luminosity reaches ≳ 10³ L_⊙ for a time-scale of ∼10⁷ yr. During this epoch, the water sublimation radius (equation 9) moves outwards to hundreds of au, far beyond the Solar Kuiper Belt. Although the comets in an OCA will generally orbit far enough away to survive, most volatile-dominated KBO analogues will likely be destroyed (Stern et al. 1990) well before the AGB phase.¹⁰ Large KBO analogues could survive sublimation in two ways: either by possession of a chemically differentiated, rocky interior, or by simply being large enough that their surface escape velocity v_esc ≫ v_th. Since a L = 6 × 10³ L_⊙ giant star will produce an equilibrium temperature T_d ≈ 330 K at |$r=50\,{\rm {\rm au}}$|⁠, the thermal velocity of water will be v_th ≈ 0.7 km s⁻¹ and sublimated ice will escape from even self-gravitating, Triton-sized objects. Therefore, chemical differentiation and production of rocky cores offers the more promising mechanism for survival of very large KBO analogues.

Although the sublimation of KBO analogues will leave behind a debris ring of gas and dust (as in the case of cometary debris described in Section 3.3), the stellar luminosity will be much higher when this occurs. If we assume a typical horizontal branch luminosity ≈6 × 10³ L_⊙, then equation (27) predicts that radiation pressure will remove all grains smaller than 2300 μm, while PR drag (equation 24) sets an even more stringent minimum size of ≈ 10⁴ μm (see also Section 4 of Bonsor & Wyatt 2010). These lower limits for grain size are orders of magnitude above the ≈ 100 μm maximum grain size inferred based on the Deep Impact results for Tempel 1 (Section 3.3). The KBO scenario may have serious difficulty producing enough dust that survives until the birth of the WD, unless surviving rocky cores have both a high enough spatial density to begin a collisional cascade, and sufficiently large total mass to reproduce observed IR luminosities.

A second alternative is the possibility that unresolved binary companions have created dust-rich discs during the WD's prior post-AGB phase (Chu et al. 2011; Bilíková et al. 2012). Analogous discs around post-AGB stars are generally associated with binarity (van Winckel et al. 2009), but, as was pointed out by Chu et al. (2011), circumstellar dust discs in post-AGB systems are typically too compact to produce an IR excess peaked at 24 μm (Taam & Ricker 2010). Although they possess roughly correct temperatures for matching the IRAC-band (3-8 μm) excess seen in 40 per cent of the Chu et al. (2011) sample, circumstellar discs generally have fractional luminosities L/L_WD ∼ 0.3 (de Ruyter et al. 2006), orders of magnitude higher than the observed IRAC-band excesses, or for that matter the 24 μm excesses. It remains unclear if the hot, very luminous discs observed around post-AGB stars are able to survive in attenuated form into and past the PN stage. During completion of this paper, new observational work appeared which indicates that 8 out of 13 CSPN stars with IR excesses also possess a binary companion (Clayton et al. 2014), strengthening the post-AGB disc scenario.

In the following subsections, we discuss the observational implications of this paper's model for IR excesses around young WDs. In addition to future IR observations of these dusty sources, we highlight more indirect implications as well, such as the atmospheric hydrogen content of old WDs, natal kick velocity distributions, and exoplanet survival during post-MS evolution.

5.1 Future observations of discs around young WDs

Future observations by JWST of debris discs around young WDs may distinguish between our cometary hypothesis and competing theories involving objects from within the plane of the ecliptic such as KBO analogues (Su et al. 2007; Dong et al. 2010). One ‘smoking gun’ discriminant would be to spatially resolve the central IR cavity to distinguish a quasi-spherical cloud from one with a disc-like concentration (see discussion at end of Section 4.3). Although a cavity of size |${\sim } 10\,{\rm {\rm au}}$| is predicted in both scenarios (due to a combination of radiation pressure blowout and PR drag), the cometary hypothesis predicts a relatively constant surface brightness profile due to projection effects (Fig. 7) rather than the complete absence of emission predicted for a disc-like configuration observed at most inclinations.

Though in principle a promising diagnostic, direct imaging of the central cavity will be challenging. For the nearest young WD debris disc (distance 129 pc; Chu et al. 2011) |$10\,{\rm {\rm au}}$| corresponds to an angular size 0.16 arcsec. This is similar to the anticipated intrinsic resolution of the Mid-Infrared Instrument (MIRI) on JSWT (0.11 arcsec pixel⁻¹), but the diffraction-limited resolution at 10 μm (20 μm) is substantially greater ≈0.45 arcsec (0.89 arcsec). It is therefore unlikely that the central cavities can be resolved, although obviously the discovery of a significantly closer (and likely fainter) young WD debris disc would improve the prospects for a resolved source substantially.

A more feasible alternative discriminant might be large-scale imaging of these dusty IR sources. Dust clouds produced by an OCA will generally appear circular in projection (unless f_rot ∼ 1, which is disfavoured for the Solar Oort cloud). However, if it is collisions between coplanar KBO analogues (or binary interactions) that produce the observed dust, then depending on the inclination angle of the debris disc, a range of ellipticities may result in projection. In Fig. 7, we show the diffraction-limited angular resolution capabilities of JWST for the two nearest debris systems around young WDs at 10 μm and 20 μm; if these systems can be imaged out to |${\sim } 50\hbox{-}100\,{\rm {\rm au}}$|⁠, then their overall geometry can serve as a test of this paper's hypothesis.

Unfortunately, at these distances, both 10 and 20 μm observations will probe the Wien tail of the dust's thermal emission. For WD0439+466 (the 129 pc system) the 10 and 20 μm values of B_ν are 4.9 × 10⁻³ and 1.4 × 10⁻¹ of the blackbody peak, respectively. For the Helix nebula (at 210 pc), these fluxes are reduced from the blackbody peak by factors of 3.6 × 10⁻³ and 1.6 × 10⁻³, respectively. Observations at 20 μm are therefore much more favourable.

The much greater spatial resolution of the atacama large millimeter array (ALMA) offers a possibly superior pathway for resolving the geometry of dusty IR sources around young WDs. However, the SEDs we plot in Figs 8–10 are only accurate for wavelengths λ ≲ 100 μm. A more accurate model for dust emission and absorption would be needed to accurately assess the observability of these debris clouds at ALMA wavelengths.

Figure 10.

Open in new tabDownload slide

Same as Fig. 8, but calculated for two different values for the strength of the WD natal kick v_k = 0.5 km s⁻¹ (dashed coloured lines) and 2 km s⁻¹ (dot–dashed coloured lines). The colours, solid coloured curves (fiducial model), and black sensitivity curves have the same meaning as in Fig. 8.

5.2 Implications for WD natal kicks

A variety of observations provide circumstantial, if not conclusive, evidence for WD birth kicks. We briefly summarize these motivating observations, before discussing how observations of WD debris discs can substantiate the existence of natal kicks.

Young WDs in globular clusters may possess a more extended radial distribution than their older counterparts. Such an effect was claimed in NGC 6397 (Davis et al. 2008), although recent HST observations have cast doubt on this conclusion (Heyl et al. 2012). Weaker evidence exists for an extended young WD radial distribution in ω Cen (Calamida et al. 2008). If correct, these results are counterintuitive: WD progenitors are relatively massive stars expected to migrate to the centre of their host cluster due to mass segregation, while WDs themselves are undermassive and hence should migrate back out of the centre after forming. Detailed dynamical modelling has found that such observations can be reproduced if a sizeable fraction of WDs receive birth kicks v_k ∼ σ (Heyl 2007a), where σ ∼ few km s⁻¹ is the cluster velocity dispersion. Fregeau et al. (2009), for example, find that ∼4 km s⁻¹ kicks could reproduce observations of NGC 6397. Importantly, the required kicks must also be impulsive with respect to the orbital time-scale of the WD about the globular cluster, i.e. to occur on a time-scale ≲ 10⁵ yr. Note that this condition is similar to that requiring the kick be impulsive with respect to the OCA (2), one of the key assumptions of our model.

Open clusters, with velocity dispersions lower than those in globular clusters, offer in some ways a more promising arena to investigate weak WD natal kicks. Open clusters are observed to possess a large depletion of WDs relative to their expected abundance from the IMF (Weidemann et al. 1992), which birth kicks v_k ∼ 1–5 km s⁻¹ could resolve (Fellhauer et al. 2003). This evidence is not conclusive, however, since alternative explanations for such WD deficits exist, most notably the hiding of WDs in binary systems.¹¹

WD kicks may also help explain the large orbital eccentricies of WD-barium star binary systems, which are significantly higher than those expected from tidal circularization of the orbit during the AGB phase of the WD progenitor (Izzard, Dermine & Church 2010). Asymmetric mass-loss and velocity kicks during the AGB superwind phase provides a natural explanation for such high eccentricities (Heyl 2007b). On the other hand, barium stars as a whole possess a large (≈99 per cent) binarity fraction, implying the kick velocity cannot be too high enough (v_k ≲ 2–3 km s⁻¹) or a large fraction of such systems would be dissociated (e.g. Izzard et al. 2010).

Finally, WD natal kicks are a possible explanation for the short birth periods of WDs (Spruit 1998), which are otherwise difficult to explain given the expected magnetic coupling between the WD-forming core and the envelope of the progenitor red giant.

Given the critical role that WD kicks play in producing observational signatures of cometary sublimation (e.g. IR emission detected by JWST), a confirmation of our interpretation of such observations would likewise provide strong evidence for kick velocities in the range 0.1 km s⁻¹ ≲ v_k ≲ 2 km s⁻¹. As we now describe, observations of older, hydrogen-deficient WDs provide an independent (and possibly contradictory) probe of WD natal kicks.

5.3 Hydrogen abundances in old WDs

Roughly half of the DB WDs in the Bergeron et al. (2011) sample¹² have so little atmospheric hydrogen that their observed spectra are incompatible with accretion of a single modestly sized comet, and are many orders of magnitude below our estimates for total accreted M_H. Even most of the DBA atmospheres in Bergeron et al. (2011) have hydrogen abundances typically one or more orders of magnitude below our fiducial M_H ∼ 5 × 10²³ g estimates.

These observed hydrogen deficiencies could, theoretically, be explained in the following ways.

• Convective dilution of a small surface hydrogen layer once the WD cools below a certain temperature.

• Extremely late-time surface burning of hydrogen.

• Loss of surface hydrogen in a line-driven wind.

• An unfavourably low OCA mass.

• A significantly small, large, or non-impulsive WD natal kick.

Although option (i) is widely believed to explain the observed DB deficiency among older, cooler WDs (Fontaine & Wesemael 1987), it manifestly cannot explain the observed hydrogen absences in the younger, hotter DB population discovered more recently (Eisenstein et al. 2006). Option (ii), the so-called born again scenario held to be responsible for Sakurai's Object (Duerbeck et al. 2000; Herwig et al. 2011), invokes a Very Last Thermal Pulse (VLTP) of hydrogen burning after the star has left the AGB branch, briefly reinflating the star to a large size. However, this VLTP is short-lived, and is generally taken to occur ∼10⁵ yr after the star leaves the AGB track (Iben et al. 1983). Although a H-burning pulse at these times could wipe out the hydrogen accreted from many, perhaps most, of the comets in our MC samples, a fraction of these comets come from semimajor axes |$a \gtrsim 10^4 \,{\rm {\rm au}}$| and therefore do not arrive at the WD until ∼10⁶ yr after the natal kick. Given the enormous value of M_H relative to Bergeron et al. (2011) limits, option (ii) seems unlikely. Option (iii) is also not promising; the high surface gravity of WDs poses a significant barrier to launching line-driven winds, and even if they can be driven at low mass-loss rates, they will not collisionally couple to hydrogen and helium (Unglaub 2008).

The remaining possibilities, options (iv) and (v), are of much astrophysical interest, as they touch on important open questions: the existence of OCAs and WD natal kicks. The formation of the Solar Oort Cloud depended critically on the existence of Solar System gas giants; without a properly configured gas giant system to scatter planetesimals, it is possible that a subset of exoplanetary systems would lack OCAs entirely, or at least see serious reductions in their mass (see Section 5.4 for further discussion).

Option (v) also has some plausibility: impulsive natal kicks v_k ≳ 4 km s⁻¹ are found to significantly reduce M_evap (Fig. 1). However, we again note that an enormous reduction in M_evap, much larger than is seen for a 4 km s⁻¹ kick, is needed to explain some of the observed DB hydrogen fractions. Furthermore, kicks much larger than this value are significantly constrained by observations of WDs in globular clusters (Section 5.2). A very small natal kick, v_k ≪ 0.1 km s⁻¹, could also shut off the hydrogen delivery mechanism proposed in this paper; a non-impulsive natal kick of larger magnitude would have the same effect. Although some combination of these possibilities appears to be the cleanest explanation of the low DB hydrogen fractions, we note that a prediction in this cases is the survival of a massive OCA. Over much longer time-scales, perturbations to this OCA by the galactic tide or passing stars could lead to the accretion of a handful of comets and the deposition of significant M_H (Wiegert & Tremaine 1999); however, calculating typical rates of comet accretion from these mechanisms is beyond the scope of this paper.

An important point raised in Section 5.2 concerns the survival of OCAs in dense stellar environments. Because WDs residing in globular or sufficiently dense open clusters will lose their OCAs over time, the floor for atmospheric hydrogen abundances will not be set by cometary accretion but instead by ISM accretion. A key prediction of our model is therefore that the ratio of extremely H-deficient DBs to the overall DB population should be greater in globular (and perhaps dense open) clusters. This argument is strengthened by the low metallicities of globular cluster stars, which likely impeded or prevented the formation of OCAs.

5.4 Implications for exoplanetary systems

The production of OCAs requires relatively massive planets orbiting at substantial distances from their host stars. The schematic picture we present in Appendix is largely based on Tremaine (1993) and describes how pericentric interactions with an interior planet will perturb the orbital energies of a coplanar planetesimal belt. Following a gradual process of energy diffusion, orbital evolution is taken over by the galactic tide, which leads to angular momentum diffusion at fixed energy, and isotropization of cometary inclinations. One notable feature of this model is that it requires a relatively restricted range of planetary parameters (mass M_xp and semimajor axis a_xp) in order to form an OCA. A sizeable (at least super-Earth in mass) exoplanet must reside at a semimajor axis |${\gtrsim } 1\,{\rm {\rm au}}$|⁠. There is therefore a reasonable chance that any young WD with evidence of a surrounding OCA would possess an exoplanet that survived its host star's post-MS evolution, particularly if the progenitor star was not too massive (Villaver & Livio 2007).

Fig. 11 shows the allowed parameter space of exoplanetary systems that will allow the production of an OCA around a solar mass star. Generally, a Neptune- or Uranus-like exoplanet is required, although at small semimajor axis, super-Earths may suffice. The allowed mass range vanishes for |$a_{\rm xp} \lesssim 1\,{\rm {\rm au}}$|⁠, as planetesimals then collide with the perturbing planet before the galactic tide can begin to alter their pericentre. Stars with a greater zero-age MS mass require larger planets; a 2.5 M_⊙ star requires Uranus- to Saturn-sized exoplanets to form an OCA.

Figure 11.

Open in new tabDownload slide

Parameter space of exoplanet properties (mass M_xp and semimajor axis a_xp) conducive to OCA formation (see Appendix for details). Solid lines represent the minimum (red) and maximum (blue) planet mass necessary to perturb a belt of planetesimals into a long-lived OCA, calculated for a solar-type star of mass M_* = M_⊙ and age t_* = 10¹⁰ yr embedded in a region of stellar density ρ_* = 0.15 M_⊙ pc⁻³ similar to that of the Solar neighbourhood. The region of parameter space capable of producing an OCA is shaded for this fiducial scenario. The dashed lines show an otherwise identical model calculated for a a more massive 2.5 M_⊙ star. Dotted and dot–dashed lines show the maximum mass calculated for a  M_⊙ star in regions of lower (ρ_* = 0.015 M_⊙ pc⁻³) and higher (ρ_* = 1.5 M_⊙ pc⁻³) ambient stellar density, respectively, than the fiducial case (the lower bound on exoplanetary mass is unchanged from the fiducial model shown with a solid line).

If our interpretation of the infrared excesses around young WDs is correct, then the ∼15 per cent detection rate of 24 μm excess (Chu et al. 2011) implies an equally large fraction of exoplanetary systems possess a middle to large-weight gas giant with a ≳ few au. WD debris discs may thus provide a unique probe of such objects, which are otherwise challenging to detect with radial velocity or transit surveys. This number is broadly compatible with inferred abundances of large-separation planets from microlensing surveys (Cassan et al. 2012), which are more sensitive to this type of exoplanet. However, given the existence of competing explanations (both dynamically active KBO analogues, and AGB wind fallback discs) for the IR excesses discussed in this paper, more followup work is needed before firm conclusions can be drawn about surviving exoplanets in such systems.

6 CONCLUSIONS

We have examined the response of extrasolar OCAs to the kicks that may accompany WD birth. Under the assumption that WDs receive a modest, ∼1 km s⁻¹ natal kick over a sufficiently short time-scale (≲10⁵ yr), then ∼10⁻³ of the mass of the associated OCA will be placed on to nearly radial orbits resulting in sublimation near pericentre. This has several important observational implications, which we list here.

• The solid debris from the sublimation process will form a roughly spherical, optically thin cloud around the WD, with dust grains absorbing and re-radiating a fraction of the WD's luminosity. The resulting IR emission (Figs 7, 8, and 10) can be fit into rough agreement with the 24 μm IR excesses observed by Spitzer around ∼15 per cent of newborn WDs. For roughly half the sample this can be accomplished with OCA masses comparable to the mass of the Solar Oort Cloud (M_OCA ∼ 10 M_⊕), but for the other half of the sample values of M_OCA > 100 M_⊕ are needed.

• Roughly half of young WD debris discs are seen to possess a second IR excess at much shorter ( ∼ 3 μm) wavelengths, which is not well fit by our model. The simple model presented in this paper is a time-independent one that excises all dust removed from the system by PR drag; we plan to investigate this steady-state, small-radius PR flow in future work to determine whether it could produce a second IR peak. It is also possible that the systems which exhibit the short-wavelength IR peak acquire both dust excesses through an alternate mechanism, e.g. an unresolved stellar binary companion.

• Our proposed mechanism for producing the observed IR excess differs from past hypotheses primarily in its geometric structure: spherical, rather than disc like. Our mechanism also employs a progenitor population residing far enough from the WD to resist sublimation during post-MS evolution; more tightly bound populations of planetesimals (e.g. Kuiper Belt analogues) are likely to sublimate before the PN phase.

• Young WDs with evidence for OCAs are also likely to possess sizeable planets orbiting at large enough semimajor axis to have survived post-MS evolution (Fig. 11). These WDs could be attractive targets for planet searches.

• The very late heavy bombardment predicted in this paper is at odds with the hydrogen fraction in extremely H-deficient DB WDs. Many DB WDs have less hydrogen in their atmospheres than would be delivered by the accretion of a single moderately sized comet, and therefore must have either (i) received a birth kick v_k ≪ 1 km s⁻¹ or v_k ≫ 1 km s⁻¹; or (ii) not have possessed an OCA.

Although this exploratory study has described the basic features of debris clouds produced around young WDs by natal kicks, several important theoretical problems remain for future work. These include time-dependent models for dust inspiraling due to PR drag; a realistic model for frequency-dependent dust absorption and emission; and a dynamical study of comet accretion on to WDs in the absence of a natal kick, which will test whether the extremely H-deficient DB WDs can only be explained by option (ii): a complete absence of OCAs.

For the nearest young WDs with IR excesses, it is possible that future imaging by JWST will provide a smoking gun test for this hypothesis, either by detecting the central cavity (due to PR drag) characteristic of disc scenarios, or by directly measuring the shape of the emitting area at larger sizes (true discs will appear non-circular due to inclination with respect to our line of sight). Careful spectroscopy may also be able to distinguish between these hypotheses. However, even with current Spitzer observations, the OCA mechanism proposed in this paper has one strong advantage over KBO scenarios: its cometary reservoirs of dust and gas survive the post-MS evolution of the parent star, which is not obviously the case for Kuiper Belt analogues. Current and future observations of debris discs around young WDs may therefore offer a rare opportunity to probe distributions of comets in extrasolar planetary systems.

We would like to thank John Debes, Jay Farihi, John Fregeau, Mike Jura, Chris Matzner, Norm Murray, Roman Rafikov, and Yanqin Wu for helpful discussions. We would also like to thank the anonymous referee for many useful suggestions that have improved this paper. AL was supported in part by NSF grant AST-1312034. BDM gratefully acknowledges support from the NSF grant AST-1410950 and the Alfred P. Sloan Foundation.

REFERENCES

APPENDIX A: OCA FORMATION

Taken together, equations (A3), (A9), and especially (A10) provide lower limits on the needed exoplanetary mass to form an OCA, while equations (A6), (A8), and (A12) provide corresponding upper limits. These constraints are plotted for a fiducial scenario in Fig. A1; in general it is equations (A3), (A6), and (A12), and at small a_xp, equation (A10) that most tightly regulate the mass of the required exoplanet.

Figure A1.

Open in new tabDownload slide

Constraints on the formation of OCAs by exoplanets of mass M_xp and orbital semimajor axis a_xp. In this fiducial case, we take a solar-type star in the solar neighbourhood of the galaxy, with M_* = M_⊙ and ρ_* = 0.15 M_⊙ pc⁻³; we plot these constraints at a time t_* = 5 × 10⁹ yr after its formation. We furthermore take the exoplanet's density ρ_xp = 3 g cm⁻³, and consider a belt of planetesimals of initial thickness Δθ = 0.1. The thin red curves represent lower bounds on the required planetary mass, while the thick blue curves represent upper bounds. The solid thin curve, solid thick curve, thick dashed curve, thin dashed curve, thin dotted curve, and thick dotted curve denote equations (A3), (A6), (A8), (A9), (A10), and (A12), respectively. The most constraining requirements are generally that orbital energy diffusion can happen in a time less than t_*, that energy diffusion does not lead to the escape of most comets, and (at small a_xp) that the comets do not directly impact the perturbing planet before beginning angular momentum diffusion.

In this section, we have neglected certain additional physical effects – for example, radial migration of a star's orbit through the galactic potential (Kaib et al. 2011), or dense stellar birth environments (Brasser et al. 2006; Kaib & Quinn 2008) – in order to provide a simple analytic picture. In general, numerical orbital integrations will be necessary to fully model the birth and evolution of OCA candidates.

APPENDIX B: IMPULSIVE mass-loss

APPENDIX C: EVAPORATIVE mass-loss OVER MULTIPLE ORBITS

APPENDIX D: ECCENTRIC PR DRAG

The solid-state debris left over after sublimation of the volatiles in a comet will travel on highly eccentric orbits, which are subject to PR drag from the central WD of luminosity L_WD. Here, we outline the dynamics of PR drag on dust grains of size b for general values of the orbital eccentricity.

Figure D1.

Open in new tabDownload slide

PR drag time-scale, t_PR, as a function of the eccentricity deficit 1 − e of the initial orbit, calculated for a WD of luminosity L_WD = 100 L_⊙ and for orbits with different pericentre radii |$r_{\rm p}=50\,{\rm {\rm au}}$| (black) and 1 au (red). Solid, dashed, and dotted lines are calculated for particles of size b = 10⁴, 10², and 1 μm, respectively. Note that t_PR represents only the ‘initial’ drag time-scale, which will in general increase (⁠|$\propto L_{\rm WD}^{-1} \propto t^{1.2}$|⁠, approximately; equation 8) as the WD cools, such that the actual drag time becomes infinite once t_PR exceeds the WD cooling time-scale. The effect of eccentricity for values 1 − e ∼ 10⁻¹–10⁻² (typical of the OCA debris streams of interest here) is to increase the drag time-scale by up to an order of magnitude as compared to the circular case (right edge of the plot).