Skip to main content
Sandbox Physics
Debris disks from grain collisions to thermal evidence

AC16 · L2

Couple a fragment mass budget to orbital transport

Sandbox Physics ·

A coupled calculation must let a collision change the particles that later make the image. Begin with one encounter and keep its accounting intact.

A nonrotating frame centered on an accelerating star

Use astronomical units, Julian years and solar masses for dynamics, with the gravitational constant converted from SI. The stellar origin accelerates toward a prescribed Kepler planet. Subtracting that acceleration gives the indirect term in the grain equation. Omitting it would make the trajectory inconsistent with the chosen origin, even when the planet’s position looks plausible.

The planet’s mean motion includes star plus planet mass. Grain gravity includes the radiation-reduced stellar term. PR drag retains physical light speed and removes kinetic energy through both radial and transverse motion. It acts slowly at these radii. The display accelerates the clock while retaining the physical drag coefficient.

r¨=−GM⋆(1−β)rr3+GMp[Rp−r∣Rp−r∣3−RpRp3]−βGM⋆cr2(vrr^+v)\ddot{\mathbf r}=-\frac{GM_\star(1-\beta)\mathbf r}{r^3}+GM_p\left[\frac{\mathbf R_p-\mathbf r}{|\mathbf R_p-\mathbf r|^3}-\frac{\mathbf R_p}{R_p^3}\right]-\frac{\beta GM_\star}{cr^2}\left(v_r\hat{\mathbf r}+\mathbf v\right)
Positions and velocities are relative to the star; the last bracket is the first-order PR force per unit mass.
aPR ⁣⋅ ⁣v=−βGM⋆cr2(vr2+v2)≤0\mathbf a_{\rm PR}\!\cdot\!\mathbf v=-\frac{\beta GM_\star}{cr^2}(v_r^2+v^2)\leq0
The work rate is nonpositive, but planetary work can still change the grain’s total orbital energy.

Reading: REBOUNDx 5.1.0 · radiation forces and physical light speed

Count the same collision once

A parcel represents a large number of physical grains of one radius and velocity. Within a finite spatial cell, two distinct populations produce an encounter rate proportional to the product of their numbers. A single population has unordered pairs; its diagonal rate needs the one-half factor and must remove two grains per event. Processing both ordered pairs would double the mass transfer.

The weighted Monte Carlo implementation considers each distinct parcel pair once. Its packet event removes the smaller available grain number from both participants, with event probability equal to the physical rate times the substep divided by that packet weight. The summed destructive hazard per grain limits the collision substep to 0.05. Parcels with identical velocities have no destructive self-collision; distinct same-size streams can collide.

Finite cells introduce another trap. Subtracting velocities at two different azimuths mistakes smooth Kepler shear for a local collision speed. The operator maps cylindrical peculiar velocities to a common collision position before comparing them. A perfectly cold circular flow then has zero destructive relative velocity. This local interpolation is a numerical approximation, and its momentum and kinetic jumps are recorded separately.

Cij=NiNjVπ(si+sj)2vij,i≠j;Cii=Ni(Ni−1)2Vπ(2si)2viiC_{ij}=\frac{N_iN_j}{V}\pi(s_i+s_j)^2v_{ij},\quad i\ne j;\qquad C_{ii}=\frac{N_i(N_i-1)}{2V}\pi(2s_i)^2v_{ii}
These are physical grain-pair rates. Off-diagonal events consume one grain from each side; diagonal events consume two from the same side.

Reading:

Integrate fragments by mass, not by bin labels

A catastrophic event combines both destroyed grain masses. The teaching threshold is 10 joules per kilogram, applied to the relative center-of-mass impact energy divided by combined mass. Less energetic encounters are omitted; there is no elastic-scattering or cratering prescription hidden behind a destructive event.

For accepted events, use a fragment number law with exponent 3.5 between 0.01 micrometres and half the combined equivalent-volume radius. Mass in a size interval is the integral of radius cubed times the number law. The part below the 0.1 micrometre resolved edge enters an explicit subgrid reservoir. Integrating each clipped bin and the unresolved interval before normalization makes the fractions close exactly.

Representative bins are coarse: an interval can overlap a parent bin even when the continuous upper cutoff lies inside its edge. The representative radius is not an assertion that every daughter has exactly that physical radius. Refining the bins is part of checking emission and transport errors.

Qimp=mimj2(mi+mj)2vij2,Qimp≥QD=10 J kg−1Q_{\rm imp}=\frac{m_im_j}{2(m_i+m_j)^2}v_{ij}^2,\qquad Q_{\rm imp}\geq Q_D=10\,\mathrm{J\,kg^{-1}}
The threshold is a declared synthetic material law, not a measured strength for an identified mineral.
fk=∫max⁡(sk,−,smin⁡)min⁡(sk,+,smax⁡)s3−qf ds∫smin⁡smax⁡s3−qf ds,fsub+∑kfk=1f_k=\frac{\int_{\max(s_{k,-},s_{\min})}^{\min(s_{k,+},s_{\max})}s^{3-q_f}\,ds}{\int_{s_{\min}}^{s_{\max}}s^{3-q_f}\,ds},\qquad f_{\rm sub}+\sum_k f_k=1
An empty overlap contributes zero. For the chosen exponent, the antiderivative is proportional to the square root of radius.

Reading: Stark & Kuchner 2009 · collisional grooming, fixed arXiv 0909.2227v1

Worked encounter: follow both mass and kinetic energy

Take two equal grains with equal and opposite velocities of 100 metres per second. The relative speed is 200 metres per second, and the impact energy is 5000 joules per kilogram of combined mass. After supplying the 10 joule-per-kilogram fracture energy, ten percent of the remaining relative energy becomes symmetric daughter kicks; the rest is assigned to internal energy. The center-of-mass velocity stays zero.

Each size contribution is launched as two opposite kicks. Resolved fragments plus unresolved material carry exactly the original mass and center-of-mass momentum, and kinetic plus fracture plus internal energy closes for this local event. If the center of mass were moving, its translational kinetic energy would remain. Globally, daughters change beta, the planet does work and PR dissipates. Spatial mapping adds separately measured numerical changes, so total mechanical energy is not conserved.

Reading:

A resonance needs an angle, not only a radius

For the exterior 2:1 commensurability, the grain completes one orbit while the planet completes about two. Equate mean motions using the grain’s radiation-reduced central mass and the planet’s total-mass clock. This predicts a nominal semimajor axis, not a guarantee of capture. At 13 astronomical units and the largest representative grain, it lies near 20.6 astronomical units.

The critical angle combines grain mean longitude, planet mean longitude and grain periapsis. If it repeatedly circulates, the encounter geometry sweeps through every phase. Bounded libration instead indicates organized phase evolution over the observed interval. The separate conservative probe starts at 20.9209 astronomical units, eccentricity 0.15 and critical angle pi; its range is 4.69188 radians over 150 planet periods. Removing the planet from the same initial state yields 20.8150 radians of unwrapped circulation. These are finite-window results for that probe, not the full collisional cloud or a permanent-stability proof.

nnp=12,a2:1=ap[M⋆(1−β)M⋆+Mp]1/322/3\frac{n}{n_p}=\frac12,\qquad a_{2:1}=a_p\left[\frac{M_\star(1-\beta)}{M_\star+M_p}\right]^{1/3}2^{2/3}
The nominal location shifts with grain radiation pressure.
ϕ=2λ−λp−ϖ\phi=2\lambda-\lambda_p-\varpi
Compute osculating elliptic elements with the reduced central mass; do not assign elliptic longitudes to unbound grains.

Reading: REBOUNDx 5.1.0 · radiation forces and physical light speed · Han et al. · ARKS II, arXiv 2601.13670v1

What is conserved, and under which switch settings?

Sum active mass with the outer-domain, inner-boundary, subgrid and numerical-stop reservoirs. It must equal initial mass plus any explicitly declared injection. Production records have zero injection. A numerical close-planet stop is not counted as physical accretion, and an outer crossing is not automatically labeled escape.

For a circular prescribed planet, fixed beta, no PR and no collisions, the rotating Hamiltonian is constant. Without a planet, the pressure-reduced two-body energy is constant. Neither test applies unchanged to an eccentric driver or a user dragging its orbit. Independent IAS 15 trajectories test the actual force beyond these limiting invariants.

M0+Minj=Mactive+Mouter+Minner+Msub+MnumericalM_0+M_{\rm inj}=M_{\rm active}+M_{\rm outer}+M_{\rm inner}+M_{\rm sub}+M_{\rm numerical}
The mass error includes the visible grains and every mass reservoir.
J=12v2−GM⋆(1−β)r−GMp∣r−Rp∣+GMpr ⁣⋅ ⁣RpRp3−np(r×v)z\mathcal J=\frac12v^2-\frac{GM_\star(1-\beta)}r-\frac{GM_p}{|\mathbf r-\mathbf R_p|}+\frac{GM_p\mathbf r\!\cdot\!\mathbf R_p}{R_p^3}-n_p(\mathbf r\times\mathbf v)_z
This form is valid for the circular, time-stationary rotating driver and conservative force switches stated above.

Reading: REBOUNDx 5.1.0 · radiation forces and physical light speed

Resampling is a numerical operation, not a physical collision

A finite parcel count limits computational work. The implementation samples within each size and corrects the velocity mean and dispersion so that mass, linear momentum and kinetic energy match before compression. It does not preserve gravitational potential, angular momentum, spatial structure or higher velocity moments. Raw and corrected jumps have different meanings and are saved separately.

The all-parent counterexample shows why this matters. Breaking a fixed mass into smaller grains can greatly increase area and collision hazard. A stiff operator or excessive optical depth is a reason to stop and report a partial record, not to slow the collisions by an arbitrary visual factor. The saved refinement compares density moments and radiative flux after changes in timestep, parcel number, cell grid and seed. A resonant-clump claim also needs stable density structure under those refinements.

Reading:

Saved calculation, readable without the Lab

Solar star, planet mass 0.003 solar masses, circular 13 AU orbit, initial normal optical depth 0.0003. Eight grain bins, 384 initial weighted parcels, RK4 at 128 steps per planet orbit, finite-cell catastrophic collisions and a 128-per-bin compression limit. The representative record is eight planet periods old.

QuantitySaved resultInterpretation
Elapsed time374.42 yr374.42\,\mathrm{yr}A real finite integration window
Weighted collision events590Statistical packet events, not individual grain count
Maximum normal optical depth0.0019440.001944Finite-cell column estimate
Largest raw kinetic compression jump3.16%3.16\%Before explicit velocity-moment correction
Conservative probe angle range4.6919 rad4.6919\,\mathrm{rad}Separate 150-period record, not the whole collisional cloud

The stochastic four-period refinement changes 850 µm flux by about 0.7–4.0% and 70 µm flux by 8.7–32%. Centroid and low-order density moments change enough that individual cloud clumps are not a converged planet inference. The independent conservative probe only establishes bounded angle behavior for its specified initial state and saved 150-period window.

Inspect actual refinement records · Inspect the all-parent finite stress counterexample

Audit an evolving collision record

  1. Load the denser release, let it evolve, and export its physical populations and mass reservoirs.
  2. Compare the conservative frozen angle reference with the evolving cloud; keep the two records distinct.

Exercises with explanations

  1. For one population of six grains, how many unordered pairs exist, and how many grains remain after two diagonal collision events?

    Read the explanation

    There are fifteen possible unordered pairs. Each diagonal event consumes two grains, so two events leave two grains. Using thirty ordered pairs would double count encounters.

  2. Why can a resampling step conserve kinetic energy exactly and still alter the image?

    Read the explanation

    The selected positions and the correlation between position, size and velocity can change. Temperature depends on distance, and sky pixels depend on projected position. Conserving a few velocity moments is insufficient to preserve a radiative map.

Sources and calculation scope

  1. Stark & Kuchner 2009 · collisional grooming, fixed arXiv 0909.2227v1

    Primary full text: sections II.1–II.4, III and IV; local rates, grid/noise tests and the limits of fragment reassignment. This weighted initial-value operator is not a reproduction of their steady-stream solver.

  2. REBOUNDx 5.1.0 · radiation forces and physical light speed

    Official force documentation and example; independent records use actual REBOUND 5.2.2 IAS15 plus REBOUNDx 5.1.0, with physical-c PR and barycentric coordinates converted to the star frame.

  3. Han et al. · ARKS II, arXiv 2601.13670v1

    Fixed 20 January 2026 full-text target, published A&A 705 A196. Close reading: sections 2 and 4.3–4.5, Figure 8 and Appendix H. No ARKS data reproduction or planet discovery is claimed.

FeedbackSomething broken, confusing, or missing on this page?Tell us
What kind of feedback?

Sent with this page’s address, its Lab ID, the site language and your device type. Please leave out personal details. Feedback is kept for up to 12 months.Privacy policy