COMETARY ORIGIN OF THE ZODIACAL CLOUD AND CARBONACEOUS MICROMETEORITES. IMPLICATIONS FOR HOT DEBRIS DISKS

Our model starts with the source population of objects. This may be either of the following: (1) individual asteroid groups such as the Karin, Veritas, or Beagle asteroid families (NVBS06; Nesvorný et al. 2008), (2) asteroid belt as a whole, (3) active JFCs defined by their debiased orbital distribution obtained from Levison & Duncan (1997, hereafter LD97) and their physical lifetime, (4) JFCs with orbital distribution similar to (3) but dynamically evolved beyond their nominal active lifespan, (5) HTCs, and (6) OCCs. These source populations are described in a more detail below.

We ignore the contribution of interstellar dust particles because thermal emission from these small particles (diameter D < 1 μm) would create diagnostic spectral features in the MIR wavelengths that are not observed (e.g., Reach et al. 2003). More specifically, observations of the MIR spectrum of the zodiacal cloud by Reach et al. (2003) suggest that the bulk of the zodiacal cloud is produced by ∼10–100 μm particles. We also ignore dust produced by disruptive collisions in the Kuiper Belt (hereafter KB dust) because Landgraf et al. (2002) and Moro-Martín & Malhotra (2003) showed that KB dust particles should represent only a minor contribution to the inner zodiacal cloud.

Results for (1) were taken from NVBS06 and Nesvorný et al. (2008) who found that the three main dust bands discovered by Low et al. (1984) originate from the Karin, Veritas, and Beagle asteroid families. According to these results (Figure 1), the three main dust bands represent only ≈9%–15% of the zodiacal cloud emission at the ecliptic and <5% overall. About a dozen other dust bands have been identified (Sykes 1990). Since these dust bands are much fainter than the three main dust bands, Nesvorný et al. results set the upper limit on the contribution of identified asteroid breakups to the zodiacal cloud.

The orbital distribution of the asteroid belt source (2) is modeled by using the observed orbital distribution of asteroids with D > 15 km. We obtained this distribution from the ASTORB catalog (Bowell et al. 1994). This sample should be complete and unbiased (Jedicke & Metcalfe 1998; Gladman et al. 2009). The orbital distribution of asteroids with D < 15 km, which may be a better proxy for the initial distribution of dust produced in main belt collisions (Sykes & Greenberg 1986), is roughly similar to that of large asteroids. Thus, the orbital distribution that we use should be a reasonable assumption for (2).

The orbital distributions obtained by numerical integrations of test particle trajectories from individual comets would suffer from the heavily biased comet catalogs. Moreover, there are hundreds of known comets, each producing dust at a variable and typically unknown rate. In addition, it is also possible that the present zodiacal dust complex contains particles from lost parents as suggested by the orphaned Type-II trails (Sykes 1990) and identification of meteoroid streams with disrupted JFCs (see Jenniskens (2008) for a review). In view of these difficulties, we resorted to the following strategy.

The orbital distribution of JFCs was taken from LD97 who followed the evolution of bodies originating in the Kuiper Belt as they are scattered by planets and evolve in small fractions into the inner solar system. Starting at the time when the comets' perihelion distance, q, first drops below 2.5 AU in the LD97 simulations, we include it as an active JFC in our list of source objects. LD97 showed that the orbital distribution of visible JFCs obtained in this way nicely approximates the observed distribution. Moreover, LD97 argued that the inclination distribution of new JFCs (reaching q < 2.5 AU for the first time) is relatively narrow. The inclination distribution widens at later times as the JFC orbits become more spread by Jupiter encounters.

By comparing the width of the model inclination distribution with that of the observed JFCs, LD97 were able to estimate the fading lifetime of JFCs, t_JFC, defined as the characteristic time between their first and last apparitions. They found that t_JFC = 12,000 yr with a 90% confidence interval 3000 < t_JFC < 30,000 yr. See Figure 3 for an illustration of the steady-state JFC orbit distribution for t_JFC = 12,000 yr. This is our initial distribution of particles produced by active JFCs. We use t_JFC as a free parameter in our model with values extending to t_JFC = 100,000 yr to account for the possibility that an important dust component may be produced by old dormant JFCs as they spontaneously disrupt.

Zoom In Zoom Out Reset image size

Our models for the orbital distributions of HTCs and OCCs are simpler than the one described for JFCs above, because we do not have in hand an appropriate numerical model that we could use with confidence. Fortunately, this is not a major limitation factor in this work because our main results described in Section 4 are not sensitive to the detailed properties of the HTCs and OCCs populations.

For HTCs, we assume that the differential distribution of the perihelion distance, dN(q), is dN(q) ∝ q dq, and set an upper limit of q at 3 AU. HTCs typically become visible/active only if they reach q < 3 AU. Similarly, the cumulative semimajor axis distribution of HTCs, N(<a), is assumed to be N(<a) ∝ a with an upper cut on a at 50 AU. The differential inclination distribution dN(i) is taken from Levison et al. (2006). This distribution can be approximately described by

$$\begin{equation} dN(i) \propto \sin (i) e^{-{1\over 2} \left({i\over \sigma }\right)^2}\, di \end{equation} \tag{ 1 }$$

with σ = 30°.

According to Francis (2005), the long-period comets have $dN(q) \propto (1 + \sqrt{q})\, dq$ for q < 2 AU. For q > 2 AU, Francis' study predicts dN(q) being flat or declining while we would expect the perihelion distribution to increase with q. It probably just shows that the distribution is not well constrained for q > 2 AU. We use dN(q) ∝ 2.41(q/2)^γ dq for q > 2 AU with 0 < γ ⩽ 1. The semimajor axis values of OCCs are set between 10,000 and 50,000 AU, which is known as the Oort spike (Wiegert & Tremaine 1999). We also use a = 1000 AU to check on the dynamical behavior of dust particles launched from the returning OCCs. The inclination distribution of OCCs is set to be uniformly random between 0° and 180°.

We tested two different methods to launch particles from their source objects. In the first method, chosen to approximate the ejection of particles from active comets, we launched particles at the perihelion when the mean anomaly M = 0. In the second method, we launched particles along the orbit with uniform distribution in M. This second method is more appropriate for the asteroidal particles and for particles produced by comet disruptions. Indeed, the identified comet disruptions do not seem to be correlated in any way with the perihelion passage (e.g., Weissman 1980). The two ejection methods produced comparable results. Given that disruptions should represent the main mass loss in comets (see Section 6), we use the second method for all results presented in Section 4. To simplify things, we neglect the ejection velocities of dust particles from their parent object and assume that they initially follow the parent object's orbit modified by the radiation pressure.

The individual comets in our model are assumed to contribute in roughly the same proportion to the circumsolar dust complex. The reasoning behind this assumption is that if an individual super-comet were the dominant source of circumsolar dust, the zodiacal cloud would not have such a smooth and symmetrical structure. In fact, the observed smooth structure of the zodiacal cloud probably implies a source population that contains numerous objects that are well mixed in the orbital space.

The orbits of particles were tracked using the Wisdom–Holman map (Wisdom & Holman 1991) modified to include effects of radiation forces (Burns et al. 1979; Bertotti et al. 2003). The acceleration $\vec{F}$ on a particle due to these forces is

$$\begin{equation} \vec{F} = \beta G {m_\odot \over R^2} \left[ \left(1 - {\dot{R} \over c} \right) {\vec{R} \over R} - {\vec{V} \over c} \right], \, \end{equation} \tag{ 2 }$$

where $\vec{R}$ is the orbital radius vector of the particle, $\vec{V}$ is its velocity, G is the gravitational constant, m_☉ is the mass of the Sun, c is the speed of light, and $\dot{R} = dR/dt$. The acceleration (2) consists of the radiation pressure and velocity-dependent PR terms. Parameter β is related to the radiation pressure coefficient, Q_pr, by

$$\begin{equation} \beta = 5.7\times 10^{-5}\ {Q_{\rm pr} \over \rho s},\, \end{equation} \tag{ 3 }$$

where radius s and density ρ of the particle are in cgs units. Pressure coefficient Q_pr can be determined using the Mie theory (Burns et al. 1979). We used Q_pr = 1 which corresponds to the geometrical optics limit where s is much larger than the incident-light wavelength. We assumed that the solar-wind drag force has the same functional form as the PR term and contributes by 30% to the total drag intensity (Gustafson 1994).

We used particles with diameter D = 2s = 10, 30, 100, 200, 300, 1000 μm and set their bulk density to ρ = 2 g cm⁻³. For comparison, Love et al. (1994) reported ρ ≈ 2 g cm⁻³ for stratospheric-collected IDPs, while McDonnell & Gardner (1998) found mean ρ = 2–2.4 g cm⁻³ from the analysis of data collected by the LDEF and Eureca satellites. On the other hand, density of 1 g cm⁻³ has been often assumed for cometary matter (e.g., Joswiak et al. 2007; Wiegert et al. 2009). Grün et al. (1985) suggested that 20%–40% of particles may have low densities whereas most meteoroids have ρ = 2–3 g cm⁻³. Our results may be easily rescaled to any ρ value and we explicitly discuss the effect of ρ wherever it is appropriate.

The particle orbits were numerically integrated with the swift_rmvs3 code (Levison & Duncan 1994) which is an efficient implementation of the Wisdom–Holman map and which, in addition, can deal with close encounters between particles and planets. The radiation pressure and drag forces were inserted into the Keplerian and kick parts of the integrator, respectively. The change to the Keplerian part was trivially done by substituting m_☉ by m_☉(1 − β), where β is given by Equation (3).

The code tracks the orbital evolution of a particle that revolves around the Sun and is subject to the gravitational perturbations of seven planets (Venus to Neptune) until the particle impacts a planet, is ejected from the solar system, or drifts to within 0.03 AU from the Sun. We removed particles that evolved to R < 0.03 AU because the orbital period for R ≲ 0.03 AU is not properly resolved by our integration time step. In Section 4.2, we also consider the effect of collisional disruption by removing particles from the simulations when they reach their assumed physical lifespan. We followed 1000–5000 particles for each D, source, and parameter value(s) that define that source.

Particles were assumed to be isothermal, rapidly rotating spheres. The absorption was assumed to occur into an effective cross section πs², and emission out of 4πs². The infrared flux density (per wavelength interval dλ) per unit surface area at distance r from a thermally radiating particle with radius s is

$$\begin{equation} F(\lambda) = \epsilon (\lambda,s) B(\lambda,T) {s^2 \over r^2}, \, \end{equation} \tag{ 4 }$$

where is the emissivity and B(λ, T) is the energy flux at (λ, λ + dλ) per surface area from a blackbody at temperature T:

$$\begin{equation} B(\lambda,T) = {2 \pi h c^2 \over \lambda ^5} \left[ e^{hc/\lambda k T} - 1 \right]^{-1}. \end{equation} \tag{ 5 }$$

In this equation, h = 6.6262 × 10⁻³⁴ J s is the Planck constant, c = 2.99792458 × 10⁸ m s⁻¹ is the speed of light, and k = 1.3807 × 10⁻²³ J K⁻¹ is the Boltzmann constant. We used = 1 which should be roughly appropriate for the large particles used in this work. See NVBS06 for a more precise treatment of (λ, s) for dust grains composed of different materials.

To determine the temperature of a particle at distance R from the Sun, we used the temperature variations with R that were proposed by different authors from spectral observations of the zodiacal cloud (e.g., Dumont & Levasseur-Regourd 1988; Renard et al. 1995; Leinert et al. 2002; Reach et al. 2003). For example, Leinert et al. (2002) proposed that T(R) = 280/R^0.36 K near R = 1 AU from ISOPHOT spectra. We used T(R) = T_1AU/R^δ K, where T_1AU is the temperature at 1 AU and δ is a power index. We varied T_1AU and δ to see how our results depend on these parameters. Values of T_1AU ≈ 280 K and δ = 0.5 correspond to the equilibrium temperature of large dark particles. Values δ < 0.5 would be expected, for example, for fluffy particles with small packing factors (e.g., Gustafson et al. 2001). Note that the power law is used here as a simple empirical parameterization of the temperature gradient. The real temperature gradient is likely to be a more complicated function of the heliocentric distance (and particle properties).

To compare our results with IRAS observations described in Section 2, we developed a code that models thermal emission from distributions of orbitally evolving particles and produces infrared fluxes that a space-borne telescope would detect depending on its location, pointing direction, and wavelength. See NVBS06 for a detailed description of the code.

In brief, we define the brightness integral along the line of sight of an infrared telescope (defined by fixed longitude l and latitude b of the pointing direction) as

$$\begin{equation} {\it \int _{a,e,i} da de di \int _{0}^{\infty } dr\ r^2 \int _D dD\, F(D,r) N(D;a,e,i) S(R,L,B),} \end{equation} \tag{ 6 }$$

where r is the distance from the telescope and F(D,r) is the infrared flux (integrated over the wavelength range of the telescope's system) per unit surface area at distance r from a thermally radiating particle with diameter D. S(R, L, B) defines the spatial density of particles as a function of the heliocentric distance, R, longitude, L, and latitude, B (all functions of r as determined by geometry from the location and pointing direction of the telescope). N(D; a, e, i) is the number of particles having effective diameter D and orbits with semimajor axis, a, eccentricity, e, and inclination, i.

We evaluate the integral in Equation (6) by numerical renormalization (see NVBS06). F(D,r) is calculated as described in Section 3.3. N(D; a, e, i) is obtained from numerical simulations in Section 3.2. S(R, L, B) uses theoretical expressions for the spatial distribution of particles with fixed a, e, and i, and randomized orbital longitudes (Kessler 1981; NVBS06).

We assume that the telescope is located at (x_t = r_tcos ϕ_t, y_t = r_tsin ϕ_t, z_t = 0) in the Sun-centered reference frame with r_t = 1 AU. Its viewing direction is defined by a unit vector with components (x_v, y_v, z_v). In Equation (6), the pointing vector can be also conveniently defined by longitude l and latitude b of the pointing direction, where x_v = cos bcos l, y_v = cos bsin l, and z_v = sin b. As described in Section 2, we fix the solar elongation l_☉ = 90° and calculate the thermal flux of various particle populations as a function of b and wavelength. The model brightness profiles at 12, 25, and 60 μm are then compared with the mean IRAS profiles shown in Figure 2.

To check our code, known as Synthetic InfraRed Telescope (SIRT), we programmed a particle version of the algorithm, which should be in many ways similar to the core algorithm in SIMUL (Dermott et al. 1988; see also Dermott et al. 2001). The particle version inputs the orbital elements of particles obtained in the orbital simulations (Section 3.2) and produces their orbit clones that are uniformly spread over 2π in mean anomaly M. Thus, every test particle is assumed to trace a cloud of real particles having the same orbit as the test particle but different angular locations along the orbit. See Vokrouhlický et al. (2008; their Section 2.5) for a technical description of the algorithm. This procedure is based on the assumption that any concentration of particles with a specific M value would be quickly dispersed by the Keplerian shear. We employ this procedure to improve the resolution. Without it, the number of integrated test particles would be too small to obtain a useful result.

To be able to compare the results of the particle algorithm with SIRT, the particle algorithm must also use smooth distributions in perihelion longitude ϖ and nodal longitude Ω. This is achieved by generating additional orbital clones with ϖ and Ω uniformly spread over 2π. Figure 4 shows examples of the results obtained from the particle algorithm and SIRT. The agreement between the two codes is excellent which gives us confidence that both codes work properly. We find that the SIRT algorithm based on the Kessler distribution is much faster than the particle one. For this reason, we use the original SIRT code in this study.

Zoom In Zoom Out Reset image size

Figure 5 shows the 25 μm flux of D = 100 μm particles produced by different source populations. Note that these profiles do not sensitively depend on the particle size (see Figures 6 and 7 for the results for different D). Instead, the main differences between results for different source populations in Figure 5 reflect the initial orbit distribution of particles in each source and their orbit evolution. Therefore, these profiles can help us to identify the source population(s) that can best explain IRAS observations.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The asteroidal particles produce a profile with a very sharp peak centered at the ecliptic. The emission from asteroidal particles near the ecliptic poles is relatively faint. The polar emission comes from the particles that evolved by PR drag from a > 2 AU to R = 1 AU. While most asteroidal particles indeed reach 1 AU, they pass too briefly near b = ±90° to produce important polar fluxes. This is why most radiation is received from b ∼ 0, where the telescope collects the thermal emission of particles over a wide distance range. A broader distribution of orbital inclinations is apparently needed to match IRAS measurements.

The profile produced by HTC particles is much broader than the observed one (Figure 5). In this case, the magnitude of the polar fluxes is ≈1/2 of that near the ecliptic. This result reflects the very broad inclination distribution of HTCs (Section 3.1). A potentially significant contribution of HTCs to the zodiacal cloud is also problematic because the two large HTCs, 109P/Swift-Tuttle and 1P/Halley, tend to librate about mean-motion resonances, causing relatively stable orbits for long periods of time. Thus, dust released by HTCs is expected to be concentrated along certain locations on the sky making it difficult to explain the smooth profile of the zodiacal dust. Note also that Altobelli et al. (2007) have not detected HTC particle impacts in the Cassini dust experiment, indicating that HTC particles are relatively sparse.

The OCC particles, which have a nearly isotropic inclination distribution, produce the MIR flux that is constant in latitude (not shown in Figure 5). Therefore, the ecliptic and polar fluxes from OCC particles are roughly the same and do not match observations. We conclude that a single-source model with either asteroidal, HTC, or OCC particles cannot match the observed profile of the zodiacal cloud.

We are left with JFCs. It is notable that the width and shape of the JFC profile in Figure 5 closely match observations. The only slight difference is apparent for large ecliptic latitudes where the model flux, shown here for D = 100 μm and t_JFC = 12,000 yr, is slightly weaker than the one measured by IRAS. We will discuss this small difference below and show that it could be explained if: (1) slightly smaller JFC particles were used, and/or (2) the zodiacal cloud has a faint isotropic component. We thus believe that the close resemblance of our model JFC profile with IRAS data is a strong indication that JFCs are the dominant source of particles in the zodiacal cloud.

Since asteroids and active JFCs have similar inclination distributions (Hahn et al. 2002), it may seem surprising that JFC particles produce substantially wider MIR flux profiles than asteroidal particles. By analyzing the results of our numerical integrations, we find that the encounters with terrestrial planets and secular resonances are apparently not powerful enough to significantly affect the inclination distribution of drifting asteroidal particles. The inclination distributions of the asteroidal particles and their source main belt asteroids are therefore essentially the same (≈10° mean i). On the other hand, we find that JFC particles are scattered by Jupiter before they are able to orbitally decouple from the planet and drift down to 1 AU. This results in a situation where the inclination distribution of JFC particles is significantly broader (≈20° mean i for R < 3 AU) than that of their source JFCs. This explains Figure 5 and shows limitations of the arguments about the zodiacal cloud origin based on the comparative analysis of sources (e.g., Hahn et al. 2002).

We will now address the question of how the MIR fluxes from the JFC particles depend on D and t_JFC. We define

$$\begin{equation} \eta ^2 = {1 \over \pi } \int _{-\pi /2}^{\pi /2} {\left[{F_{\rm IRAS}(b)-F_{\rm model}(b)}\right]^2 \over \sigma ^2(b)}\ db,\, \end{equation} \tag{ 8 }$$

where F_IRAS is the mean IRAS flux, σ is the standard deviation of F_IRAS determined from the spread of representative IRAS scans for each b (Section 2), and F_model = F_JFC (i.e., α_JFC = 1 and α_AST = α_HTC = α_OCC = 0 in Equation (7)). Note that the integration in Equation (8) is set to avoid the intervals in b with strong galactic emission.

While the definition of η² in Equation (8) is similar to the usual χ² statistic (e.g., Press et al. 1992), we will not assign a rigorous probabilistic meaning to the η² values obtained from Equation (8). This is mainly because it is not clear whether the σ values computed in Section 2 from the IRAS data can adequately represent the measurement errors. Instead, we will use Equation (8) only as an indication of whether a particular model is more reasonable than another one. Models with η² ≲ 1 will be given priority. For a reference, the JFC model in Figure 5 gives η² = 5.1.

We calculated η² as a function of D and t_JFC to determine which values of these parameters fit IRAS observations best. We found that the best fits with η² < 10 occur for 30 ⩽ D ⩽ 100 μm and t_JFC ⩽ 30,000 yr.

Figure 7 illustrates how the shape of the 25 μm profile produced by JFC particles depends on D and t_JFC. The profiles become wider with increasing D and t_JFC values. For t_JFC = 12,000 yr, the best results were obtained with D = 30 μm and D = 100 μm (η² = 3.2 and 5.1, respectively). The profiles for D = 10 μm are too narrow and clearly do not fit data well (η² = 70), while those for D = 1000 μm are slightly too wide (η² = 58). We also found that there are no really good fits with t_JFC > 30,000 yr, because the profiles are too broad near the ecliptic independently of the particle size.

The best single-source fits discussed above have η² > 1 which is not ideal. According to our additional tests this is probably not due to the coarse resolution and studied range of D and t_JFC. Instead, this may point to (1) a subtle problem with our JFC model, and/or (2) the possibility that additional minor sources (such as asteroids, OCCs, or HTCs) should be included in the model. Option (1) is difficult to test unless a better model of the JFC population becomes available. Here, we concentrate on (2) because an ample evidence exists (e.g., for ∼10% near-ecliptic asteroid contribution from asteroid dust band modeling) that these additional sources may be relatively important.

We start by discussing the results obtained by assuming that the zodiacal cloud has two sources. The motivation for considering the two-source model was the following. First, we wanted see whether a combination of two sources could successfully fit the observed profile. Second, we attempted to place upper limits on the contributions of asteroid, OCC, and HTC sources. While it is obvious that models with more than two sources can be tuned to fit the data better, it is not clear whether more than two sources are actually required. Our two-source models were used to test these issues.

In the first test, we used the two-source model with asteroids and OCCs (i.e., α_JFC = α_HTC = 0 in Equation (7)). We found that this particular model produces unsatisfactory fits (η² > 10) to IRAS observations for all particle sizes considered here (10–1000 μm) (Figure 8(a)). The model profile is significantly narrower near the ecliptic, where the asteroid component prevails, and is too wide for |b| ≳ 40°, where the OCC component prevails. This happens mainly due to the fact that the asteroid dust is confined to the ecliptic plane and produces a very narrow profile near the ecliptic (Figure 5). We also find it unlikely that two so distinct populations of objects, such as the main belt asteroids and OCCs, would have comparable dust production rates. Thus, we believe that the two-source model with asteroid and OCC dust can be dismissed.

Zoom In Zoom Out Reset image size

In the second test, we set α_OCC = α_HTC = 0 and considered models of the zodiacal cloud with the JFC and asteroid components. We found that a small contribution of asteroid dust can improve the fits. For example, η² = 0.92 for the D = 30 μm JFC particles with α_JFC = 0.9 and t_JFC = 12,000 yr, and D = 100 μm asteroidal particles with α_AST = 0.1. This represents a significant improvement from η² = 3.2 that we obtained for a single-source model with JFC particles only. Values α_AST ≳ 0.3 are clearly unsatisfactory because η² > 10 for α_AST > 0.3. Also, η² > 3.1 for α_AST > 0.2 which shows that the fit does not improve when we add a ≳20% asteroid contribution. These results suggest that a very large asteroid contribution to the zodiacal cloud can be ruled out. This agrees with the conclusions of NVBS06 who found that α_AST ≲ 0.15 from modeling of the main asteroid dust bands.

Finally, the two-component model with the JFC and isotropic OCC sources can fit the data very well (Figure 8(b)). With D = 100 μm and t_JFC = 12,000 yr, corresponding to one of our best single-source fits with JFCs, α_JFC = 0.97 and α_OCC = 0.03, we find that η² = 0.36, by far the best fit obtained with any two-source model. As Figure 8(b) shows, the fit is excellent. We may thus find an evidence for a small contribution of OCC particles to the zodiacal cloud. A much larger OCC contribution is not supported by the data because the fit gets significantly worse for α_OCC > 0.1. For example, η² > 10 for α_OCC > 0.15 which is clearly unsatisfactory. A large contribution of OCC particles can therefore be rejected.

We propose based on the results described above that the zodiacal cloud has a large JFC component (α_JFC ≳ 0.9), and small asteroid/OCC components (α_AST ≲ 0.2 and α_OCC ≲ 0.1). To verify this conclusion, we considered three-component models with α_HTC = 0 and used α_JFC, α_AST, and α_OCC as free parameters. We found that the best two-source fit with α_JFC = 0.97 and α_OCC = 0.03 cannot be significantly improved by including a small asteroid contribution. Similarly, the fit with α_JFC = 0.9 and α_AST = 0.1 cannot be improved by adding a small OCC contribution. Thus, the asteroid/OCC contributions cannot be constrained independently because their effects on the combined profiles can be compensated by adjusting the D and t_JFC values of the dominant JFC particles.

If we set the parameters of the dominant JFC particles to be D = 100 μm and t_JFC = 12,000 yr, however, the α_AST and α_OCC values can be constrained much better (Figure 9). For example, models with η² < 1 require that α_AST < 0.22 and α_OCC < 0.13. Thus, under reasonable assumptions, the contribution of asteroid particles to the near-ecliptic IRAS fluxes is probably <10%–20%, in agreement with the results obtained in NVBS06 from modeling of the asteroid dust bands. This means that asteroid dust contributes only by <10% to the overall zodiacal dust emission at MIR wavelengths. The thermal emission of OCC particles can constitute as much as ∼10% of the near-ecliptic emission with ≈5% providing the best fits (Figure 9). When integrated over latitude, the overall OCC component in the zodiacal cloud is likely to be ≲10%.⁶

Zoom In Zoom Out Reset image size

For the sake of consistency with the results suggested from modeling of the asteroid dust bands (see Section 3.1; NVBS06), we impose a small asteroid contribution in the JFC/OCC model. Figure 10 shows our preferred fit at different IRAS wavelengths. The η² values of this fit in different wavelengths are the following: 0.29 at 12 μm, 0.35 at 25 μm, and 0.06 at 60 μm. This is very satisfactory. Since our model does not include detailed emissivity properties of dust grains at different wavelengths (Section 3.3), we set the emissivity at 25 μm to be 1 and fit for the emissivities at 12 and 60 μm. We found that the relative emissivities at 12 and 60 μm that match the data best are 0.76 and 0.87, respectively. Such a variability of MIR emissivity values at different wavelengths is expected for D ≈ 100 μm silicate particles with some carbon content (NVBS06). Note also that our preferred D values (D ≈ 100 μm) are within the range of dominant sizes of particles at 1 AU as determined from spacecraft impact experiments (D = 20–200 μm; Grün et al. 1985).

Zoom In Zoom Out Reset image size

The observational evidence for collisional disruption of interplanetary particles is undeniable (see, e.g., Grün et al. 1985), yet it is very difficult to model the full collisional cascade in a computer code as each disrupted dust grain produces numerous fragments. The exponentially increasing number of particle fragments, which in reality exceeds 10²⁵ for D > 30 μm (NVBS06), renders any full N-body integration impossible. To circumvent this problem, the N-body integration of a smaller number of "tracer" particles can be coupled with a Monte Carlo model for collisions as in NVBS06. This method is not ideal. Also, any model for the collisional cascade would suffer from our lack of detailed understanding of particle properties and their fragmentation during impacts.

Here, we opt for a very simple approximation of the effect of disruptive collisions. We assume that the collisional lifetime of particles is t_col(D) and stop the N-body integration of diameter D particles when t = t_col(D). Thus, particles keep the same D for t < t_col(D) and vanish at t = t_col(D). This is a very crude approximation of the real collisional cascade, in which particles can be eroded by small collisions and do not vanish upon disruptive impacts (but produce a range of new particle sizes). Also, t_col(D) should be a function of R while we assume here that it is not. Still, as we show below, our simple model should be able to capture the main effects of particle collisions.

Our choice of t_col(D) is motivated by the published estimates of the collisional lifetime of particles based on satellite impact rates and meteor observations. For example, Grün et al. (1985) argued that the collisional lifetime of D = 1 mm particles at 2.5 AU should be ∼10⁴ yr (see also Jacchia & Whipple 1961). This relatively short lifetime is a consequence of the dominant population of D ≈ 100–300 μm particles in the inner solar system (e.g., Love & Brownlee 1993) that are capable of disrupting mm-sized particles upon impacts.

For comparison, the approximate PR drag timescale of particles to spiral down from 2.5 AU to 1 AU is t_PR = 2500 yr × ρs, which for ρ = 2 g cm⁻³ and s = 500 μm gives t_PR = 2.5 × 10⁶ yr. Thus, the PR drag lifetime of these large particles is significantly longer than t_col, indicating that they must disrupt before they can significantly evolve by PR drag. Using this assumption in the model, we found that the profiles produced by large JFC particles with t_col ∼ 10⁴ yr are much narrower in latitude than the ones we obtained in Section 4.1. This is because large particles die before they can evolve to R ∼ 1 AU, where they could contribute to polar fluxes. The zodiacal cloud cross section therefore cannot be dominated by large JFC grains. The large grains are important to explain radar and optical observations of meteors (see Section 5.3; Taylor & Elford 1998; Jenniskens 2006; Wiegert et al. 2009).

On the other side of the size spectrum, D < 10 μm particles have t_PR ≪ t_col due to the lack of small D < 1 μm impactor particles that are blown out of the solar system by radiation pressure, and because the PR drag timescale is short for small D (see, e.g., Dermott et al. 2001). Thus, the small particles are expected to spiral down by PR drag from their initial orbits to R < 1 AU without being disrupted. Our original results described in Section 4.1 are therefore correct for small particles. We showed in Section 4.1 that D < 30 μm JFC particles do not fit IRAS observations well.

Since t_PR ≪ t_col for small particles and t_PR ≫ t_col for large ones, there must exist an intermediate particle size for which t_PR ∼ t_col. These intermediate-size particles are expected to be very abundant in the zodiacal cloud simply because they have the longest lifetimes. Grün et al. (1985) and others argued that the transition from the PR drag to collision-dominated regimes must happen near 100 μm. This is consistent with the LDEF measurements which imply that the D ≈ 200 μm particles represent the dominant mass fraction at R = 1 AU (Love & Brownlee 1993).

The question is therefore how to model collisional effects for D ∼ 100 μm. This is not a simple problem because the effects of the full collisional cascade, including gradual erosion of particles and their supply from breakups of the large ones, should be particularly important in this transition regime. It would be incorrect, for example, to take the Grün et al. estimates of t_col at their face value and remove particles when t > t_col. In reality, each particle can accumulate PR drift during previous stages of evolution when it is still attached to its (slightly) larger precursor particles.

To test these issues, we assumed a wide range of effective t_col and calculated model JFC profiles for each of these cases. Figure 11 shows that the profiles obtained with t_col ⩽ 5 × 10⁵ yr are significantly narrower than the IRAS profiles, even if we tried to compensate for part of the apparent discrepancy by OCC particles (Figure 11(b)). On the other hand, profiles with t_col ⩾ 6 × 10⁵ yr are almost indistinguishable from the original results that we obtained in Section 4.1 with t_col = ∞. The transition between 5 × 10⁵ yr and 6 × 10⁵ yr occurs near the mean PR drag lifetime of D = 100 μm JFC particles in our model. It clearly makes a large difference whether particles are allowed to drift down to R = 1 AU or not. The main lesson we learn from this exercise is that IRAS observations imply that the zodiacal cloud particles have been significantly affected by PR drag.

Zoom In Zoom Out Reset image size

Additional constraints on the micrometeoroid environment near 1 AU are provided by radar and optical observations of meteors. For example, Hunt et al. (2004) determined the meteor entry speeds from the high-gain ALTAIR radar. For 30 km s⁻¹, the minimum detectable mass is 10⁻⁶ g (corresponding to D = 100 μm for ρ = 2 g cm⁻³), while only mm-sized and larger meteoroids can be detected by the ALTAIR radar for <20 km s⁻¹. The ALTAIR measurements represent a significant improvement in sensitivity relative to that of previous radar programs (e.g., Taylor 1995; Taylor & Elford 1998). For example, Taylor (1995) cited the minimum detectable mass of 10⁻⁴ g at 30 km s⁻¹ for the Harvard meteor radar, corresponding to D ≳ 500 μm particles.

In Section 5.3, we estimate that the mean atmospheric entry speed of D ∼ 100 μm zodiacal cloud particles is ≈14 km s⁻¹, and that >90% impact at <20 km s⁻¹. Thus, the ALTAIR and Harvard radars cannot detect the bulk of small zodiacal cloud particles impacting Earth at low speeds. These measurements are instead sensitive to large meteoroids, which carry relatively little total mass and cross section, have short t_col, and are expected to impact on JFC-like orbits. This explains why radar observation show little evidence for populations of small particles with orbits strongly affected by PR drag. In Figure 12, we compare the atmospheric entry speeds of D = 1 mm JFC particles with t_col = 10⁴ yr with the Harvard radar data. This figure documents the dominant role of large JFC particles in meteor radar observations.

Zoom In Zoom Out Reset image size

The spatial distribution of sporadic meteors shows several concentrations on the sky known as the helion/antihelion, north/south apex, and north/south toroidal sources (e.g., Younger et al. 2009, and the references therein). Wiegert et al. (2009) have developed a dynamical model to explain these observations. Their main result concerns the prominent helion/antihelion sources for which the large JFC particles clearly provide the best match. Our model for large JFC particles is in many ways similar to that of Wiegert et al. (2009). It should therefore be consistent with the observed relative strength of the helion/antihelion sources. The more recent high-gain antenna observations show that smaller meteoroids appear to show a weaker helion/antihelion source of eccentric short-period orbits (Mathews et al. 2001; Hunt et al. 2004; Galligan & Baggaley 2004). This implies that orbits of smaller particles should be more affected by PR drag (in agreement with Section 4.2).

The motion of interplanetary particles can be probed by high-resolution spectral observations of the zodiacal cloud. Reynolds et al. (2004) measured the profile of the scattered solar Mg i λ5184 Fraunhofer line in the zodiacal cloud. The measurements indicate a significant population of dust on eccentric and/or inclined orbits. In particular, the inferred inclination distribution is broad extending up to about 30°–40° with respect to the ecliptic plane. The absence of pronounced asymmetries in the shape of the profiles limits the retrograde population of particles to less than 10% of the prograde population.

These results are in a broad agreement with our model. As we discussed in Section 4.1, small JFC particles are scattered by Jupiter before they are able to orbitally decouple from the planet and drift down to 1 AU. This results in a situation where the inclination distribution of JFC particles is broad and extends beyond 20°. The model eccentricities of JFC particles show a broad range of values with most having e = 0.1–0.5 (see Section 5.3). This is in a nice agreement with the analysis of Ipatov et al. (2008) who found that e ∼ 0.3 best fits the Reynolds et al. data.

According to our preferred model with the dominant contribution of JFC particles to the zodiacal cloud, the inner circumsolar dust complex has the total cross section area of (2.0 ± 0.5) × 10¹¹ km². This is a factor of ∼10 larger than the cross section of asteroidal particles in the main asteroid dust bands (NVBS06). The uncertainty of the total cross section was estimated from the range of values obtained for models with η² < 3. Also about 40% of the total estimated cross section of the zodiacal cloud, or ≈8 × 10¹⁰ km², is in particles that reside inside Jupiter's orbit (i.e., with R < 5 AU).

The estimated values are comparable to the effective emitting area of the zodiacal cloud defined as 1 ZODY in Gaidos (1999) (1 ZODY = 1.7 × 10¹¹ km², assuming blackbody emission at 260 K and a bolometric luminosity of 8 × 10⁻⁸ L_☉, where L_☉ is the Sun's value; Reach et al. 1996; Good et al. 1986). Note, however, that we estimate in Section 5.5 that the bolometric luminosity of the inner zodiacal cloud is ∼2.5 times larger than the one assumed by Gaidos (1999).

The total mass of the zodiacal cloud is a function of the unknown particle density and loosely constrained dominant particle size. With ρ = 2 g cm⁻³ and D = 200 μm, we estimate that the total mass is m_ZODY = 5.2 × 10¹⁹ g, which is roughly equivalent to that of a 37 km diameter body. The zodiacal cloud therefore currently contains relatively little mass. Note that these estimates apply to the inner part of the circumsolar dust complex that is detectable at MIR wavelengths. The outer circumsolar dust complex beyond Jupiter is likely more massive due to the contribution from KB particles (e.g., Landgraf et al. 2002; Moro-Martín & Malhotra 2003). According to Greaves et al. (2004), the KB dust disk may represent up to ∼1.2 × 10²³ g. This is up to ∼3 × 10³ times the mass of the inner zodiacal cloud estimated here. Note that this is an upper bound only and that the real KB dust disk can be much less massive.

Our mass estimate is at least a factor of ∼2 uncertain. For example, if ρ = 1 g cm⁻³ or D = 100 μm, we find that m_ZODY = 2.6 × 10¹⁹ g. These values are a factor of ∼2–4 lower than the mass of the zodiacal cloud suggested by NVBS06 from modeling of the asteroid dust bands. NVBS06 assumed that the radial distribution of zodiacal particles is similar to that of the asteroid dust bands, which is incorrect if JFCs are the dominant source. On the other hand, NVBS06 determined the realistic size distribution of zodiacal particles by tracking the collisional evolution, while we used the single-size distributions here.

We estimate that ≳80% of the total mass at <5 AU should be contained in JFC particles. Since these particles can efficiently decouple from Jupiter by PR drag, a large fraction of the total mass is distributed relatively close to the Sun. (For reference, we find that 53%, 19%, and 3.7% of D = 10, 100, and 1000 μm particles released by JFCs, respectively, can decouple from Jupiter.) Figure 13 shows the mass fraction of JFC particles contained in a sphere of radius R around the Sun. The distribution is steep for R < 5 AU and shallower for R > 5 AU reflecting the orbital properties of our model JFC population. About 30% of JFC particles, or about 1.6 × 10¹⁹ g in total mass (for ρ = 2 g cm⁻³ and D = 200 μm), are located within R < 4 AU. Also, ≈10%, or about 5.2 × 10¹⁸ g, has R < 2 AU.

Zoom In Zoom Out Reset image size

For a comparison, assuming that the asteroidal particles with D = 200 μm and ρ = 2 g cm⁻³ contribute by 15% to the near-ecliptic MIR fluxes, we find that the total masses in asteroidal particles with R < 2 AU and R < 4 AU are 5.3 × 10¹⁷ g and 1.3 × 10¹⁸ g, respectively. Thus, the total mass (or number) ratio of JFC to asteroidal particles in the inner solar system is ≲10. Note that this estimate applies as far as the size distributions of JFC and asteroidal particles in the zodiacal cloud are similar, which is expected because both particle populations live in the common collisional environment and have similar PR drag timescales.

We used the Öpik algorithm (Öpik 1951; Wetherill 1967) to estimate the terrestrial accretion rate of JFC particles expected from our model. For 30 < D < 300 μm, the average impact probability of JFC particles on the Earth is ∼5 × 10⁻⁹ yr⁻¹ per one particle in the zodiacal cloud. A similar value is obtained if the impact probability is estimated from the number of direct impacts recorded by the N-body integrator. Thus, in a steady state with ∼2 × 10¹⁹ g in the zodiacal cloud, we estimate that ∼10⁵ tons of JFC particles should be accreted by the Earth annually. This is larger than the nominal Earth's accretion rate of 20,000–60,000 tons yr⁻¹ as determined from LDEF (Love & Brownlee 1993) and the antarctic MM record (Taylor et al. 1996).

This may imply that the real Earth's accretion rate is somewhat larger than the LDEF values. Alternatively, the LDEF constraints may imply that the real mass of the zodiacal cloud is lower than the one estimated here. As we discussed above, the mass of the zodiacal cloud estimated here from the IRAS data is at least a factor of ∼2 uncertain. It is thus plausible that m_ZODY ∼ 10¹⁹ g (e.g., if ρ = 1 g cm⁻³), which would bring our results into a better agreement with LDEF. Additional uncertainty in these estimates is related to the effects of collisional disruption of particles and continuous size distribution.

For comparison, if we assumed that D = 200 μm asteroidal particles are producing the full near-ecliptic MIR flux measured by IRAS, the estimated terrestrial accretion rate of asteroidal particles would be 1.5 × 10⁵ tons yr⁻¹. According to NVBS06 and the results obtained here, however, the asteroidal particles contribute by only ∼10% of the near-ecliptic MIR flux. Thus, we find that the asteroid particle accretion rate should be ∼15,000 tons yr⁻¹, or only ∼15% of the JFC particle accretion rate. The asteroidal particles should therefore represent a relatively minor fraction of IDPs and MMs in our collections. This explains the paucity of the ordinary chondritic material in the analyzed samples (see, e.g., Genge 2006).

Using the same assumptions, we estimate from our model that 16,000 tons yr⁻¹ and 1600 tons yr⁻¹ of JFC particles should be accreted by Mars and the Moon, respectively. The accretion rate of JFC particles on the Moon is thus only about ∼2% of Earth's accretion rate. This corresponds to a smaller physical cross section and smaller focusing factor of the Moon. The mass influx on Mars is ∼20% of Earth's accretion rate. For a comparison, 1600 tons and 100 tons of asteroidal particles are expected to fall on Mars and the Moon annually.

Love & Brownlee (1993) found from the LDEF impact record that D ≈ 200 μm particles should carry most of the mass of zodiacal particles near 1 AU, while we find here that D ≈ 100 μm provides the best fit to IRAS observations. This slight difference may be related to some of the limitations of our model. It can also be real, however, because the LDEF size distribution computed by Love & Brownlee (1993) is bending from the steep slope at D > 300 μm to the shallow slope at D < 50 μm. The cross section area should therefore be dominated by smaller particles than the mass. From Figure 4 in Love & Brownlee (1993), we estimate that D ≈ 100 μm particles should indeed dominate the total cross section area of the zodiacal cloud at 1 AU.

These results have implications for the origin of MMs. MMs are usually classified according to the extent of atmospheric heating they endure (e.g., Engrand & Maurette 1998). Cosmic spherules are fully melted objects. Scoriaceous MMs are unmelted but thermally metamorphosed objects. The fine-grained MMs and coarse-grained MMs are unmelted objects which can be distinguished on the basis of their grain size. Based on bulk composition, carbon content, and the composition of isolated olivine and pyroxene grains, fine-grained MMs and scoriaceous MMs, which appear to be thermally metamorphosed fine-grained MMs, are likely related to carbonaceous chondrites. It has been estimated that the ratio of carbonaceous to ordinary chondrite MMs is ∼6:1 or larger (see, e.g., Levison et al. 2009). This stands in stark contrast to the terrestrial meteorite collection, which is dominated by ordinary chondrites.

A possible solution to this discrepancy is that a large fraction of the collected MMs are particles from the JFCs. This possibility has to be seriously considered because we find here that the carbonaceous JFC grains should prevail, by a large factor, in the terrestrial accretion rate of micrometeoroids. It has been suggested in the past that a possible problem with this solution is that the cometary particles should encounter the Earth at large velocities (e.g., Flynn 1995), so that they either burn up in the atmosphere or are converted into cosmic spherules. Thus, while cometary particles could produce fully melted objects such as the cosmic spherules it was not clear whether the less thermally processed carbonaceous MMs, such as the fine-grained and scoriaceous MMs, may represent cometary material.

By assuming t_col ≳ 5 × 10⁵ yr for D = 100 μm as required by IRAS observations (see Section 4.2), we find from our model that the mean impact speed of D ∼ 200 μm JFC particles on Earth is ≈14.5 km s⁻¹ (Figure 14; see Section 4.3 for a discussion of the size dependence of impact speed and its relevance to meteor observations). This value is only slightly higher than that of the asteroidal particles (≈12.5 km s⁻¹). The comparable impact speeds of JFC and asteroidal particles in our model are a consequence of PR drag which efficiently circularizes the orbits before they can reach 1 AU (Figure 15). We thus find that the impact speeds of the JFC particles are low and do not pose a serious problem. Based on this result and the high terrestrial accretion rate of JFC particles on Earth (Section 5.2), we propose that the carbonaceous MMs in our collections are grains from the JFCs. A large contribution from primitive material that may have been embedded into the main asteroid belt according to Levison et al. (2009) is probably not needed.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Dobrica et al. (2010) compared antarctic MMs with particles that were collected from JFC 81P/Wild (known also as Wild 2) by the Stardust spacecraft. They found that the bulk composition (CI-like) and petrography (presence of chondrules and calcium–aluminium-rich inclusions; CAIs) is similar. Moreover, MMs and 81P/Wild particles have strikingly similar mineralogy. For example, the pyroxene and olivine have similar abundances in both samples, which is unusual for most meteorite groups which typically have more olivine (Lauretta & McSween 2008). The same three types of sulfides (troilite, pyrrhotite, and pentlandite) were found both in MMs and Wild-2 particles. In addition, trace elements such as chromium and manganese have similar concentrations in the collected samples. Though phyllosilicates, present in some MMs (Engrand & Maurette 1998), were not found in Wild-2 particles, they were detected in JFC Tempel 1 (Lisse et al. 2006). Finally, Wild-2 particles and MMs have the same oxygen isotopic composition of chondrules and CAIs's components (Dobrica et al. 2010).

These results indicate that MMs have physical properties similar to those of cometary particles, which is consistent with most MMs being derived from comets. This provides independent evidence (to that obtained from our dynamical model) for cometary origin of MMs.

It is believed that the main source of JFCs is the scattered trans-Neptunian disk, which should have decayed by a factor of ∼100 over the past ∼4 Gyr (LD97; Dones et al. 2004). If the JFC population decayed proportionally, we can estimate that the ecliptic component of the zodiacal dust should have been ∼100 times brighter initially that it is now. This corresponds to the near-ecliptic 25 μm flux of about 7 × 10³ MJy sr⁻¹.

A different insight into the historical brightness of the zodiacal cloud can be obtained in the framework of the Nice model (Tsiganis et al. 2005), which is the most complete model currently available for the early evolution of the outer solar system. In the Nice model, the giant planets are assumed to have formed in a compact configuration (all were located at 5–18 AU). Slow migration was induced in these planets by gravitational interaction with planetesimals leaking out of a massive primordial trans-planetary disk. After a long period of time, most likely some 700 Myr after formation of the giant planets (Gomes et al. 2005), planets crossed a major mean-motion resonance. This event triggered a global instability that led to a violent reorganization of the outer solar system. Uranus and Neptune penetrated the trans-planetary disk, scattering its inhabitants throughout the solar system. Finally, the interaction between the ice giants and the planetesimals damped the orbits of these planets, leading them to evolve onto nearly circular orbits at their current locations.

The Nice model is compelling because it can explain many of the characteristics of the outer solar system (Tsiganis et al. 2005; Morbidelli et al. 2005; Nesvorný et al. 2007; Levison et al. 2008; Nesvorný & Vokrouhlický 2009). In addition, the Nice model can also provide an explanation for the Late Heavy Bombardment (LHB) of the Moon (Tera et al. 1974; Chapman et al. 2007) because the scattered inhabitants of the planetesimal disk, and main belt asteroids destabilized by planetary migration, would provide prodigious numbers of impactors in the inner solar system (Levison et al. 2001; Gomes et al. 2005).

Assuming that the historical brightness of the zodiacal cloud was proportional to the number of primitive objects that were scattered into the inner solar system on JFC-like orbits, we can estimate how it changed over time. In the pre-LHB stage in the Nice model, the leakage rate from the planetesimal disk beyond 15 AU was likely not significant relative to that at LHB. We thus expect that the MIR emission from the inner zodiacal cloud at R < 5 AU should have been relatively faint, except if a massive population of particles was sustained by collisions in the pre-LHB asteroid belt. Here, we focus on the LHB and post-LHB stages.

According to Wyatt et al. (2007), the asteroidal debris disk is expected to decay by orders of magnitude from the time of Jupiter's formation, which marked the start of the fragmentation-dominated regime in the asteroid belt (e.g., Bottke et al. 2005), to LHB. It thus seems unlikely that a massive population of debris could be sustained over 700 Myr by the collisional grinding of main belt asteroids. Instead, it has been suggested that the collisional grinding in the massive trans-planetary disk at R > 15 AU should have produced strong MIR emission peaking at ∼100 μm (Booth et al. 2009, hereafter B09). Being more distant the trans-planetary disk probably decayed more slowly by collisions than the asteroid belt. Thus, in the pre-LHB stage, the Wien side of the trans-planetary disk emission may have exceeded the one from the inner zodiacal cloud down to ∼20 μm (B09).

During the LHB, as defined by the Nice model, large numbers of outer disk planetesimals were scattered into the inner solar system and the inner zodiacal cloud could have become orders of magnitude brighter than it is now. To estimate how bright it actually was, we used simulations of the Nice model from Nesvorný & Vokrouhlický (2009, hereafter NV09). NV09 numerically tracked the orbital evolution of four outer planets and 27,029 objects representing the outer planetesimal disk. The mass of the disk was set to be 35 Earth masses. In total, NV09 performed 90 different numerical integrations of the Nice model, only some of which ended with the correct orbits of the outer planets.

We used these successful simulations to determine the number of scattered objects with JFC-like orbits as a fraction of the total initial number of planetesimals in the trans-planetary disk. Figure 16 shows how this fraction changed over time in one of the NV09 successful simulations. Consistently with the estimated physical lifetime of modern JFCs (LD97), we assumed that the physical lifetime of planetesimals after reaching q < 2.5 AU for the first time was 10⁴ yr. Objects past their physical lifetime did not contribute to the statistic.

Zoom In Zoom Out Reset image size

Immediately after the planetary instability occurred in the Nice model, the estimated fraction of planetesimals having JFC-like orbits was ≈7 × 10⁻⁵ (Figure 16). It then decayed to ≈10⁻⁶ at 50 Myr after the start of LHB. Even though the NV09 simulations gradually lose resolution at later times due to the insufficient number of tracked particles, we can still estimate that the fraction was ∼10⁻⁸ at 500 Myr, or about 3.4 Gyr ago in absolute chronology.

Charnoz et al. (2009) and Morbidelli et al. (2009) argued, using the crater record on Iapetus and the current size distribution of Jupiter's Trojans, that the total number of D > 2 km planetesimals in the pre-LHB trans-planetary disk was ∼10¹⁰–10¹². Using this value and Figure 16, we find that there were ∼7 × 10⁶ JFCs with D > 2 km at time of the LHB peak, t_LHB, and ∼10⁵ JFCs at t_LHB + 50 Myr. These estimates are at least an order of magnitude uncertain mainly due to the poorly known size distribution of small planetesimals in the trans-planetary disk.

For comparison, Di Sisto et al. (2009) found, in a good agreement with the previous estimates of LD97, that there are ≈100 JFCs with D > 2 km and q < 2.5 AU in the current solar system (with about a factor of 50% uncertainty in this value). Therefore, if the inner zodiacal cloud brightness reflects variations in the size of the historical JFC population, we find that it has been ∼7 × 10⁴ brighter at t_LHB and 2 × 10³ brighter at t_LHB + 50 Myr than it is now. This would correspond to the near-ecliptic 25 μm fluxes of 5 × 10⁶ and 10⁵ MJy sr⁻¹, respectively. These values largely exceed those expected from dust particles that were scattered from the trans-planetary disk (B09). Most of the action was apparently over by t_LHB + 500 Myr, when our model suggests that the inner zodiacal cloud was only ∼10 times brighter than it is now.⁷

Figure 17 shows how the present zodiacal cloud would look to a distant observer. If seen from the side, the brightest inner part of the zodiacal cloud has a disk-like shape with a ≈1.6 ratio between the ecliptic and polar dimensions. Similar shapes have been reported by Hahn et al. (2002) from Clementine observations of scattered light. At a larger distance from the Sun, the shape of the zodiacal cloud is oblate and shows cusps at the ecliptic. The axial ratio becomes ≈1.3 at R = 5 AU.

Zoom In Zoom Out Reset image size

The radial brightness profiles in Figure 17 show a steep dimming of the zodiacal cloud with R. For R < 1 AU, a factor of ∼10 in brightness is lost per 1 AU. For 1 < R < 5 AU, factor ∼10 is lost per 2 AU. These profiles are approximate because we ignored the effect of collisions in our model, which should be especially important for R ≲ 1 AU. It is unclear how the shape of the zodiacal cloud would look for R > 5 AU because we did not model the contribution from KB dust.

Figure 18 shows the spectral energy distribution (SED) for distant unresolved observations of the zodiacal cloud. At a distance of 10 pc from the Sun, SED of the present inner zodiacal cloud is 1.4 × 10⁻⁴ Jy at 24 μm and 5.5 × 10⁻⁵ Jy at 70 μm, corresponding to the excesses over the Sun's photospheric emission at these wavelengths of about 3.4 × 10⁻⁴ and 1.1 × 10⁻³, respectively. For comparison, the approximate 3σ excess detection limits of Spitzer telescope observations of Sun-like stars are 0.054 at 24 μm and 0.55 at 70 μm (Carpenter et al. 2009). The MIR emission of the present inner zodiacal cloud is therefore undetectable by distant unresolved observations with a Spitzer-class telescope. Specifically, the detectable emission levels are ≈160 and ≈500 larger at 24 and 70 μm, respectively, than those of the present inner zodiacal cloud.

Zoom In Zoom Out Reset image size

When the flux is integrated over wavelengths, we find that the fractional bolometric luminosity of the inner zodiacal cloud, L_ZODY, relative to that of the Sun, L_☉ = 3.839 × 10²⁶ W, is L_ZODY/L_☉ ∼ 2 × 10⁻⁷. This is a larger value than 10⁻⁸–10⁻⁷ suggested by Dermott et al. (2002) and perhaps comparable to that of KB dust at >30 AU (Stern 1996). The effective blackbody temperature of the zodiacal cloud can be estimated from T_eff = 5100/λ_max, where λ_max is the wavelength of the SED maximum in microns. With λ_max = 18.5 μm, this gives this gives T_eff = 276 K.

B09 studied how the MIR excess of the solar system debris disk varied with time. According to them, the main source of the pre-LHB MIR emission should have been the population of dust particles produced by collisions in the massive trans-planetary disk at R > 20 AU. In Figure 18, we show the model SED produced by the B09 trans-planetary disk. Being dominated by collisions (as opposed to PR drag regime; see Wyatt 2005), the trans-planetary particles are destroyed before they could evolve to R < 20 AU. The SED emission therefore peaks at longer wavelengths (≈100 μm) than the SED of the present zodiacal cloud (≈20 μm). Also, with ∼35 Earth masses in the pre-LHB trans-planetary disk, its estimated MIR emission is strong and produces excesses of ∼0.1 at 24 μm and ∼50 at 70 μm over Sun's photospheric emission at these wavelengths. These values are comparable to those of observed exozodiacal debris disks (B09).

The trans-planetary disk objects, including small dust particles, became scattered all around the solar system during LHB. This led to a significant depletion of the trans-planetary particle population which could not have been compensated by the collisional cascade because collisions became increasingly rare in the depleted disk. The MIR excess should have thus dropped by orders of magnitude within several hundred Myr after the LHB start. B09 estimated that the 24 μm excess of dust particles scattered from the trans-planetary disk should have dropped to ∼3 × 10⁻⁵ at the present time. This is about an order-of-magnitude lower value than the 24 μm excess estimated by us for the current inner zodiacal cloud. Thus, there must have been a transition epoch some time after LHB when the 24 μm excess stopped being dominated by dust particles scattered from the trans-planetary disk and became dominated by particles produced by JFCs.

Figure 19 illustrates how the MIR emission of the inner zodiacal cloud should have varied with time during LHB. The size of the JFC population was estimated by using the methods described in Section 5.4. We then scaled up the MIR emission of the present inner zodiacal cloud by the appropriate factor (see Section 5.4). We found that the 24 μm excess reached F_ZODY/F_Sun(24 μm) ∼ 20 at the LHB peak and stayed for about 100 Myr above the Spitzer's 3σ detection limit. It dropped down to ∼10 times the value of the present zodiacal cloud at t_LHB + 500 Myr. We were unable to determine how this trend continues after t_LHB + 500 Myr because of the resolution issues with the NV09 simulations (Section 5.4). We expect that F_ZODY/F_Sun(24 μm) should have decayed by an additional factor of ∼10 from t₀ + 500 Myr to the present time.

Zoom In Zoom Out Reset image size

The 70 μm excess behaves similarly (Figure 19(b)). It reaches F_ZODY/F_Sun(70 μm) ∼ 70 at the LHB peak and decays at later times. Given the tighter detection limit of Spitzer at 70 μm, the 70 μm excess would remain detectable by Spitzer for 50 Myr, which is roughly half of the interval during which the 24 μm excess could be detected. Also, when these values are compared to the ones estimated by B09 for the pre-LHB trans-planetary disk, we find that the 70 μm excess expected for JFC particles at the LHB peak is comparable to that of the pre-LHB excess produced by collisions in the trans-planetary disk. We would thus expect that the solar system's LHB has not produced any significant increase of the 70 μm emission.

Conversely, the 24 μm excess raises a factor of ∼100 above the B09 pre-LHB level, indicating that the solar system became significantly brighter during LHB at these shorter wavelengths. During LHB, the emission from JFC particles should have exceeded that of the scattered trans-planetary particles by ∼20 (compare Figure 19(a) with Figure 5 in B09). Thus, the system does not need to have a significant cold dust disk at the same time as the hot dust disk to provide material for the hot disk. In fact, we find here that the hot disk can be fed by D ≳ 1 km objects, which have little total cross section area to be detected in the cold disk, but large mass to sustain the hot disk upon their disintegration at <10 AU (see Section 6). Since the decay rates of both populations after LHB should have been similar, their ratio should have remained roughly constant over time suggesting that the trans-planetary dust did not represent any significant contribution to the post-LHB emission of the zodiacal cloud at 24 μm.

These results have interesting implications for our understanding of hot debris disks observed within 10 AU around mature stars (Trilling et al. 2008; Meyer et al. 2008; Carpenter et al. 2009; see also a review by Wyatt 2008). It has been argued that some of the observed brightest hot disks, such as HD 69830, HD 72905, HD 23514, η Corvi, and BD+20307, cannot be explained by assuming that they are produced by the collisional grinding of the local population of asteroids (Wyatt et al. 2007). Specifically, Wyatt et al. (2007) pointed out that the emission from a locally produced population of debris is expected to be much weaker than the observed emission because disks become depleted over time by collisions. This problem cannot be resolved by assuming a more massive initial population because the massive population would decay faster. Instead, Wyatt et al. proposed that the bright hot debris disks can be seen around stars with planetary systems that are undergoing the LHB instability akin to that invoked in the Nice model.

In Section 5.6, we estimated how the 24 μm and 70 μm excesses varied during the solar system's LHB. We found that the 24 μm excess should have rapidly risen by a large factor from the pre-LHB value and then gradually decayed. It would remain detectable by a Spitzer-class telescope for about 100 Myr after the LHB start, or ∼2% of Sun's current age. The 24 μm excess reached values ≳10 at LHB, which is comparable to those of the brightest known hot disks (Wyatt et al. 2007). Conversely, the solar system's LHB has not produced a significant increase of the 70 μm excess relative to the pre-LHB level. The 70 μm excess decayed after LHB and became undetectable by a Spitzer-class telescope after ∼50 Myr. Thus, if the timing is right, debris disks may show the 24 μm excess but not the 70 μm excess. This could be relevant for systems such as HD 69830, which shows a large excess emission at 8–35 μm (Beichman et al. 2005), but lacks 70 μm emission, and HD 101259 (Trilling et al. 2008).

In a broader context, our study of the zodiacal cloud implies that (1) the populations of small debris particles can be generated by processes that do not involve, at least initially, disruptive collisions (see Section 6) and (2) observed hot dust around mature stars may not be produced from a population of objects that is native to <10 AU. Instead, in the solar system, most particles located within the orbit of Jupiter are fragments of planetesimals that formed at >15 AU. These icy objects are transported to <5 AU by gravitational encounters with the outer planets and disintegrate into small particles by disruptive splitting events (thought to occur due to processes such as the pressure buildup from heated volatiles or nucleus spin-up; see Section 6). If these processes are common around stars harboring planets, the collisional paradigm in which debris disks are explained by collisions of local populations of planetesimals (e.g., "exo-asteroids") may not be as universal as thought before.

If our basic assumptions are correct, large quantities of dust should have been accreted by the Moon, Earth, and other terrestrial planets during LHB. For example, assuming 10¹⁴ g yr⁻¹ mean accretion rate over 100 Myr we estimate that ∼10²² g of extraterrestrial material should have fallen on the Earth at the time of LHB (with ∼50% of this mass accumulating in the first 10 Myr). This is ∼50 times more mass than the quantity accumulated by the Earth at its current accretion rate over 4 Gyr. The Moon should have accreted about 2% of Earth's value during LHB, or ∼2 × 10²⁰ g in total over 100 Myr.

These estimates are at least 1 order of magnitude uncertain. They were obtained by scaling up the present accretion rates by the ratio between the estimated population of JFCs at LHB and their present number. The contribution from planetesimals larger than the largest presently observed JFCs was ignored. This contribution may have been important if most mass of the trans-planetary disk were in D ≳ 20 km planetesimals (Morbidelli et al. 2009). Assuming that all planetesimals reaching q < 2.5 AU during LHB became disrupted and dissolved into small particles, we estimate that up to ∼10²⁴ g of this material could have been accreted by the Earth.⁸

In reality, significant mass losses must have occurred because many planetesimals or their large fragments were removed by Jupiter before they could fully disintegrate. Also, small dust grains could have been disrupted and blown away by the radiation pressure before they would reach 1 AU. The real mass accreted by the Earth during LHB is thus probably significantly lower than ∼10²⁴ g. This also places an upper mass limit, ∼2 × 10²² g, on the LHB accretion of dust particles by the Moon. For comparison, the mass of large impactors estimated from the number and size distribution of lunar basins is 6 × 10²¹ g (Hartmann et al. 2000). Thus, the total mass of carbonaceous dust deposited on the Moon during LHB could have been as low as ∼1/30 of that of the large impactors, or up to 3 times larger. These larger values may be inconsistent with the lunar geological record.

The LHB is of fundamental interest in studies of the origin of life because it immediately precedes the oldest evidence for a biosphere (Awramik et al. 1983; Schidlowski 1988; Mojzsis et al. 1996). The significance of our results in this context is that JFC dust grains can bring in unaltered primitive material from the outer solar system. They could potentially be the source of the earliest organic material that gave rise to life on Earth (e.g., Jenniskens 2001; Jenniskens et al. 2004). Note that the high abundance of organics is a defining property of comets (Bockelée-Morvan et al. 2004). Both MMs and Wild-2 particles contain two important classes of organic molecules: amino acids and PAHs (Sandford et al. 2006; Engrand & Maurette 1998).