Experiments and the theory behind them
Each experiment has an ID, the question it answers, a config under experiments/En_<slug>/config.toml, a script that writes outputs with a run manifest, and one or more figures. Sizes given are the study sizes; CI runs a reduced version (see 03-verification.md).
1. Estimator theory (basis for E5 to E13)
Fix the state N (one evaluation of the RHS). The target is the vector
F_k = Σ_(i,j) K_kij N_i N_j = Σ_(i,j) a_ij v^(ij)_k, with a_ij = N_i N_j R_ij and v^(ij)_k = X_kij,
where v^(ij) is the per-collision recipe of pair (i, j) (its gain/loss vector, ADR 0005). Write T_ij = a_ij v^(ij) for the term vector.
Thesis estimator (without replacement, no weights). Drawing a set 𝒮 of distinct pairs with inclusion probabilities π_ij and summing T_ij over 𝒮 gives E[F̂] = Σ π_ij T_ij. With π_ij < 1 every process is slowed down by its own factor π_ij, so the bias is B = −Σ (1 − π_ij) T_ij. This is finding F1.
Horvitz-Thompson (HT). F̂ = Σ_{(i,j)∈𝒮} T_ij / π_ij is unbiased for any design with π_ij > 0 wherever T_ij ≠ 0.
Importance sampling with replacement (IS). Draw S pairs i.i.d. from p, F̂ = (1/S) Σ_s T_{(ij)s} / p{(ij)_s}. Unbiased; covariance (1/S)(Σ T T^⊤/p − F F^⊤). This is Dullemond’s Eq. 1.70 with p_ij = F_i F_j.
Variance-optimal proposal. For a weighted norm ‖x‖²_w = Σ_k w_k x_k², E‖F̂ − F‖²_w = (1/S)(Σ_ij ‖T_ij‖²_w / p_ij − ‖F‖²_w). Minimizing over p with Σp = 1 (Lagrange, Cauchy-Schwarz) gives
p*_ij ∝ ‖T_ij‖_w = N_i N_j R_ij ‖v^(ij)‖_w.
With the mass weighting w_k = m_k² (error in mass per bin) and replacing L2 by L1 (a common robust choice) this is p_ij ∝ N_i N_j Σ_k m_k |K_kij| = N_i N_j ΔM^net_ij: the relevance measure the thesis wrote down, without the extra m_i m_j factor (F3). Cost per draw differs between pairs (recipe length L_ij is 3 to 4 for coagulation, up to N_m for shattering), so the cost-aware optimum is p_ij ∝ ‖T_ij‖ / sqrt(L_ij). E7 tests this.
Hybrid estimator. Sort pairs by p*_ij; sum the top set 𝒟 holding a fraction q of the proposal mass exactly, sample the rest. Unbiased, variance from the tail only. For peaked distributions (the typical case, notes section 1.3.1) a small 𝒟 removes most variance. E12.
Implicit sampled form. Replacing a_ij = N_i N_j R_ij by a_ij evaluated at the unknown end-of-step state, but with weights from the start-of-step proposal (Dullemond Eq. 1.71), gives a sparse quadratic kernel K̃ with one entry block per drawn pair. Unbiased at the start of the step; first-order bias in Δt that vanishes as the step shrinks. E10 measures it.
Detailed balance per sample. Because each draw contributes a full recipe vector v^(ij) whose mass moment is zero, every sampled RHS conserves mass exactly, independently of which pairs were drawn. Mass conservation does not depend on estimator quality. ADR 0006.
2. Experiment catalogue
Physics and baseline
| ID | Question | Setup | Output |
|---|---|---|---|
| E1 | Do our disk profiles match the thesis and diskmodel.py? What does the viscosity fix (F6) change? | Thesis disk, r ∈ [10⁻², 10²] au; ν_mol vs ν_t | Profiles 2×2; Δv_RD with both choices |
| E2 | Which Δv component dominates where, and how do toggles change the kernel? | N_m = 100, m ∈ [10⁻¹⁵, 10⁻⁸] kg; each component on/off | Δv maps per component, C_ij, P_frag, P_bounce maps |
| E3 | Kernel structure and conservation | coag, frag (pulverization, MRN up to m_i+m_j, independent shattering, cratering), bouncing; with/without cancellation handling; f64 vs double-double | Kernel slices (thesis Figs. 101, 102), ΔX^rel maps, sparsity vs N_m |
| E4 | Integrator correctness and cost on deterministic kernels | Analytic kernels; thesis full model; all integrators; f64 vs double-double long-term mass error; DustPy steady state at matched settings | Convergence plots; work-precision; steady state; stiffness spectrum |
| E5 | Faithful reproduction of thesis chapter 4 | Thesis estimator with scalar weights (as run) and with documented weights | Thesis Figs. 104, 106, 107, 109 regenerated |
Estimators on a fixed state (no time integration)
| ID | Question | Setup | Output |
|---|---|---|---|
| E6 | RQ1: bias vs variance of the thesis estimator | Snapshots N(t) from the full-kernel run at 6 times; R = 10⁴ replicas per estimator; ρ_sample ∈ {0.05, 0.1, 0.25, 0.5} | Bias and SE per bin; decomposition of MSE into bias² and variance; same for HT and IS |
| E7 | RQ2: proposal comparison | IS with uniform, ρ_iρ_j, F_iF_j (W = 1, m, m²), ΔM^net N_iN_j, p*, cost-aware p* | η vs S for each proposal and snapshot |
| E8 | RQ3: shared draws for all k vs independent draws per k | IS with p*; same total cost | MSE and mass conservation per variant |
Estimators inside time integration
| ID | Question | Setup | Output |
|---|---|---|---|
| E9 | RQ4: sample handling | Frozen per macro step (thesis), frozen per integrator step, batch-averaged B draws, resampled per step | e_acc(t), lag ⟨ln m⟩(t), negativity events, step statistics |
| E10 | RQ4: implicit sampled form | Eq. 1.71 vs explicit sampled RHS inside backward Euler, Rosenbrock, Radau | Bias vs Δt, stability at large Δt |
| E11 | RQ5: noise at steady state | Full model to t = 10⁴ yr; S scan; Δt scan | Noise floor vs S and Δt (test H4); replica spread of ρ_k |
| E12 | RQ5: hybrid and balance-preserving estimators | q scan for the hybrid; stratification by outcome type (coag vs frag) | η and noise floor vs q; small-grain tail stability |
| E13 | RQ6: Jacobian strategies | Dense exact, sparse sampled, dominant-pair approximate, JFNK; with W-method and simplified-Newton Radau | Cost split RHS/Jacobian/solve; work-precision |
Scaling and comparison
| ID | Question | Setup | Output |
|---|---|---|---|
| E14 | RQ7, p = 1 | N_m ∈ {50, 100, 200, 400, 800}; with and without fragmentation and bouncing; best estimator and integrator from E7 to E13 | Work-precision diagrams; break-even N_m vs accuracy |
| E15 | RQ7, p = 2 | N_m × N_φ ∈ {50, 100} × {10, 20, 40}; full 2D sum as reference where affordable | Break-even surface; test of H5 |
| E16 | RQ8, comparison | Literature numbers and reproduced test problems for low-rank/TT and DG; RPMC via Beutel’s code if Dullemond can share it; no own implementations unless a referee asks | Comparison table and plot at matched accuracy |
3. Experiment hygiene
-
Headline experiments (E9 to E15) run in two physics setups (corrected thesis model, DustPy-like), each with and without bouncing.
-
Every run writes a manifest (config hash, git revision, seed, machine, timings). No figure in the paper without a manifest.
-
Experiments that use replicas use independent RNG streams per replica (ADR 0009) so results are reproducible and parallel.
-
Heavy runs go through
heavy(global rules) and are logged inexperiments/LOG.mdwith date, command, duration, and outcome.