© The Authors 2026

1 Introduction

To date, over 6000 exoplanets have been discovered around main-sequence stars. Most solar-type stars in the Milky Way host at least one planet on average (Zhu & Dong 2021). As most known host stars of a planetary system end their life as white dwarfs (WDs), these post-main-sequence stars offer a unique window for understanding the composition and evolution of exoplanetary systems (McDonald & Veras 2023).

Owing to the high surface gravities of WDs, heavy elements (metals) undergo rapid diffusive settling in their photospheres on timescales that are short compared with typical WD cooling ages (Paquette et al. 1986). Consequently, metals accreted from external reservoirs are readily detectable against otherwise pure hydrogen- or helium-dominated WD atmospheres; this phenomenon is referred to as WD pollution (Veras 2021). Spectroscopic surveys show that approximately 25-50% of the surveyed WDs display signatures of metal pollution (Zuckerman et al. 2003; Zuckerman et al. 2010; Koester et al. 2014; Coutu et al. 2019), which implies that a substantial fraction of WDs are currently accreting external material.

Composition analyses of polluted WDs indicate that the abundances of accreted elements are broadly consistent with the bulk compositions of Solar System asteroids and comets (Jura 2006; Jura et al. 2009; Zuckerman et al. 2007, 2011; Swan et al. 2019; Swan et al. 2023; Xu et al. 2017, e.g.,). These planetesimal-derived pollutants are driven within the Roche sphere of the WD, where they undergo tidal disruption, and their fragments are subsequently accreted onto the stellar surface (Jura 2003; Debes et al. 2012; Chen et al. 2019, e.g.,). The delivery mechanisms that bring material to this state strongly depend on size: Small particles (sizes of ~10⁻⁶−10³m) are transported inward primarily by radiative forces such as Poynting-Robertson (PR) drag and the Yarkovsky effect (Veras et al. 2022), whereas large bodies (~10³−10⁶m) must be driven to very high eccentricities (this is typically done through planetary perturbations or interactions with stellar/companion bodies) in order to approach the Roche radius (Debes & Sigurdsson 2002; Frewen & Hansen 2014; McDonald & Veras 2023; Veras & Rosengren 2023, e.g.,).

Planets orbiting WDs are widely invoked as dynamical agents that can deliver small bodies into star-grazing orbits and thereby contribute to photospheric pollution (Frewen & Hansen 2014; Veras & Rosengren 2023, e.g.,). Nevertheless, only a handful of planetary mass objects have been confirmed in post-main-sequence systems to date (e.g., PSR B1620-26(AB)b: Thorsett et al. (1993); Sigurdsson et al. (2003), WD 0806-661b: Luhman et al. (2011), WD J0914+1914b: Gänsicke et al. (2019), WD 1856+534b: Vanderburg et al. (2020), MOA-2010-BLG-477Lb: Blackman et al. (2021), and PHL 5038AB: Casewell et al. (2024), with several additional candidate detections recently reported from James Webb Space Telescope (JWST) (Mullally et al. 2024; Limbach et al. 2024; Debes et al. 2025). While WDs with planetary systems are generally more susceptible to pollution, the well-studied system WD 1856+534 (hosting a close-in giant at ~0.02 au) shows no detectable photospheric metal pollution despite the presence of massive short-period companions (Vanderburg et al. 2020; Xu et al. 2021). These observations motivate the hypothesis that sufficiently massive close-in planets act as a dynamical barrier, either by ejecting, capturing, or redirecting incoming bodies, thereby suppressing the delivery of pollutant material to the WD photosphere.

Previous studies have examined multiple dynamical mechanisms by which planetary companions can reduce WD photospheric pollution originating from diverse sources. For meter-to kilometer-scale particles, Veras (2020) showed that close-in giant planets can impede inward migration of radiatively driven material, such as that caused by Poynting-Robertson drag and the Yarkovsky effect, by trapping particles into mean-motion resonances. Regarding cometary contributions, O’Connor et al. (2023) and Pham & Rein (2024) estimated that planets on wide orbits (~10-100 au) can reduce cometary delivery rates through a combination of direct scattering and companion-induced secular precession. Despite these advances, an important gap persists concerning the planet’s suppressive effect on the dominant pollution channel: dynamically excited, high-eccentricity asteroids that ultimately intersect the WD Roche radius (Fröhlich & Regály 2024). While specific case studies (e.g., Rogers et al. 2024) have shown that individual planets can eject infalling bodies, a systematic quantification of whether these gravitational interactions statistically reduce the overall pollution rate across a broad planetary parameter space remains absent.

We perform short-term three-body simulations to quantify how a close-in planet affects near-parabolic asteroids¹. We then statistically evaluate the efficiency with which various planetary companions inhibit the injection of pollutant material into the WD, that is, their effectiveness as dynamical protectors of the host star. In Section 2 we describe our physical model and numerical setup, specify the ranges of planetary and asteroid parameters we explored, and define the quantitative criteria we used to measure planetary protection. Section 3 presents representative case studies and a statistical analysis of population outcomes for different planetary masses and orbital separations. In Section 4 we examine several factors that affect our results, discuss limitations and caveats, and compare our findings with previous studies. Finally, we summarize our main conclusions and implications in Section 5.

2 Methods

2.1 Simulation setup

We considered a system in which a planet moves on a circular orbit around a central WD, while an asteroid approaches the star on a highly eccentric trajectory. Because the mass of the asteroid is negligible relative to the planet, we treated the asteroid as a massless test particle throughout the simulations. The schematic model is shown in Fig. 1.

We set the central WD mass to M_WD = 0.518 M_⊙, adopting the best-fit mass of WD 1856+534 reported by Vanderburg et al. (2020). We also adopted the corresponding WD radius of 0.0131 R_⊙ to identify asteroids that collided with the WD. Although multiple studies reported marginally different stellar parameters for this system (Xu et al. 2021; Alonso et al. 2021; Limbach et al. 2025), tests we performed indicated that variations within the quoted uncertainty in M_WD did not produce statistically significant changes in our results (see Section 4.1 for detailed statistics).

The mass of WD1856+534b remains observationally unconstrained; theoretical models of brown dwarf evolution and atmospheres provide an upper limit of 13.8 M_Jup (95% confidence) (Vanderburg et al. 2020). More recent JWST observations (Limbach et al. 2025) constrained the upper limit to approximately 6 M_jup. For the primary parameter study, we adopted the more conservative upper limit 13.8 M_Jup in order to quantify the maximum potential protective effect of a planet; we then explored a broad range of planetary masses from one Earth mass to 14 Jupiter masses to assess how mass affects our results. We set the planet semimajor axis (a_p) to 0.02 au for the benchmark case and additionally selected several values of a_p between 0.02 au and 1 au to study how the orbital separation affects our outcomes. The results of the parameter study appear in Section 3.4.

Because the physical radius of a giant planet is significantly larger than that of a WD, asteroid-planet collisions represent a critical sink for potential polluters. We defined a collision event whenever the distance between a particle and a massive body fell below its physical radius. For our benchmark case, we adopted the planet radius R_p = 10.4 R_⊕ based on the observed radius of WD 1856+534b (Vanderburg et al. 2020). For other planets in our explored mass range (1 M_⊕−14 M_Jup), we estimated the physical radius using the empirical mass-radius relations from Müller et al. (2024), $$R_p=\{\begin{array}{ll}(1.02\pm 0.03)M_p^{(0.27\pm 0.04)},& M_p<(4.37\pm 0.72),\\ (0.56\pm 0.03)M_p^{(0.67\pm 0.05)},& (4.37\pm 0.72)<M_p<(127\pm 7),\\ (18.6\pm 6.7)M_p^{(-0.06\pm 0.07)},& M_p>(127\pm 7),\end{array}$$(1)

where R_p and M_p are expressed in Earth units (R_⊕ and M_⊕).

We sampled asteroid semimajor axes a_ast uniformly within [1, 10] au. This range was chosen to represent the inner-system reservoir where the majority of rocky polluters are expected to originate, while ensuring a_ast > a_p for all planetary configurations (a_p ≤ 1 au), thereby guaranteeing that the initial asteroid orbits can cross that of the planet through subsequent eccentricity excitation. We note that the dynamical effect of the asteroid initial semimajor axis depends on the planet properties. For massive close-in planets, this 1-10 au range typically sits above the threshold a_in (see Section 4.1). For lower-mass or more distant planets, however, a_in can shift toward higher values, meaning that their shielding effectiveness might be more sensitive to the specific spatial distribution of the asteroid belt. By adopting a fixed 1-10 au range, we provide a consistent benchmark to compare the protective capacity across the entire planetary parameter space, while acknowledging that the exact protection fractions for the least massive planets might vary when more distant reservoirs are considered. More realistic asteroid semimajor axis distributions can be set by specific preexcitation mechanisms, which are beyond the scope of this work.

We assumed that each asteroid initially possessed a perihelion within the WD Roche limit. Specifically, we sampled the initial pericenter q_ast uniformly in the range (0.5, 1.0) r_roche. This ensured that the asteroids were on near-disruptive trajectories while providing sufficient variation to sample planet-induced gravitational deflections. The Roche radius is defined as $$r_{\mathrm{R}\mathrm{o}\mathrm{c}\mathrm{h}\mathrm{e}}={\left(\frac{3M_{\mathrm{W}\mathrm{D}}}{2\pi {\rho }_{\mathrm{a}\mathrm{s}\mathrm{t}}}\right)}^{1/3},$$(2)

where we adopted an average asteroid density of ρ_ast = 2 g/cm³ (Carry 2012), yielding r_Roche ~ 0.004 au. To further explore the sensitivity of our results, we also considered a population of asteroids perturbed onto orbits just exterior to the Roche limit (Li et al. 2025) by conducting additional tests with q_ast sampled in (1, 2) r_Roche.

In each simulation, the asteroid was initialized at its apoastron (true anomaly f = π). The argument of pericenter ω was drawn uniformly at random from [0,2π). Inclination i and longitude of the ascending node Ω were sampled in two separate suites to explore both coplanar and noncoplanar configurations: in the coplanar suite, we set i = Ω = 0; in the noncoplanar suite, we drew i and Ω randomly from [0, π] and [0,2π), respectively. The dynamical differences between these configurations are discussed in Section 3.2.

We set the simulation time in units of the initial asteroid orbital period T_ast, determined by the initial semimajor axis a_ast of the asteroid. An ideal mechanism for a planet to prevent pollution is to eject the asteroid before it is completely disrupted. Crucially, an asteroid entering the Roche limit on a highly eccentric orbit is expected to undergo significant tidal fragmentation and complete disintegration within approximately ten periapsis passages (Malamud & Perets 2020). Since we treated the asteroid as a single intact body to track its orbital evolution, we limited our primary simulations to 10 T_ast. For the majority of the cases we studied, most ejections and collisions occurred within these first few orbits. While some planets might eventually eject asteroids over much longer timescales (as discussed in Section 3.2), bodies that remained bound and within the disruption zone after 10 T_ast were assumed to be totally disrupted. Beyond this point, further orbital evolution would be dominated by the dynamics of debris streams and nongravitational effects such as Poynting-Robertson drag (e.g., Li et al. 2021), which are beyond the scope of our current N-body model.

Our simulations primarily focused on resolving the dynamics during periapsis. Although the Roche region occupies only a small fraction of the total orbit, the highly eccentric nature of the trajectories demands exceptional temporal resolution. This necessitates substantial computational resources, as extremely small time steps are required to maintain precision during close encounters. We conducted N-body simulations using the REBOUND package (Rein & Liu 2012) and integrated them with the IAS15 integrator (Rein & Spiegel 2015), which uses a high-order adaptive time-stepping scheme. The integrator therefore automatically shortens time steps during periapsis passages of high-eccentricity trajectories to resolve close approaches with high fidelity while maintaining a rigorous constraint on global integration errors throughout the simulation.

Fig. 1 Simplified model of the WD system. The planet orbits the central WD on a circular path, while an asteroid initially follows a highly eccentric trajectory toward the WD. The figure is not to scale.

2.2 Definition of protection

We ran each set of simulations until one of the following stopping criteria was met: (1) the asteroid was either ejected from the inner system or collided with one of the massive bodies (the WD or the planet); or (2) the simulation time reached the upper limit of 10 T_ast. When multiple criteria might be satisfied simultaneously, the earliest interaction was recorded. Each asteroid was then classified into one of five mutually exclusive outcome classes based on three diagnostic variables: (i) the minimum distance between the asteroid and the WD, d_min; (ii) the final orbital energy of the asteroid, E_final (negative values indicate bound orbits); and (iii) the minimum distance between the asteroid and the planet, d_p,min. The outcome classes are defined as listed below.

Ejected-OutRoche (Ejected-OR) (E_final > 0, d_min > r_Roche, d_p,min > R_p): ejected before entering the Roche radius and without collision with the planet.

Ejected-InRoche (Ejected-IR) (E_final > 0, R_WD < d_min < r_Roche, d_p,min > R_p): ejected after passing inside the Roche radius, but prior to any impact.

Collide (d_p,min < R_p): collides with the close-in planet.

Crash (d_min < R_WD, d_p,min > R_p): crash into the WD.

Bound (E_final < 0, d_min > R_WD, d_p,min > R_p): remains gravitationally bound to the WD at the end of the integration and does not collide with the planet or the WD.

Without planets in close-in orbits, all bodies whose initial orbits cross the Roche limit would pollute the WD. The existence of planets will change the asteroid fate, and we quantified the protective effect of planets with the proportion of asteroids that were finally ejected at the end of the simulation and that collided with the planet. We defined the protective fraction $$F_{\mathrm{P}\mathrm{r}\mathrm{o}\mathrm{t}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}}=F_{\mathrm{E}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}-\mathrm{O}\mathrm{R}}+F_{\mathrm{E}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}-\mathrm{I}\mathrm{R}}+F_{\mathrm{C}\mathrm{o}\mathrm{l}\mathrm{l}\mathrm{i}\mathrm{d}\mathrm{e}},$$(3)

where F_Ejected-OR and F_Ejected-IR are the fractions of asteroids in the Ejected-OR and Ejected-IR classes, respectively, and F_Collide is the fraction in the Collide class.

3 Results

3.1 Typical cases

In the previous section, we established criteria to assess the planetary effect on small bodies and classified the outcomes into five distinct classes. Below, we analyze each class and present representative examples.

Fig. 2a illustrates the first class (Ejected-OR), in which planetary perturbations deflect asteroid trajectories before they reach the Roche sphere. This deflection confines these asteroids exterior to the planet orbit, preventing close approaches to the WD and representing the strongest form of planetary protection.

In the second outcome class (Ejected-IR), asteroids can be expelled either by a single close planetary encounter or via multiple encounters that progressively increase their orbital energy, ultimately producing unbound trajectories. Representative examples are shown in Figs. 2b,c. Although classified as the same category, asteroids with different initial conditions exist within the Roche limit of a WD for significantly varying durations. Because tidal disruption is not instantaneous upon entering the outer Roche limit (Malamud & Perets 2020), these residencetime differences produce divergent evolutionary outcomes. Bodies with short Roche-limit transits (∆t ≲ τ_disrupt) might be ejected before significant tidal fragmentation, thereby avoiding accretion onto the WD. Conversely, prolonged residence (∆t ≫ τ_disrupt) leads to full tidal disintegration, producing debris that can ultimately pollute the WD atmosphere. This outcome represents an intermediate imperfect protective effect (termed weak protection) in which the planet reduces but does not eliminate the risk of WD pollution. The population of this category can provide a wide reference constraint on the protection effect of planets. We address a critical extension below in Section 3.3: whether planets can effectively clear debris prior to accretion when asteroids are shattered inside the Roche limit.

When an asteroid trajectory is strongly perturbed by planetary gravity, it can undergo a close encounter that leads to a direct impact with the planet. Fig. 2d illustrates a typical example of such a collision. The orbit shows a sudden termination at the planetary radius. Once captured by the planetary atmosphere, these pollutants are effectively removed from the system, preventing any further interaction with the WD. This direct interception constitutes an additional efficient mechanism for planetary protection.

Planetary scattering does not always redirect small bodies outward; in some cases, it instead promotes inward scattering that can lead to WD accretion. Fig. 2e shows a representative trajectory in which a small body is scattered on a nearly radial path toward the WD. We did not model post-Roche processes (e.g., tidal fragmentation, collision cascades, and subsequent debris evolution) in our simulations. However, including tidal-disruption physics might alter this outcome: fragments produced inside the Roche limit might still be cleared by planetary perturbations before accretion, thereby reducing the pollution. We present the results below in Section 3.3.

Fig. 2f presents a representative bound case (E_final < 0). These orbits often arise from planetary perturbations that significantly increase orbital periods, thereby delaying the asteroid return to the WD vicinity. Unless a planet ejects a body in a single scattering event, the object typically undergoes multiple close encounters, each of which can lead to either ejection or inward scattering. While the cumulative probability of ejection increases with each encounter (O’Connor et al. 2023), asteroids in systems with low-mass planets can persist for thousands of years, potentially crossing the Roche limit numerous times. Since an asteroid might be tidally disrupted and accreted before ejection occurs, we classified these temporarily bound objects as a separate category and excluded them from the protection fraction calculation to ensure a conservative estimate of the defensive role of the planet.

Fig. 2 Typical asteroid trajectories of four classifications in Section 2.2 on a logarithmic scale, with units in au. The line colors denote classifications, and the diamond markers indicate the starting points of the asteroids. (a) Ejected-OR: asteroid was ejected before it entered the Roche limit. (b) and (c) Ejected-IR: both are ejected after entering the Roche limit, the difference is the duration they spent within the Roche limit. (d) Collide case: asteroid collided with the planet when it crossed the planetary orbit (we depict the planet at the moment of collision (navy circles); note that the planetary radii do not represent the actual planet sizes.). (e) Crash case: asteroid was scattered directly to the WD. (f) Bound: asteroids still orbit the WD until the end of simulation.

3.2 Population analysis

While the previous section categorized individual asteroid trajectories, this section provides a statistical quantification of these outcomes across a broad parameter space. We evaluated the effectiveness of planets as dynamical shields by analyzing the relative proportions of each outcome class. The comprehensive statistical results for various planetary mass and orbital configurations are summarized in Table 2.

3.2.1 Two dimensions

Although the Ejected-OR category represents the most effective form of planetary protection in which asteroids are removed before they can ever threaten the WD, it is not the dominant outcome. In our benchmark case (see Table 1), only 12% of the asteroids are deflected while remaining outside the Roche limit; these are visualized as the blue points within the black shaded region in the left panel of Fig. 3. The statistics in Table 2 reveal that at a fixed semimajor axis of 0.02 au, the fraction of asteroids achieving this fate decreases significantly as planetary mass declines, reflecting a reduced scattering cross section. Conversely, for a fixed mass of 13.8 M_Jup, the proportion of this category increases slightly as the planetary semimajor axis is extended beyond 0.1 au, exceeding 20%. This suggests that more distant planets, while having a smaller overall impact area, might more gently nudge asteroids into unbound orbits before they reach the inner system.

In contrast, the Ejected-IR class constitutes the most populous outcome, where asteroids are scattered out of the system only after penetrating the Roche limit. In our benchmark case, this class accounts for approximately 70% of the total population. This high fraction implies that while planets are efficient at eventually clearing eccentric asteroids, they often fail to prevent the initial penetration into the immediate WD vicinity. This proportion shows a nonmonotonic dependence on planetary mass: it peaks at nearly 80% as the mass decreases slightly from the benchmark, but subsequently declines for planets significantly smaller than one Jupiter mass, where the scattering efficiency drops. Furthermore, the fraction of this class decreases as the planetary semimajor axis increases because more distant planets exert a weaker dynamical effect on the inner regions near the Roche limit. Since this category includes asteroids that may have completed multiple periastron passages, we further quantify the planet effect on the resulting tidal fragments in Section 3.3.

Asteroids remaining bound to the system at the end of our simulations are the primary potential sources for WD pollution. For planets located at 0.02 au, we find that when the mass exceeds 5 M_Jup, the bound fraction remains below 2%, indicating a high clearing efficiency. However, F_Bound increases significantly for lighter or more distant planets. To illustrate the long-term clearing process, Fig. 4 shows the temporal evolution of the bound fraction for a 1 M_Jup planet at 0.1 au. The planet requires approximately 100 initial orbital periods (roughly 2000 years) to fully clear this population. This delay is critical, as it provides a window for asteroids to undergo tidal disruption or radiative circularization (Li et al. 2021), potentially leading to accretion before dynamical ejection can take place.

A small but distinct subset of asteroids undergoes a more violent dynamical evolution, leading to the Crash outcome. In our benchmark case, approximately 7% of the planet-affected asteroids are scattered directly onto the WD surface. This fraction is highly sensitive to the planetary parameters: at a fixed semimajor axis, it drops to 0% for a 1 M_Jup planet, but climbs to 20% as the planet distance increases beyond 0.5 au. While these direct impacts represent the most immediate form of pollution, they differ from the gradual accretion typically inferred from debris disks. A high Crash fraction suggests a regime in which planetary perturbations are strong enough to bypass tidal disruption entirely. If a reservoir of eccentric objects is continuously replenished, this direct channel (even at an efficiency of 10-20%) might contribute to the observed steady-state pollution, albeit with different temporal characteristics than the disk-fed accretion model.

Beyond gravitational scattering, the physical size of giant planets allows them to act as a direct sink for infalling debris. In our benchmark case, approximately 10% of the asteroids end their evolution by colliding with the planet. Interestingly, due to the mass-radius relation for giant planets (see Eq. (1)), planets with masses around 0.1-1 M_Jup possess larger physical cross sections than their more massive counterparts. Consequently, these intermediate-mass planets exhibit higher collision fractions, making them more effective at sequestering potential polluters. This suggests that the chemical signatures of pollution in some WD systems might be mirrored by the enrichment of their close-in giant planets, offering a potential planetary shield signature for future atmospheric characterization. However, this protective mechanism becomes inefficient as the planetary semimajor axis increases beyond 0.5 au, where the collision fraction falls below 5%. This decline highlights that for distant planets, protection is achieved almost entirely through dynamical ejection rather than physical interception.

For asteroids that do not initially cross the Roche limit, the presence of a close-in planet provides an additional layer of security. Our results show that the majority of these objects are ejected from the system before they can even reach the Roche zone. Specifically, in the benchmark case, only about a quarter of the asteroids are driven into the Roche limit by planetary perturbations. Of these, 17.7% are subsequently ejected after their inward excursion (Ejected-IR), and only a small fraction (5%) result in a direct Crash onto the WD. Notably, the collision fraction remains largely insensitive to the initial perihelion locations of the asteroids. This suggests that the planetary shield is not merely a local obstacle, but extends its effect well beyond the immediate tidal disruption zone, efficiently cleaning the system of potential polluters before any significant inward orbital migration can occur, such as that driven by radiative effects.

Fig. 3 Distribution of the final asteroid energy and minimum distance to the WD during the entire simulation. The system parameters selected for this simulation are the benchmark parameters of WD1856+534b system given in Table 1. The red dots in the figure represent the initial state of asteroids. The black dots represent the final state of each asteroid group at the end of the simulation; the color of the region in which they land denotes their ultimate fate. The inverted black triangles indicate those that collided with planets.

Fig. 4 Proportion of bound asteroids to all simulated asteroids as a function of simulation time. The two different colors represent the trends in different time units.

3.2.2 Three dimensions

Although planets and small bodies in the Solar System are predominantly coplanar, many exoplanetary systems host planets with high eccentricities and significant mutual inclinations, enabling pollutant delivery from a wide range of directions. In addition, contaminants originating in the outer reservoir (e.g., exo-Oort-cloud comets) can approach the inner system from arbitrary inclinations. To examine the difference between planetary effect under coplanar and noncoplanar circumstances, we randomly selected 1000 simulations and plotted the distribution of asteroids in the energy-distance plane (Fig. 3).

The probability of extremely close planet-asteroid encounters is higher for coplanar geometries than for three-dimensional configurations, producing stronger instantaneous perturbations in the coplanar case. As shown in Fig. 3, the distribution of outcomes in the energy-distance plane is broader for the coplanar sample than for the three-dimensional sample. Additionally, we quantified the difference in coplanar and noncoplanar scenarios by the average minimum gravitational interaction distance between the asteroid and the planet. In the noncoplanar case, this distance is 2.8 R_Hill, while in the coplanar case, it is 2.1 R_Hill. At this interaction distance, the effect of the planet to alter the fates of the asteroids in three dimensions remains significant, although the magnitude of individual changes in orbital direction decreases. Accordingly, the fraction of Ejected-OR objects (F_Ejected-OR) falls to about 8% in the three-dimensional experiments. Simultaneously, the Ejected-IR fraction (F_Ejected-IR) increases to roughly 88%, so the ejection fraction in three dimensions exceeds that found for coplanar runs. The fraction of direct collisions with either the WD or planets has fallen to about 1%. Consequently, the total protective fraction, F_Protection, rises to approximately 97% for isotropic pollutant sources.

This superior protective effect in three dimensions is highly sensitive to planetary parameters, however. While the shield remains stronger in 3D for slight increases in semimajor axis or small decreases in mass, it is not universally more effective than the coplanar shield. As the planetary mass continues to decrease or the semimajor axis increases, the effectiveness of the 3D protection degrades much more precipitously than its twodimensional counterpart. Consequently, the protection efficiency in the three-dimensional scenario suffers a far more dramatic reduction than in the two-dimensional case. For instance, when the planetary mass is reduced to 0.1 M_Jup at 0.02 au, the threedimensional protective effect nearly vanishes, dropping to only 5.1%, while the 2D protection remains robust at 82.2%.

3.3 Effect of fragmentation

In the preceding simulations, asteroids were treated as intact particles even after crossing the Roche limit. However, a more realistic and observationally motivated scenario is that these bodies undergo tidal disruption upon close approach, forming a debris stream or disk that eventually accretes onto the WD (Veras et al. 2024). Since the single-body approximation cannot capture the complex dynamical evolution that occurs within the Roche limit, we extended our model to incorporate tidal fragmentation. This allowed us to evaluate whether the divergent trajectories of fragments, under the combined effect of the WD and the planetary perturbations, significantly alter the overall protection efficiency.

Following the analytical framework of Malamud & Perets (2020), we adopted a simplified disruption model based on the following assumptions: (1) The parent asteroid was assumed to be spherical and dynamically undisturbed prior to entering the Roche limit. (2) Upon the first crossing of the Roche limit, the body was disrupted instantaneously into ten fragments, each of which was subsequently treated as a massless test particle. (3) These fragments inherited the instantaneous center-of-mass velocity of the parent body at the moment of disruption, while their initial positions were distributed isotropically at a distance r_frag from the parent center. In this prescription, the divergent evolution of fragments is driven solely by their spatial offsets, which translates into slight differences in their orbital energy and angular momentum relative to the WD.

The spatial offset r_frag, which scales with the physical asteroid radius, introduces a degree of dynamical diversity among fragments. Specifically, there exists a critical radius r_crit above which fragments with r_frag > r_crit become gravitationally unbound from the WD immediately upon disruption. Using Equation (2) from Malamud & Perets (2020), we calculated the critical value for our system $$r_{\mathrm{c}\mathrm{r}\mathrm{i}\mathrm{t}}=\frac{r_{\mathrm{R}\mathrm{o}\mathrm{c}\mathrm{h}\mathrm{e}}^2}{2a_{\mathrm{a}\mathrm{s}\mathrm{t}}-r_{\mathrm{R}\mathrm{o}\mathrm{c}\mathrm{h}\mathrm{e}}^2}\sim 127\mathrm{k}\mathrm{m}.$$(4)

In the absence of a planet, fragments that become unbound due to tidal forces alone would be unable to pollute the WD. Because our definition of planetary protection (F_Protection) is predicated on a 100% baseline pollution rate, however, we primarily stipulated r_frag < r_crit to ensure that any observed protection was attributable to planetary effects. Nevertheless, considering that realistic polluters might be larger (Veras et al. 2024), we also explored cases where r_frag > r_crit as a reference. Consequently, we tested r_frag values of 10, 50, 100, and 500 km, with the corresponding statistical outcomes presented in Table 2.

We conducted these fragmentation simulations using exactly the same parameters as the two-dimensional benchmark case to isolate the effect of tidal disruption on the resulting statistics.

Notably, the fraction of asteroids ejected prior to reaching the Roche limit (Ejected-OR) remained identical to the single-body model, as only those parent bodies that penetrate r_Roche were subject to fragmentation. Upon disruption, however, fragments from the same parent can follow markedly divergent dynamical trajectories. To illustrate this stochasticity, Fig. 5 presents a representative case comparing the evolution of a parent asteroid and its child fragments. In the single-body integration (upper panel), the parent body is directly accreted by the WD. In contrast, the lower panels reveal that while three out of ten fragments (30%) follow the parent fate and accrete onto the WD, the remaining seven fragments (70%) are eventually ejected from the system due to subsequent planetary perturbations. This example highlights how tidal disruption can transform a potential accretion event into a series of ejections, diversifying the possible outcomes for the original pollutant mass.

As discussed in Section 3.1, asteroids in the Ejected-IR category exhibit a wide diversity in their encounter histories prior to ejection. We refined our analysis of this variation by tracing the ultimate fates of individual fragments. By assigning each fragment a mass fraction of 0.1 and adopting a fragment-weighted accounting scheme, we determined the ejection and accretion fractions with higher precision than by treating parent bodies as indivisible units. For instance, as illustrated in Fig. 5, a parent body previously categorized as a Crash (count = 1) is now recorded as contributing 0.7 to the Ejected-IR fraction and 0.3 to the Crash fraction. This refined approach effectively captures how the stochastic history of close planetary encounters dictates the divergent evolutionary paths of an asteroid’s constituent debris.

Under this fragment-weighted accounting scheme, the overall protective fraction is F_protection ~ 91.5%. This value is remarkably consistent with the F_protection = 92.2% reported for the single-asteroid model (Table 2, first row), indicating that fragmentation does not fundamentally undermine the planetary shield. Crucially, the initial radius of the asteroid has a negligible effect on the final outcome. The insensitivity of the protection fraction to the asteroid radius (from 10 to 500 km) suggests that the phase space for ejection is dominated by the strong planetary scattering. The initial energy spread caused by tidal disruption even for a 500 km body is insufficient to nudge fragments out of the vast gravitational reach of the planet or to significantly alter the statistical probability of accretion versus ejection. Although this simplified fragment model does not capture all real-world complexities, it provides a useful first-order reference: Eq. (3) yields a reasonable estimate of the planetary protection efficiency.

Fig. 5 Evolution of the distance between the particles and the WD over time. The top panel shows the evolution of an asteroid without tidal disruption. The middle and bottom panels show the evolution of the fragments considering the tidal disruption, with different line styles representing different fragments.

3.4 Effect of planet mass and orbital parameters

In this subsection, we investigate the effect of planetary mass and semimajor axis on the protection efficiency, identifying the conditions under which planets provide significant pollution shielding. Our simulations sampled planetary masses M_p from 1 M_Earth to 14 M_Jup and semimajor axes α_p from 0.02 to 1 au. For each (M_p, a_p) pair, we performed 100 simulations to obtain outcome distributions analogous to those in Section 3.2. Table 2 summarizes representative results for selected (M_p, a_p) combinations, and Fig. 6 illustrates the dependence of F_protection on planetary mass for five fixed semimajor axes.

The key findings are as follows. First, protective efficiency scales positively with planetary mass, approaching a maximum of over 90%. For planets at 0.02 au, masses ≳1 M_Jup prevent more than 90% of potential asteroid incursions, demonstrating the formidable shielding power of close-in gas giants. Second, the overall shielding effect of the planet diminishes with increasing semimajor axis. At a_p = 0.5 au, the protective fraction drops by ≳20% across the explored mass range. Beyond a_p = 1 au, even the most massive planets in our sample eject only ~60% of asteroids. This declining trend with distance stems from a competition between scattering strength and encounter frequency. Although the ratio of the surface escape velocity of a planet to the local Keplerian velocity increases at larger distances (which theoretically favors ejection during a single close encounter), this is outweighed by the longer orbital periods and the significantly reduced frequency of encounters. Consequently, within a finite integration timescale, the cumulative probability of a successful ejection is lower for more distant planets, even when they are intrinsically stronger scatterers per event.

The direct scattering effect of a planet can be estimated analytically based on its mass and orbital parameters. If the kinetic energy imparted during a planetary encounter exceeds the particle gravitational binding energy, the particle becomes unbound and is prevented from accreting onto the WD (O’Connor et al. 2023). To evaluate this protective capability, we adopted a qualitative analytical approach. While O’Connor et al. (2023) primarily addressed perturbations from distant encounters, close-in planets are more likely to undergo strong, deep encounters with asteroids. Our simulations indicate that the typical closest-approach separation scales as approximately twice the Hill radius (2 R_Hill). By adopting this characteristic encounter distance, we modified the scattering parameter Λ (Eq. (39) of O’Connor et al. 2023) to define a close-encounter ratio Λ_c, $$\begin{array}{rl}{\mathrm{\Lambda }}_c& =\left(\frac{G{M}_p}{2R_{\mathrm{H}\mathrm{i}\mathrm{l}\mathrm{l}}}\right)/\left(\frac{G{M}_{\star }}{a_{\mathrm{a}\mathrm{s}\mathrm{t}}}\right)=\left(\frac{M_p}{M_{\star }}\right)\frac{a_{\mathrm{a}\mathrm{s}\mathrm{t}}}{2a_p{\left(\frac{M_p}{3M_{\star }}\right)}^{\frac{1}{3}}}\\ & \sim {\left[\frac{M_p}{M_{Jup}}\right]}^{\frac{2}{3}}{\left[\frac{M_{\star }}{M_{\odot }}\right]}^{\frac{2}{3}}\left[\frac{a_{\mathrm{a}\mathrm{s}\mathrm{t}}}{10\mathrm{a}\mathrm{u}}\right]{\left[\frac{a_p}{0.1\mathrm{a}\mathrm{u}}\right]}^{-1}.\end{array}$$(5)

The modified ratio of close encounter Λ_c of a planet depends on the semimajor axis a_ast of the interacting asteroid. Since our sampled asteroids cover a range of initial orbits, we adopted a representative midpoint of a_ast = 5 au for this estimation. The resulting values of Λ_c for our simulated (M_p, a_p) pairs are displayed in the right panel of Fig. 6. Theoretically, planets with Λ_c ≳ 1 are able to produce strong direct-scattering effects on asteroids. We highlight three representative cases in the two panels of Fig. 6 to illustrate this correspondence. A filled circle marks the benchmark case discussed previously, where Λ_c > 10. Two triangles indicate cases with Λ_c ~ 0.1 and Λ_c ~ 1, which correspond to protection fractions of ~80% and ~25% in our simulations, respectively. The simulation outcomes strong agree with the theoretical expectations derived from the Safronov-like criterion overall, validating Λc as a useful predictor of planetary shielding efficiency.

Fig. 6 Left panel: protection rates (Fprotection) given by Eq. (3) as a function of mass, with different line styles representing distinct semimajor-axis (ap) values. Right panel: mapping the modified close encounter ratio (Eq. (5)) in the planetary mass-semimajor-axis diagram.

4 Discussion

4.1 Effect of other parameters

We first examined the effect of WD mass by performing simulations using various mass estimates for WD 1856+534 reported in the literature. In these tests, planetary parameters were kept constant, and the integration timescale remained unchanged to ensure a controlled comparison. As shown in Table 3, the resulting asteroid outcome fractions across different WD masses vary by no more than 0.5%. This consistency confirms that the uncertainty in the WD mass within the observationally constrained range does not significantly alter the dynamical evolution or the ultimate protective efficiency of the system.

While the results presented thus far are based on a specific set of initial orbits, it is crucial to verify whether the inferred protective efficiency is sensitive to these assumptions. We therefore examined how varying the asteroid initial semimajor axis (a_ast) and periapsis (q_ast) distributions alters the statistical outcomes. By employing a double-scan strategy (systematically varying a_ast while holding q_ast constant, and vice versa), we isolated the dynamical effect of each parameter. As shown in Figure 7, as the initial q_ast decreases, F_Ejected-OR declines, while F_Ejected-IR increases. This shift indicates that asteroids with smaller initial periapses are more likely to bypass the outer guard of the planet and reach the Roche limit. However, the combined F_Protection remains remarkably stable across the sampled q_ast range, confirming that the choice of initial periapsis has a negligible effect on the overall protective capacity.

The dependence on the initial semimajor axis a_ast can be understood through Eq. (5). Theoretically, there exists a threshold a_in above which the scattering ratio Λ_c > 1, signaling highly efficient ejection. For our benchmark case (13.8 M_Jup at 0.02 au), a_in ~ 0.6, which is well below our sampled range of 1-10 au. This explains why, for massive close-in planets, F_Protection shows little sensitivity to the exact value of a_ast. For lower-mass or more distant planets, however, a_in shifts to higher values. In these regimes, our fixed 1-10 au range serves as a conservative benchmark. While these planets are less effective at shielding the inner system, Eq. (5) implies they would still efficiently eject bodies originating from more distant reservoirs (e.g., a Kuiper-belt-like analog). By maintaining a consistent [1, 10] au range as justified in Section 2.1, we ensured a uniform comparison across the planetary parameter space, while acknowledging that the absolute protection rate for the smallest planets may increase when more distant polluters are considered.

Fig. 7 Relation between the fractions of asteroids within classifications and their initial pericenter.

4.2 Comparison with previous works

Our finding that close-in giant planets serve as a pollution filter is consistent with recent dynamical studies of specific WD systems. For instance, Rogers et al. (2024) investigated the WD 0141-675 system (prior to the retraction of its planetary candidate), finding that for asteroids on high-eccentricity orbits (e = 0.75), a massive planet (9.26 M_Jup) at 0.172 au could eject approximately 52% of them while intercepting 13% via collisions. Similarly, our simulations for a 13.8 M_Jup planet at 0.1 au showed that a total of 77.4% of the asteroids are removed via ejection, with a planetary collision fraction of 3.7%. The differences in these branching ratios stem from two primary methodological factors. First, the higher collision rates reported by Rogers et al. (2024) likely result from their specific setup of injecting particles directly into the planetary chaotic zone. Second, while they adopted the Roche radius as the collision criterion, our model allows asteroids to survive multiple Roche-limit passages, significantly increasing their probability of being eventually ejected by the planet. Despite these technical differences, both studies underscore that close-in massive planets act as efficient barriers, significantly suppressing the pollution rate by either ejecting asteroids or intercepting them before they can reach the WD.

Consequently, our results suggest that this planetary protection might be even more potent than previously estimated for other reservoirs. For example, Pham & Rein (2024) also investigated the effect of a WD1856+534 b like planet on the pollution rate from exo-Oort-cloud comets. They reported that a Jupiter-mass planet at 0.02au reduces the accretion rate by at most 50%, a limit primarily derived from the geometric capacity of the planet to intercept comets via physical collisions. In contrast, our simulations showed that a similar planet can remove roughly 90% of small bodies, indicating that the dynamical scattering barrier is more efficient than simple physical interception. This discrepancy arises because the planetary effect extends well beyond its physical cross section, especially for bodies that survive their initial Roche-limit passage. As summarized in Table 2, for a 1 M_Jup planet at 0.02au, approximately 70% of the total population is eventually scattered out of the system only after penetrating the Roche limit. Our analysis of individual trajectories (e.g., Figs. 2b,c) confirms that these highly eccentric bodies often remain trapped in the gravitational domain of the planet post-passage. Even in cases of tidal disruption, the resulting fragments are more likely to be ejected than accreted. Therefore, we conclude that the protective role of close-in planets is dominated by this sustained dynamical clearing rather than mere physical collision, making the planetary filter significantly more robust than estimates based solely on physical interception.

5 Conclusion

We used N-body simulations to quantify the pollution shielding provided by a planet orbiting a WD. By statistically tracking the outcomes of high-eccentricity asteroids in the inner system, we evaluated how a companion modifies their routes to the WD and thereby reduces contamination. Our main conclusions are listed below.

• Close-in planets can reduce the fraction of asteroids that reach the Roche limit. For a close-in companion with mass 13.8 M_Jup at α_p = 0.02 au, over 10% of the asteroids are barricaded out of the WD Roche radius compared to the noplanet case; the planet alters the asteroid orbital evolution and lowers the incidence of direct tidal disruption.

• Most bodies that enter the Roche region are subsequently removed from the system. Even when the asteroids are scattered into the Roche limit or are tidally disrupted, the resulting fragments still have a high probability of being dynamically ejected rather than accreted onto the WD.

• A WD1856+534b-like planet is capable of reducing pollution to one-tenth of its original level. The scattering shielding is even more effective when the incoming pollutant population is not coplanar with the planet.

• The effectiveness of planetary protection decreases with increasing planet semimajor axis and decreasing planet mass. At a_p ≃ 0.02 au, a planet of only ~1 M_Jup can eject ~90% of incoming asteroids. At larger separations, substantially higher masses are needed to achieve comparable protection.

We conclude that planets and companions can play an important role in regulating WD pollution by both preventing tidal encounters and removing disrupted fragments. The absence of detectable pollution on WD1856+534 is consistent with predictions of planetary protection. With the current limited sample, however, we still cannot statistically distinguish whether this is caused by planetary protection or is merely an observational coincidence. Future systematic discoveries of more similar systems will enable quantitative tests of these theoretical expectations.

Acknowledgements

We would like to thank the referee for their constructive comments, and Siyi Xu for helpful feedback that improved this manuscript. We thank Ke-Ting Shin and Hui-Gen Liu for their valuable discussions on methods, analyses, and figures. This work is supported by the National Key R&D Program of China (Grant No. 2024YFA1611803), the National Natural Science Foundation of China (Grant Nos. 12273011, 12150009 and 12403071), and the Fundamental Research Funds for the Central Universities. J.-W.X. acknowledges support from the National Youth Talent Support Program.

References

All Tables

All Figures

	Fig. 1 Simplified model of the WD system. The planet orbits the central WD on a circular path, while an asteroid initially follows a highly eccentric trajectory toward the WD. The figure is not to scale.
In the text

Fig. 2 Typical asteroid trajectories of four classifications in Section 2.2 on a logarithmic scale, with units in au. The line colors denote classifications, and the diamond markers indicate the starting points of the asteroids. (a) Ejected-OR: asteroid was ejected before it entered the Roche limit. (b) and (c) Ejected-IR: both are ejected after entering the Roche limit, the difference is the duration they spent within the Roche limit. (d) Collide case: asteroid collided with the planet when it crossed the planetary orbit (we depict the planet at the moment of collision (navy circles); note that the planetary radii do not represent the actual planet sizes.). (e) Crash case: asteroid was scattered directly to the WD. (f) Bound: asteroids still orbit the WD until the end of simulation.

In the text

Fig. 3 Distribution of the final asteroid energy and minimum distance to the WD during the entire simulation. The system parameters selected for this simulation are the benchmark parameters of WD1856+534b system given in Table 1. The red dots in the figure represent the initial state of asteroids. The black dots represent the final state of each asteroid group at the end of the simulation; the color of the region in which they land denotes their ultimate fate. The inverted black triangles indicate those that collided with planets.

In the text

	Fig. 4 Proportion of bound asteroids to all simulated asteroids as a function of simulation time. The two different colors represent the trends in different time units.
In the text

	Fig. 5 Evolution of the distance between the particles and the WD over time. The top panel shows the evolution of an asteroid without tidal disruption. The middle and bottom panels show the evolution of the fragments considering the tidal disruption, with different line styles representing different fragments.
In the text

	Fig. 6 Left panel: protection rates (Fprotection) given by Eq. (3) as a function of mass, with different line styles representing distinct semimajor-axis (ap) values. Right panel: mapping the modified close encounter ratio (Eq. (5)) in the planetary mass-semimajor-axis diagram.
In the text

	Fig. 7 Relation between the fractions of asteroids within classifications and their initial pericenter.
In the text