From e395b5b0918a56cb274f7fb77ba0be073dd229e2 Mon Sep 17 00:00:00 2001 From: Meta Vens Date: Thu, 15 Oct 2026 21:43:03 +0200 Subject: [PATCH] Research: smooth (C-infinity) gravity wells, seam jumps, conserving halo and split-deepening probe Eyes pass eyes-ag3-smooth-wells (aetheria-gravity, self-2026-10-15-ag3). Extends gravity-tail-taming: smoothness class per law, measured A&S K0/K1/E1 seam jumps, log-ratio, screened Plummer, Lamb-Oseen (with and without halo), Plummer log and conserving-halo families on the same probe, plus zone far-field against merged mass. Co-Authored-By: Claude Opus 5.5 --- docs/research/gravity-smooth-wells.md | 129 ++++++++++++ .../gravity-smooth-wells/check-output.txt | 8 + .../gravity-smooth-wells/check_laws2.py | 32 +++ docs/research/gravity-smooth-wells/laws2.py | 66 ++++++ .../gravity-smooth-wells/probe2-output.txt | 195 ++++++++++++++++++ docs/research/gravity-smooth-wells/probe2.py | 150 ++++++++++++++ .../gravity-smooth-wells/seams-output.txt | 20 ++ docs/research/gravity-smooth-wells/seams.py | 35 ++++ 8 files changed, 635 insertions(+) create mode 100644 docs/research/gravity-smooth-wells.md create mode 100644 docs/research/gravity-smooth-wells/check-output.txt create mode 100644 docs/research/gravity-smooth-wells/check_laws2.py create mode 100644 docs/research/gravity-smooth-wells/laws2.py create mode 100644 docs/research/gravity-smooth-wells/probe2-output.txt create mode 100644 docs/research/gravity-smooth-wells/probe2.py create mode 100644 docs/research/gravity-smooth-wells/seams-output.txt create mode 100644 docs/research/gravity-smooth-wells/seams.py diff --git a/docs/research/gravity-smooth-wells.md b/docs/research/gravity-smooth-wells.md new file mode 100644 index 000000000..0c3398cd4 --- /dev/null +++ b/docs/research/gravity-smooth-wells.md @@ -0,0 +1,129 @@ +# Smooth wells: C-infinity laws for the grid and the rails + +Eyes pass `eyes-ag3-smooth-wells`, campaign aetheria-gravity, session self-2026-10-15-ag3, 2026-10-15, read at master 39ebf7f6. Facts with evidence pointers; no recommendation. Extends `gravity-tail-taming.md` (same data, metrics, phases and table format) under question `aetheria-gravity:question:tail-taming-law` and ruling `aetheria-gravity:ruling:tame-tails-parent-only-rails`. + +Operator, verbatim: "Since both orbits and well shape are negotiable, is there a solution available that's more mathematically elegant? Having all the derivatives shake out nicely isn't just an affectation, the player can feel higher-order discontinuities especially in a game where moving with the grid is such a huge part of the experience." + +Scripts in `gravity-smooth-wells/`: +- `laws2.py`: the laws (force and potential). It imports c, the core radius, today's law and the A&S K0/K1 from `../gravity-tail-taming/laws.py`. +- `check_laws2.py`: force = dφ/dr for each law, the shader form of Ein against mpmath, and circular-speed shapes. Output: `check-output.txt` (first lines) and the text below. +- `seams.py`: jumps in φ…φ''' at each piecewise seam, with each branch differentiated analytically in mpmath at 50 digits. Output: `seams-output.txt`. +- `probe2.py`: `../gravity-tail-taming/probe.py` with whole-law pulls. Output: `probe2-output.txt`. + +Every law has c = 45.08-calibrated strength (c_sun = 304.9 at the median sun) and ε = today's Roche core radius 7.32·m^0.25. Law 0 of `probe2.py` reproduces the tail-taming baseline exactly (Sun ratio p50 2.01, mid-zone 1.59 … 0.51). Each law ran in 9-16 s on Starfire, one process per law. + +## 1. Smoothness class + +A radial potential φ(r) is C^∞ on the plane exactly when it is a smooth function of r². Each law below is a function of r² with no seam, so it is real-analytic everywhere, the core included. The A&S polynomial forms are piecewise, so "the HLSL form" is a separate column. + +| Law | φ on the plane (exact) | HLSL/C# form | Seam location | +|---|---|---|---| +| Today: harmonic core + log | C^1 | C^1 | r = R. φ'' jumps by 2c/R² (inside +c/R², outside −c/R²); φ''' jumps by 2c/R³ | +| Tail-taming K0 rows (harmonic core + K0) | C^1 (core seam) | C^1 at R, plus the A&S seam | r = R, and r = 2λ | +| Tail-taming cap | C^1 at R and at L (φ'' jumps by c/L² at L) | same | | +| Tail-taming force switch | C^1 at R; C^2 at L and 2L | same | | +| (a) log-ratio | C^ω | C^∞: elementary (log, div), no branch | none | +| (b) screened Plummer | C^ω | C^0 in the formula: A&S piecewise | √(r²+ε²) = 2λ | +| (c) Lamb-Oseen, with or without halo | C^ω (Ein is entire) | φ: A&S E1 is piecewise. Force c(1−e^(−r²/ε²))/r is elementary and C^∞ | φ seam at r = ε (u = r²/ε² = 1) | +| (d) Plummer log | C^ω | C^∞: one log, no branch | none | +| (e) conserving halo | C^ω | C^∞: two logs, no branch | none | + +**Measured seam jumps** (`seams-output.txt`; jump = high branch − low branch; "rel" = |jump| / |true derivative|): + +| Seam | φ | φ' | φ'' | φ''' | +|---|---|---|---|---| +| A&S K0 at x = 2, in K0^(n) (multiply by −c/λ^n) | 7.6e-9 (rel 6.7e-8) | 2.5e-7 (rel 1.8e-6) | 8.4e-7 (rel 4.6e-6) | 7.7e-5 (rel 2.9e-4) | +| A&S K1 at x = 2 (force path), in K1^(n) (multiply by c/λ^(n+1)) | 1.5e-9 (rel 1.1e-8) | 1.9e-7 (rel 1.0e-6) | 5.1e-6 (rel 1.9e-5) | 4.6e-5 (rel 1.0e-4) | +| A&S E1 → Lamb-Oseen φ at r = ε, units c/ε^n | 2.2e-9 (rel 5.5e-9) | 2.2e-6 (rel 3.5e-6) | 1.3e-4 (rel 1.2e-3) | 5.1e-3 (rel 5.4e-3) | +| Today's core at r = R, units c/R^n | 0 | 0 | 2.0 (rel 2.0) | 2.0 (rel 1.0) | + +On the same units (c/ε², with ε = R), the Lamb-Oseen shader seam's φ'' jump is 6.4e-5 of today's. The K0 seam's φ'' jump is 8.4e-7·c/λ². With λ ≥ 0.5·D ≫ R, that is below 1e-6 of today's 2c/R². The K0 value jump at the seam (6.7e-8 relative) is about one float32 ulp of K0(2) = 0.114 (float32 ε = 6e-8). + +## 2. Probe (same measures as `gravity-tail-taming.md` Question 4; 6 draws, rng seed 7) + +**Family (e), conserving halo** (added on the coordinator's brief, after the operator's observation, verbatim: "this is all a consequence of exaggerating the contributions of small bodies... a well getting deeper the more its bodies are split by the zone generator. Two systems, same mass, one of them spent it on planets, that one goes waaaay deeper."): +- φ_i = (c_near/2)·ln(r²+ε²) − ((c_near − c_far)/2)·ln(r²+L²). +- c_near is today's compressed 45.08-calibrated m_i^0.2 strength. +- c_far,i = c(M_F)·m_i/M_F, where F is the body's parent family (the parent and all its descendants; for the root, the whole zone) and only non-belt, non-empty bodies carry mass. +- L = k·ext, where ext is the body's own family extent: the longest summed-distance path from the body to any descendant. A body with no descendants falls back to D (Z for a root). +- c_far/c_near = (m_i/M_F)^0.8 ≤ 1, so the halo strength is never negative. + +**Far-field measure** (every law): the inward pull at 8 points on a circle of radius 2Z (and 10Z) round the zone centre, times the radius, averaged over the 6 draws. It is divided by c(M_zone), the strength of the zone's body mass merged into one body under the same calibration. + +"outer" = the body's outermost direct-child distance (any child type, belts and empties included); a body with no children falls back to D (Z for a root). The Sun class in "Ratio" and in the rail-speed line is **Sun-class children, i.e. suns orbiting a sun** (`probe.py` keys by the child's class). `gravity-tail-taming.md` labels the same numbers "children of root and rosette suns", and that label does not match its code. + +| # | Law | Ratio p50 / p90: Sun · GasGiant · Planet · Planetoid · Belt | Mid-zone inward pull, Z = 1000 / 2125 / 4182 / 6846 / 10000 (u/s²) | Planets round a median sun, v p10/p50/p90 (u/s); share in 12-18 | Sun-class children's parent-only rail v p10/p50/p90 | Children with no parent pull | Children with parent pull < 0.5·c/D | Zone far-field c / c(M_zone), p50/p90 at 2Z; at 10Z | +|---|---|---|---|---|---|---|---|---| +| 0 | today (harmonic core + log) | 2.01/2.46 · 1.62/2.21 · 1.24/1.71 · 1.11/1.43 · 1.18/1.55 | 1.59 / 1.27 / 0.91 / 0.67 / 0.51 | 17.1/17.5/18.0; 0.92 | 22.1/23.6/24.6 | 0 | 46 | 5.54/7.98; 5.54/7.98 | +| d | Plummer log, no taming | 2.01/2.46 · 1.63/2.22 · 1.25/1.72 · 1.13/1.44 · 1.19/1.55 | 1.57 / 1.26 / 0.91 / 0.67 / 0.51 | 16.1/17.2/17.8; 0.99 | 22.0/23.6/24.6 | 0 | 234 | 5.54/7.98; 5.54/7.98 | +| a1 | log-ratio L = 0.5·D | 1.17/1.57 · 1.04/1.37 · 0.99/1.19 · 0.99/1.16 · 0.99/1.24 | 0.34 / 0.24 / 0.15 / 0.11 / 0.08 | 14.4/16.3/17.1; 0.96 | 12.3/15.4/19.5 | 18 (L ≤ ε for 285 bodies: repulsive) | 5838 | 0.10/0.11; 0.00/0.00 | +| a2 | log-ratio L = 1·D | 1.26/1.55 · 1.12/1.43 · 1.02/1.26 · 1.01/1.20 · 1.02/1.26 | 0.64 / 0.45 / 0.29 / 0.20 / 0.15 | 15.5/16.9/17.5; 1.00 | 17.7/19.9/22.7 | 0 (L ≤ ε for 34 bodies) | 407 | 0.36/0.40; 0.02/0.02 | +| a3 | log-ratio L = 1.5·outer | 1.17/1.41 · 1.12/1.39 · 1.05/1.31 · 1.03/1.27 · 1.05/1.32 | 0.74 / 0.46 / 0.27 / 0.18 / 0.14 | 14.1/15.2/16.7; 1.00 | 18.4/20.6/23.1 | 0 (7 bodies) | 637 | 0.34/0.48; 0.02/0.02 | +| a4 | log-ratio L = 2·outer | 1.22/1.46 · 1.16/1.44 · 1.07/1.35 · 1.04/1.28 · 1.08/1.35 | 0.85 / 0.54 / 0.32 / 0.22 / 0.16 | 15.1/15.9/17.0; 1.00 | 19.8/22.0/23.7 | 0 (7 bodies) | 326 | 0.52/0.72; 0.03/0.04 | +| b1 | screened Plummer λ = 0.5·D | 1.16/1.49 · 1.04/1.34 · 0.99/1.19 · 0.99/1.16 · 1.00/1.22 | 0.41 / 0.29 / 0.18 / 0.13 / 0.09 | 14.7/16.3/17.0; 0.98 | 14.3/17.2/20.5 | 0 | 3610 | 0.07/0.08; 0.00/0.00 | +| b2 | screened Plummer λ = 1·D | 1.28/1.56 · 1.15/1.46 · 1.04/1.29 · 1.02/1.22 · 1.04/1.28 | 0.70 / 0.50 / 0.32 / 0.23 / 0.16 | 15.6/16.8/17.5; 1.00 | 18.7/20.8/22.9 | 0 | 297 | 0.45/0.50; 0.00/0.00 | +| c0 | Lamb-Oseen core, log tail | 2.01/2.46 · 1.62/2.21 · 1.24/1.71 · 1.11/1.43 · 1.18/1.54 | 1.59 / 1.27 / 0.91 / 0.67 / 0.51 | 17.0/17.4/18.0; 0.93 | 22.1/23.6/24.6 | 0 | 121 | 5.54/7.98; 5.54/7.98 | +| c1 | Lamb-Oseen + halo L = 0.5·D | 1.17/1.57 · 1.04/1.37 · 0.99/1.19 · 0.99/1.15 · 0.99/1.22 | 0.35 / 0.24 / 0.16 / 0.11 / 0.08 | 14.8/16.6/17.5; 0.97 | 12.3/15.4/19.6 | 1 | 4944 | 0.10/0.12; 0.00/0.00 | +| c2 | Lamb-Oseen + halo L = 1·D | 1.26/1.55 · 1.12/1.44 · 1.02/1.26 · 1.00/1.19 · 1.02/1.26 | 0.66 / 0.46 / 0.29 / 0.21 / 0.15 | 16.5/17.2/17.8; 0.99 | 17.8/20.0/22.8 | 0 | 146 | 0.36/0.40; 0.02/0.02 | +| c3 | Lamb-Oseen + halo L = 1.5·outer | 1.17/1.41 · 1.11/1.39 · 1.04/1.31 · 1.02/1.26 · 1.05/1.32 | 0.75 / 0.46 / 0.28 / 0.18 / 0.14 | 14.3/15.8/17.2; 1.00 | 18.4/20.6/23.1 | 0 | 136 | 0.34/0.48; 0.02/0.02 | +| c4 | Lamb-Oseen + halo L = 2·outer | 1.22/1.46 · 1.15/1.44 · 1.06/1.34 · 1.04/1.27 · 1.07/1.34 | 0.87 / 0.54 / 0.32 / 0.22 / 0.16 | 15.4/16.4/17.5; 1.00 | 19.8/22.0/23.7 | 0 | 122 | 0.52/0.73; 0.03/0.04 | +| e1 | conserving halo L = 1.5·ext | 1.25/1.50 · 1.17/1.44 · 1.06/1.33 · 1.04/1.28 · 1.07/1.34 | 0.86 / 0.58 / 0.37 / 0.25 / 0.19 | 14.5/15.9/16.9; 1.00 | 21.7/23.4/24.5 | 0 (7 bodies) | 431 | 1.42/1.99; 1.28/1.84 | +| e2 | conserving halo L = 2·ext | 1.30/1.55 · 1.20/1.49 · 1.09/1.37 · 1.05/1.30 · 1.09/1.36 | 0.95 / 0.64 / 0.40 / 0.27 / 0.21 | 15.3/16.3/17.1; 1.00 | 21.8/23.5/24.5 | 0 (7 bodies) | 270 | 1.52/2.07; 1.28/1.85 | + +Notes from `probe2-output.txt`: +- *Core alone (rows d, c0 against 0)*: replacing the harmonic core with a C^∞ core moves no ratio or mid-zone figure by more than 0.02. It lowers rail speeds for children near their parent's core. Planetoids' p10 drops 11.5 → 11.0 (Plummer) and 11.4 (Lamb-Oseen). The count of children with a parent pull below half the log value rises from 46 to 234 (Plummer) and 121 (Lamb-Oseen). +- *Screening with the smooth core (b) against the tail-taming K0 rows (harmonic core)*: λ = 0.5·D ratio p50s 1.16/1.04/0.99/0.99/1.00, against 1.15/1.04/0.99/0.99/0.99. λ = 1·D gives 1.28/1.15/1.04/1.02/1.04, against 1.28/1.14/1.03/1.01/1.03. The mid-zone figures differ by ≤ 0.03. The median-sun planets go from 16.7 to 16.3 (0.5·D) and from 17.2 to 16.8 (1·D). +- *L = 0.5·D*: L ≤ ε for 285 of 27,165 bodies, where the log-ratio force c·r·(L²−ε²)/(…) is repulsive. 18 children then have no inward parent pull, and the ratio n drops by 2-12 per class. +- *Split-deepening (far-field column)*: + - *Today, Plummer and Lamb-Oseen*: the zone's far field is Σc_i = 5.54× (p50) and 7.98× (p90) the merged body's c(M_zone), the same at 2Z and 10Z. This is the n^(1−p) effect the operator names. + - *Screened and halo laws (a, b, c')*: at 2Z the zone's far field is 0.07-0.52× c(M_zone), and at 10Z ≤ 0.04×. Their far field is neutral by construction. + - *(e)*: 1.42×/1.99× at 2Z and 1.28×/1.84× at 10Z (p50/p90), not 1.0. The static sum Σc_far,i / c(M_zone) over each zone gives the same 1.02/1.27/1.84 (p10/50/90; one root in each of 1500 zones). The excess comes from nesting. A moon's c_far is its share of its planet family's compressed c(M_F), and the compressed strength per unit mass is higher for a smaller family. So the shares sum to c(M_F) only within one family level, not across the zone. If every body instead took its share of the zone total, c(M_zone)·m_i/M_zone, the shares would sum to exactly c(M_zone). +- *(e) rails*: + - It keeps today's Sun-class rails (21.7/23.4/24.5 against 22.1/23.6/24.6), because c_near governs inside L. + - It cuts the Sun-class ratio p50 from 2.01 to 1.25 (1.5·ext) or 1.30 (2·ext), and the mid-zone pull from 1.59 to 0.86-0.95 at Z = 1000. + - Planets round a median sun: 15.9 and 16.3 p50, all in 12-18. +- *Planets round all suns*, share in 12-18 u/s: 0.51 today, 0.53-0.64 for the smooth laws (highest: log-ratio 1.5·outer at 0.64). + +## 3. Closed forms, circular speed and shader cost + +Pull F is toward the body; v² = r·F; s = √(r²+ε²); u = r²/ε². Cost is per well per texel. The force vector is F·Δ/r, with Δ the offset to the well, so "no sqrt" below means the vector form needs no square root. + +| Law | φ | F | v²(r) | Plateau v ≈ √c | φ cost (grid) | Force-vector cost | +|---|---|---|---|---|---|---| +| Today | c r²/2R² (r= 0.9] + print(f'{name:18s} v/sqrt(c): peak {v[i]:.3f} at r = {rr[i]:.0f} ({rr[i]/eps:.2f} eps); v >= 0.9 sqrt(c) on r in ' + f'[{band.min():.0f}, {band.max():.0f}]' if len(band) else f'{name:18s} peak {v[i]:.3f} at r = {rr[i]:.0f}; never 0.9') + for x in (eps, 2 * eps, 5 * eps): + pass +print('v/sqrt(c) at r = eps, 2eps, 3eps:') +for name, F in [('today', F_today), ('plummer', F_plummer), ('lamb', F_lamb)]: + print(f' {name:8s}', [round(float(np.sqrt(x * F(c, eps, x, 0) / c)), 3) for x in (eps, 2 * eps, 3 * eps)]) +print('max pull (units c/eps):', {n: round(float(np.max(F(c, eps, rr, 0)) * eps / c), 3) for n, F in [('today', F_today), ('plummer', F_plummer), ('lamb', F_lamb)]}) diff --git a/docs/research/gravity-smooth-wells/laws2.py b/docs/research/gravity-smooth-wells/laws2.py new file mode 100644 index 000000000..ae5a3ff64 --- /dev/null +++ b/docs/research/gravity-smooth-wells/laws2.py @@ -0,0 +1,66 @@ +# Smooth (C-infinity) one-well laws (aetheria-gravity, Eyes pass eyes-ag3-smooth-wells, 2026-10-15). +# Each law: pull F(c, eps, r, s) toward the body (u/s^2) for ALL r >= 0 (the core is part of the law, no seam), +# and potential phi(c, eps, r, s) with dphi/dr = F (phi up to a constant). eps = today's Roche core radius 7.32*m^0.25. +# Reuses ../gravity-tail-taming/laws.py for c, eps, the A&S K0/K1 and today's piecewise law. +import os, sys +import numpy as np +sys.path.insert(0, os.path.join(os.path.dirname(os.path.abspath(__file__)), '..', 'gravity-tail-taming')) +from laws import k0, k1, strength, core_radius, pull, tail_log, P0 # noqa: F401 + +EULER = 0.5772156649015329 + +# ---- exponential integral E1, A&S 5.1.53 (0 < u <= 1) and 5.1.54 (u >= 1): the shader form, piecewise at u = 1 ---- +def e1_lo(u): # E1(u) = -ln u + poly(u), |err| < 2e-7 + return -np.log(u) + (-0.57721566 + u*(0.99999193 + u*(-0.24991055 + u*(0.05519968 + u*(-0.00976004 + u*0.00107857))))) +def e1_hi(u): # u e^u E1(u) = P4/Q4, |err| < 2e-8 + p = 0.2677737343 + u*(8.6347608925 + u*(18.0590169730 + u*(8.5733287401 + u))) + q = 3.9584969228 + u*(21.0996530827 + u*(25.6329561486 + u*(9.5733223454 + u))) + return np.exp(-u) / u * p / q +def ein(u): + """Ein(u) = integral_0^u (1 - e^-t)/t dt = E1(u) + ln u + gamma; entire. Shader form: A&S pieces, seam at u = 1.""" + u = np.asarray(u, float) + ul = np.clip(u, 1e-300, 1.0); uh = np.maximum(u, 1.0) + lo = e1_lo(ul) + np.log(ul) + EULER # the -ln u cancels: lo is a pure polynomial in u + hi = e1_hi(uh) + np.log(uh) + EULER + return np.where(u <= 1, lo, hi) + +# ---- laws: F(c, eps, r, s) and phi(c, eps, r, s); s is the per-body scale (L or lambda), unused where noted ---- +def F_today(c, eps, r, s): # harmonic core to eps, c/r beyond (C1 potential at eps) + return pull(tail_log, c, eps, r, s) +def phi_today(c, eps, r, s): + r = np.asarray(r, float) + return np.where(r < eps, c * r * r / (2 * eps * eps), c * (np.log(np.maximum(r, 1e-300) / eps) + 0.5)) + +def F_plummer(c, eps, r, s): # (d) Plummer-softened log, no taming + return c * r / (r * r + eps * eps) +def phi_plummer(c, eps, r, s): + return 0.5 * c * np.log(r * r + eps * eps) + +def F_logratio(c, eps, r, L): # (a) neutralised Plummer: Plummer core minus Plummer halo of width L, same strength + r2 = r * r + return c * r * (L * L - eps * eps) / ((r2 + eps * eps) * (r2 + L * L)) +def phi_logratio(c, eps, r, L): + return 0.5 * c * np.log((r * r + eps * eps) / (r * r + L * L)) + +def F_splummer(c, eps, r, lam): # (b) screened Plummer: phi = -c K0(sqrt(r^2+eps^2)/lam) + q = np.sqrt(r * r + eps * eps) + return c / lam * k1(q / lam) * r / q +def phi_splummer(c, eps, r, lam): + return -c * k0(np.sqrt(r * r + eps * eps) / lam) + +def F_lamb(c, eps, r, s): # (c) Lamb-Oseen core (Gaussian source), log tail + r = np.asarray(r, float); rs = np.maximum(r, 1e-12) + return c * (-np.expm1(-(rs / eps) ** 2)) / rs +def phi_lamb(c, eps, r, s): + return 0.5 * c * ein((np.asarray(r, float) / eps) ** 2) + +def F_lamb_halo(c, eps, r, L): # (c') Lamb-Oseen core minus the log-ratio's Plummer halo at L + return F_lamb(c, eps, r, L) - c * r / (r * r + L * L) +def phi_lamb_halo(c, eps, r, L): + return phi_lamb(c, eps, r, L) - 0.5 * c * np.log(r * r + L * L) + +def F_conserve(c, eps, r, L, cf): # (e) conserving halo: compressed c near, mass-additive c_far beyond L + r2 = r * r + return c * r / (r2 + eps * eps) - (c - cf) * r / (r2 + L * L) +def phi_conserve(c, eps, r, L, cf): + return 0.5 * c * np.log(r * r + eps * eps) - 0.5 * (c - cf) * np.log(r * r + L * L) diff --git a/docs/research/gravity-smooth-wells/probe2-output.txt b/docs/research/gravity-smooth-wells/probe2-output.txt new file mode 100644 index 000000000..828413b13 --- /dev/null +++ b/docs/research/gravity-smooth-wells/probe2-output.txt @@ -0,0 +1,195 @@ +== [0] today: harmonic core + log (10s, 6 draws) + Sun ratio p10/50/90/99 [1.44, 2.01, 2.46, 2.81] parent-only v p10/50/90 [22.06, 23.62, 24.64] n 2211 + GasGiant ratio p10/50/90/99 [1.1, 1.62, 2.21, 2.67] parent-only v p10/50/90 [18.27, 22.2, 24.57] n 3566 + Planet ratio p10/50/90/99 [0.97, 1.24, 1.71, 2.16] parent-only v p10/50/90 [15.56, 17.38, 23.61] n 7445 + Planetoid ratio p10/50/90/99 [0.95, 1.11, 1.43, 1.76] parent-only v p10/50/90 [11.51, 14.48, 20.27] n 10200 + Belt ratio p10/50/90/99 [0.95, 1.18, 1.55, 2.01] parent-only v p10/50/90 [11.41, 14.74, 20.29] n 3134 + planets round suns: parent-only v p10/50/90 [16.3, 17.88, 23.65]; share in 12-18 u/s 0.51 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [17.06, 17.48, 17.98]; share in 12-18 0.92; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law today: harmonic core + log): 27165 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 42, 'Belt': 3, 'Planet': 1} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [3.69, 5.54, 7.98], at 10Z [3.69, 5.54, 7.98] + mid-zone inward pull by Z: 1000: 1.59, 2125: 1.27, 4182: 0.91, 6846: 0.67, 10000: 0.51 +== [1] (d) Plummer log, no taming (10s, 6 draws) + Sun ratio p10/50/90/99 [1.44, 2.01, 2.46, 2.81] parent-only v p10/50/90 [21.99, 23.6, 24.62] n 2211 + GasGiant ratio p10/50/90/99 [1.11, 1.63, 2.22, 2.67] parent-only v p10/50/90 [17.58, 22.08, 24.48] n 3566 + Planet ratio p10/50/90/99 [0.97, 1.25, 1.72, 2.15] parent-only v p10/50/90 [14.9, 16.9, 23.2] n 7445 + Planetoid ratio p10/50/90/99 [0.95, 1.13, 1.44, 1.77] parent-only v p10/50/90 [11.01, 13.75, 18.14] n 10200 + Belt ratio p10/50/90/99 [0.97, 1.19, 1.55, 2.01] parent-only v p10/50/90 [11.05, 14.24, 19.45] n 3134 + planets round suns: parent-only v p10/50/90 [15.76, 17.36, 23.33]; share in 12-18 u/s 0.58 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [16.07, 17.17, 17.78]; share in 12-18 0.99; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (d) Plummer log, no taming): 27165 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 214, 'Belt': 11, 'Planet': 9} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [3.69, 5.54, 7.98], at 10Z [3.69, 5.54, 7.98] + mid-zone inward pull by Z: 1000: 1.57, 2125: 1.26, 4182: 0.91, 6846: 0.67, 10000: 0.51 +== [2] (a) log-ratio L=0.5*D (9s, 6 draws) + Sun ratio p10/50/90/99 [0.89, 1.17, 1.57, 2.01] parent-only v p10/50/90 [12.26, 15.37, 19.45] n 2211 + GasGiant ratio p10/50/90/99 [0.86, 1.04, 1.37, 1.87] parent-only v p10/50/90 [12.55, 17.03, 21.52] n 3566 + Planet ratio p10/50/90/99 [0.89, 0.99, 1.19, 1.58] parent-only v p10/50/90 [11.11, 15.78, 21.98] n 7443 + Planetoid ratio p10/50/90/99 [0.87, 0.99, 1.16, 1.98] parent-only v p10/50/90 [7.68, 12.64, 17.63] n 10188 + Belt ratio p10/50/90/99 [0.87, 0.99, 1.24, 2.04] parent-only v p10/50/90 [7.29, 12.73, 17.86] n 3130 + planets round suns: parent-only v p10/50/90 [11.39, 16.23, 22.23]; share in 12-18 u/s 0.60 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [14.35, 16.25, 17.08]; share in 12-18 0.96; n 1116 + children with no parent pull at their distance: {'Planetoid': 12, 'Belt': 4, 'Planet': 2} + bodies with scale s <= eps (law (a) log-ratio L=0.5*D): 285 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 1902, 'Planet': 955, 'Belt': 643, 'GasGiant': 905, 'Sun': 1433} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.07, 0.1, 0.11], at 10Z [0.0, 0.0, 0.0] + mid-zone inward pull by Z: 1000: 0.34, 2125: 0.24, 4182: 0.15, 6846: 0.11, 10000: 0.08 +== [3] (a) log-ratio L=1*D (9s, 6 draws) + Sun ratio p10/50/90/99 [0.98, 1.26, 1.55, 1.82] parent-only v p10/50/90 [17.73, 19.93, 22.72] n 2211 + GasGiant ratio p10/50/90/99 [0.92, 1.12, 1.43, 1.79] parent-only v p10/50/90 [16.19, 19.44, 23.25] n 3566 + Planet ratio p10/50/90/99 [0.91, 1.02, 1.26, 1.58] parent-only v p10/50/90 [13.89, 16.58, 22.83] n 7445 + Planetoid ratio p10/50/90/99 [0.9, 1.01, 1.2, 1.61] parent-only v p10/50/90 [9.99, 13.32, 18.04] n 10200 + Belt ratio p10/50/90/99 [0.9, 1.02, 1.26, 1.75] parent-only v p10/50/90 [9.83, 13.7, 18.99] n 3134 + planets round suns: parent-only v p10/50/90 [14.31, 17.05, 22.93]; share in 12-18 u/s 0.61 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [15.52, 16.86, 17.49]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (a) log-ratio L=1*D): 34 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 329, 'Belt': 47, 'Planet': 31} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.25, 0.36, 0.4], at 10Z [0.01, 0.02, 0.02] + mid-zone inward pull by Z: 1000: 0.64, 2125: 0.45, 4182: 0.29, 6846: 0.20, 10000: 0.15 +== [4] (a) log-ratio L=1.5*outer (9s, 6 draws) + Sun ratio p10/50/90/99 [0.98, 1.17, 1.41, 1.7] parent-only v p10/50/90 [18.39, 20.55, 23.07] n 2211 + GasGiant ratio p10/50/90/99 [0.93, 1.12, 1.39, 1.83] parent-only v p10/50/90 [15.02, 19.76, 23.51] n 3566 + Planet ratio p10/50/90/99 [0.92, 1.05, 1.31, 1.62] parent-only v p10/50/90 [12.84, 15.44, 22.92] n 7445 + Planetoid ratio p10/50/90/99 [0.89, 1.03, 1.27, 1.58] parent-only v p10/50/90 [9.19, 12.35, 18.07] n 10200 + Belt ratio p10/50/90/99 [0.89, 1.05, 1.32, 1.72] parent-only v p10/50/90 [9.1, 12.48, 19.08] n 3134 + planets round suns: parent-only v p10/50/90 [13.69, 15.9, 23.0]; share in 12-18 u/s 0.64 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [14.13, 15.16, 16.65]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (a) log-ratio L=1.5*outer): 7 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 485, 'Belt': 98, 'GasGiant': 9, 'Planet': 45} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.26, 0.34, 0.48], at 10Z [0.01, 0.02, 0.02] + mid-zone inward pull by Z: 1000: 0.74, 2125: 0.46, 4182: 0.27, 6846: 0.18, 10000: 0.14 +== [5] (a) log-ratio L=2*outer (9s, 6 draws) + Sun ratio p10/50/90/99 [1.01, 1.22, 1.46, 1.74] parent-only v p10/50/90 [19.75, 21.97, 23.67] n 2211 + GasGiant ratio p10/50/90/99 [0.94, 1.16, 1.44, 1.85] parent-only v p10/50/90 [16.0, 20.82, 23.9] n 3566 + Planet ratio p10/50/90/99 [0.93, 1.07, 1.35, 1.64] parent-only v p10/50/90 [13.71, 15.89, 23.05] n 7445 + Planetoid ratio p10/50/90/99 [0.9, 1.04, 1.28, 1.57] parent-only v p10/50/90 [9.91, 12.89, 18.11] n 10200 + Belt ratio p10/50/90/99 [0.91, 1.08, 1.35, 1.71] parent-only v p10/50/90 [9.84, 13.14, 19.26] n 3134 + planets round suns: parent-only v p10/50/90 [14.61, 16.38, 23.13]; share in 12-18 u/s 0.63 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [15.13, 15.94, 16.98]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (a) log-ratio L=2*outer): 7 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 277, 'Belt': 36, 'Planet': 13} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.4, 0.52, 0.72], at 10Z [0.02, 0.03, 0.04] + mid-zone inward pull by Z: 1000: 0.85, 2125: 0.54, 4182: 0.32, 6846: 0.22, 10000: 0.16 +== [6] (b) screened Plummer lam=0.5*D (16s, 6 draws) + Sun ratio p10/50/90/99 [0.9, 1.16, 1.49, 1.85] parent-only v p10/50/90 [14.3, 17.2, 20.5] n 2211 + GasGiant ratio p10/50/90/99 [0.86, 1.04, 1.34, 1.78] parent-only v p10/50/90 [14.3, 17.74, 21.89] n 3566 + Planet ratio p10/50/90/99 [0.88, 0.99, 1.19, 1.54] parent-only v p10/50/90 [12.24, 15.92, 22.05] n 7445 + Planetoid ratio p10/50/90/99 [0.88, 0.99, 1.16, 1.75] parent-only v p10/50/90 [8.7, 12.75, 17.65] n 10200 + Belt ratio p10/50/90/99 [0.87, 1.0, 1.22, 1.85] parent-only v p10/50/90 [8.43, 12.9, 18.01] n 3134 + planets round suns: parent-only v p10/50/90 [12.61, 16.37, 22.22]; share in 12-18 u/s 0.63 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [14.71, 16.3, 17.03]; share in 12-18 0.98; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (b) screened Plummer lam=0.5*D): 285 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 1326, 'Planet': 530, 'Belt': 420, 'GasGiant': 442, 'Sun': 892} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.02, 0.07, 0.08], at 10Z [0.0, 0.0, 0.0] + mid-zone inward pull by Z: 1000: 0.41, 2125: 0.29, 4182: 0.18, 6846: 0.13, 10000: 0.09 +== [7] (b) screened Plummer lam=1*D (15s, 6 draws) + Sun ratio p10/50/90/99 [1.0, 1.28, 1.56, 1.81] parent-only v p10/50/90 [18.68, 20.76, 22.88] n 2211 + GasGiant ratio p10/50/90/99 [0.93, 1.15, 1.46, 1.8] parent-only v p10/50/90 [16.95, 19.89, 23.31] n 3566 + Planet ratio p10/50/90/99 [0.92, 1.04, 1.29, 1.6] parent-only v p10/50/90 [14.02, 16.57, 22.77] n 7445 + Planetoid ratio p10/50/90/99 [0.91, 1.02, 1.22, 1.57] parent-only v p10/50/90 [10.15, 13.32, 18.0] n 10200 + Belt ratio p10/50/90/99 [0.91, 1.04, 1.28, 1.7] parent-only v p10/50/90 [10.08, 13.75, 18.94] n 3134 + planets round suns: parent-only v p10/50/90 [14.58, 17.04, 22.86]; share in 12-18 u/s 0.61 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [15.57, 16.81, 17.45]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (b) screened Plummer lam=1*D): 34 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 265, 'Belt': 21, 'Planet': 11} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.22, 0.45, 0.5], at 10Z [0.0, 0.0, 0.0] + mid-zone inward pull by Z: 1000: 0.70, 2125: 0.50, 4182: 0.32, 6846: 0.23, 10000: 0.16 +== [8] (c) Lamb-Oseen core, log tail (10s, 6 draws) + Sun ratio p10/50/90/99 [1.44, 2.01, 2.46, 2.81] parent-only v p10/50/90 [22.06, 23.62, 24.64] n 2211 + GasGiant ratio p10/50/90/99 [1.1, 1.62, 2.21, 2.66] parent-only v p10/50/90 [18.2, 22.2, 24.57] n 3566 + Planet ratio p10/50/90/99 [0.97, 1.24, 1.71, 2.15] parent-only v p10/50/90 [15.5, 17.34, 23.59] n 7445 + Planetoid ratio p10/50/90/99 [0.95, 1.11, 1.43, 1.76] parent-only v p10/50/90 [11.42, 14.34, 20.11] n 10200 + Belt ratio p10/50/90/99 [0.96, 1.18, 1.54, 2.0] parent-only v p10/50/90 [11.39, 14.69, 20.25] n 3134 + planets round suns: parent-only v p10/50/90 [16.26, 17.82, 23.64]; share in 12-18 u/s 0.53 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [17.03, 17.44, 17.96]; share in 12-18 0.93; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (c) Lamb-Oseen core, log tail): 27165 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 116, 'Belt': 3, 'Planet': 2} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [3.69, 5.54, 7.98], at 10Z [3.69, 5.54, 7.98] + mid-zone inward pull by Z: 1000: 1.59, 2125: 1.27, 4182: 0.91, 6846: 0.67, 10000: 0.51 +== [9] (c) Lamb-Oseen + halo L=0.5*D (10s, 6 draws) + Sun ratio p10/50/90/99 [0.89, 1.17, 1.57, 2.02] parent-only v p10/50/90 [12.31, 15.4, 19.64] n 2211 + GasGiant ratio p10/50/90/99 [0.85, 1.04, 1.37, 1.89] parent-only v p10/50/90 [12.69, 17.53, 21.7] n 3566 + Planet ratio p10/50/90/99 [0.89, 0.99, 1.19, 1.6] parent-only v p10/50/90 [11.57, 16.25, 22.57] n 7445 + Planetoid ratio p10/50/90/99 [0.86, 0.99, 1.15, 1.97] parent-only v p10/50/90 [8.38, 13.27, 19.14] n 10199 + Belt ratio p10/50/90/99 [0.86, 0.99, 1.22, 2.06] parent-only v p10/50/90 [7.85, 13.19, 19.18] n 3134 + planets round suns: parent-only v p10/50/90 [11.79, 16.67, 22.86]; share in 12-18 u/s 0.55 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [14.76, 16.63, 17.49]; share in 12-18 0.97; n 1116 + children with no parent pull at their distance: {'Planetoid': 1} + bodies with scale s <= eps (law (c) Lamb-Oseen + halo L=0.5*D): 285 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 1321, 'Planet': 807, 'Belt': 512, 'GasGiant': 881, 'Sun': 1423} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.07, 0.1, 0.12], at 10Z [0.0, 0.0, 0.0] + mid-zone inward pull by Z: 1000: 0.35, 2125: 0.24, 4182: 0.16, 6846: 0.11, 10000: 0.08 +== [10] (c) Lamb-Oseen + halo L=1*D (11s, 6 draws) + Sun ratio p10/50/90/99 [0.98, 1.26, 1.55, 1.83] parent-only v p10/50/90 [17.77, 20.04, 22.76] n 2211 + GasGiant ratio p10/50/90/99 [0.91, 1.12, 1.44, 1.8] parent-only v p10/50/90 [16.88, 19.7, 23.44] n 3566 + Planet ratio p10/50/90/99 [0.91, 1.02, 1.26, 1.59] parent-only v p10/50/90 [14.35, 17.0, 23.35] n 7445 + Planetoid ratio p10/50/90/99 [0.89, 1.0, 1.19, 1.63] parent-only v p10/50/90 [10.37, 13.93, 19.83] n 10200 + Belt ratio p10/50/90/99 [0.9, 1.02, 1.26, 1.77] parent-only v p10/50/90 [10.18, 14.16, 19.85] n 3134 + planets round suns: parent-only v p10/50/90 [14.83, 17.51, 23.44]; share in 12-18 u/s 0.56 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [16.49, 17.16, 17.75]; share in 12-18 0.99; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (c) Lamb-Oseen + halo L=1*D): 34 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 131, 'Belt': 8, 'Planet': 7} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.25, 0.36, 0.4], at 10Z [0.01, 0.02, 0.02] + mid-zone inward pull by Z: 1000: 0.66, 2125: 0.46, 4182: 0.29, 6846: 0.21, 10000: 0.15 +== [11] (c) Lamb-Oseen + halo L=1.5*outer (10s, 6 draws) + Sun ratio p10/50/90/99 [0.98, 1.17, 1.41, 1.71] parent-only v p10/50/90 [18.44, 20.57, 23.1] n 2211 + GasGiant ratio p10/50/90/99 [0.93, 1.11, 1.39, 1.86] parent-only v p10/50/90 [15.63, 19.92, 23.6] n 3566 + Planet ratio p10/50/90/99 [0.92, 1.04, 1.31, 1.63] parent-only v p10/50/90 [13.32, 15.95, 23.4] n 7445 + Planetoid ratio p10/50/90/99 [0.88, 1.02, 1.26, 1.6] parent-only v p10/50/90 [9.56, 13.01, 19.9] n 10200 + Belt ratio p10/50/90/99 [0.88, 1.05, 1.32, 1.72] parent-only v p10/50/90 [9.47, 12.88, 19.93] n 3134 + planets round suns: parent-only v p10/50/90 [13.86, 16.58, 23.48]; share in 12-18 u/s 0.60 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [14.29, 15.83, 17.22]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (c) Lamb-Oseen + halo L=1.5*outer): 7 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 129, 'Belt': 5, 'Planet': 2} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.26, 0.34, 0.48], at 10Z [0.01, 0.02, 0.02] + mid-zone inward pull by Z: 1000: 0.75, 2125: 0.46, 4182: 0.28, 6846: 0.18, 10000: 0.14 +== [12] (c) Lamb-Oseen + halo L=2*outer (10s, 6 draws) + Sun ratio p10/50/90/99 [1.01, 1.22, 1.46, 1.74] parent-only v p10/50/90 [19.8, 21.98, 23.7] n 2211 + GasGiant ratio p10/50/90/99 [0.94, 1.15, 1.44, 1.9] parent-only v p10/50/90 [16.76, 20.91, 23.99] n 3566 + Planet ratio p10/50/90/99 [0.92, 1.06, 1.34, 1.65] parent-only v p10/50/90 [14.25, 16.36, 23.46] n 7445 + Planetoid ratio p10/50/90/99 [0.9, 1.04, 1.27, 1.58] parent-only v p10/50/90 [10.27, 13.51, 19.99] n 10200 + Belt ratio p10/50/90/99 [0.9, 1.07, 1.34, 1.73] parent-only v p10/50/90 [10.18, 13.59, 20.05] n 3134 + planets round suns: parent-only v p10/50/90 [14.88, 16.96, 23.53]; share in 12-18 u/s 0.59 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [15.36, 16.4, 17.45]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (c) Lamb-Oseen + halo L=2*outer): 7 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 117, 'Belt': 3, 'Planet': 2} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [0.4, 0.52, 0.73], at 10Z [0.02, 0.03, 0.04] + mid-zone inward pull by Z: 1000: 0.87, 2125: 0.54, 4182: 0.32, 6846: 0.22, 10000: 0.16 +== [13] (e) conserving halo L=1.5*ext (9s, 6 draws) + Sun ratio p10/50/90/99 [1.03, 1.25, 1.5, 1.73] parent-only v p10/50/90 [21.66, 23.43, 24.48] n 2211 + GasGiant ratio p10/50/90/99 [0.96, 1.17, 1.44, 1.88] parent-only v p10/50/90 [16.23, 21.96, 24.4] n 3566 + Planet ratio p10/50/90/99 [0.93, 1.06, 1.33, 1.63] parent-only v p10/50/90 [13.42, 16.32, 23.18] n 7445 + Planetoid ratio p10/50/90/99 [0.9, 1.04, 1.28, 1.6] parent-only v p10/50/90 [9.36, 12.54, 18.14] n 10200 + Belt ratio p10/50/90/99 [0.91, 1.07, 1.34, 1.75] parent-only v p10/50/90 [9.3, 12.64, 19.43] n 3134 + planets round suns: parent-only v p10/50/90 [14.29, 16.45, 23.32]; share in 12-18 u/s 0.63 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [14.52, 15.86, 16.92]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (e) conserving halo L=1.5*ext): 7 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 356, 'Belt': 54, 'Planet': 20, 'GasGiant': 1} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [1.18, 1.42, 1.99], at 10Z [1.03, 1.28, 1.84] + mid-zone inward pull by Z: 1000: 0.86, 2125: 0.58, 4182: 0.37, 6846: 0.25, 10000: 0.19 +== [14] (e) conserving halo L=2*ext (9s, 6 draws) + Sun ratio p10/50/90/99 [1.06, 1.3, 1.55, 1.79] parent-only v p10/50/90 [21.82, 23.49, 24.53] n 2211 + GasGiant ratio p10/50/90/99 [0.97, 1.2, 1.49, 1.92] parent-only v p10/50/90 [16.79, 22.01, 24.44] n 3566 + Planet ratio p10/50/90/99 [0.94, 1.09, 1.37, 1.66] parent-only v p10/50/90 [14.07, 16.43, 23.2] n 7445 + Planetoid ratio p10/50/90/99 [0.91, 1.05, 1.3, 1.58] parent-only v p10/50/90 [10.01, 13.0, 18.14] n 10200 + Belt ratio p10/50/90/99 [0.92, 1.09, 1.36, 1.75] parent-only v p10/50/90 [9.97, 13.26, 19.44] n 3134 + planets round suns: parent-only v p10/50/90 [14.91, 16.67, 23.32]; share in 12-18 u/s 0.62 + planets round a median sun (mass p40-p60, 10607-20254): v p10/50/90 [15.32, 16.31, 17.14]; share in 12-18 1.00; n 1116 + children with no parent pull at their distance: {} + bodies with scale s <= eps (law (e) conserving halo L=2*ext): 7 of 27165 + children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): {'Planetoid': 239, 'Belt': 22, 'Planet': 9} + zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z [1.24, 1.52, 2.07], at 10Z [1.03, 1.28, 1.85] + mid-zone inward pull by Z: 1000: 0.95, 2125: 0.64, 4182: 0.40, 6846: 0.27, 10000: 0.21 diff --git a/docs/research/gravity-smooth-wells/probe2.py b/docs/research/gravity-smooth-wells/probe2.py new file mode 100644 index 000000000..0c47c7bb9 --- /dev/null +++ b/docs/research/gravity-smooth-wells/probe2.py @@ -0,0 +1,150 @@ +# Smooth-well probe: ../gravity-tail-taming/probe.py with whole-law pulls F(c, eps, r, s) from laws2.py (core inside the +# law, eps = Roche core radius) and a per-body 'outer' scale (outermost direct child distance; falls back to D). +# Same data, metrics, phases (rng seed 7, 6 draws) and output format; law 0 must reproduce the tail-taming baseline. +# Usage: python probe2.py bodies.csv [draws] [law-index,...] +import csv, sys, os, collections, time +import numpy as np +sys.path.insert(0, os.path.dirname(os.path.abspath(__file__))) +from laws2 import * + +rng = np.random.default_rng(7) +DRAWS = int(sys.argv[2]) if len(sys.argv) > 2 else 6 +rows = [r for r in csv.DictReader(open(sys.argv[1])) if r.get('type')] +zones = collections.defaultdict(list) +for r in rows: zones[(r['density'], r['seed'])].append(r) + +def hill2(D, c, cpar): # 2 x B25's tidal radius d_H = D sqrt(c / (2 c_parent)); a root or bodyless-parent child: D + ok = cpar > 0 + return np.where(ok, 2 * D * np.sqrt(c / (2 * np.where(ok, cpar, 1.0))), D) + +W = lambda F: (lambda c, eps, r, s, cf: F(c, eps, r, s)) # laws without a far strength ignore cf +ONE = lambda m, D, Z, c, cp, o, e: np.ones_like(m) +# (name, F(c, eps, r, s), p, scale(m, D, Z, c, c_parent, outer, own family extent) -> per-body s, subtract disc) +LAWS = [ + ('today: harmonic core + log', W(F_today), 0.2, ONE, False), + ('(d) Plummer log, no taming', W(F_plummer), 0.2, ONE, False), + ('(a) log-ratio L=0.5*D', W(F_logratio), 0.2, lambda m, D, Z, c, cp, o, e: 0.5 * D, False), + ('(a) log-ratio L=1*D', W(F_logratio), 0.2, lambda m, D, Z, c, cp, o, e: 1.0 * D, False), + ('(a) log-ratio L=1.5*outer', W(F_logratio), 0.2, lambda m, D, Z, c, cp, o, e: 1.5 * o, False), + ('(a) log-ratio L=2*outer', W(F_logratio), 0.2, lambda m, D, Z, c, cp, o, e: 2.0 * o, False), + ('(b) screened Plummer lam=0.5*D', W(F_splummer), 0.2, lambda m, D, Z, c, cp, o, e: 0.5 * D, False), + ('(b) screened Plummer lam=1*D', W(F_splummer), 0.2, lambda m, D, Z, c, cp, o, e: 1.0 * D, False), + ('(c) Lamb-Oseen core, log tail', W(F_lamb), 0.2, ONE, False), + ('(c) Lamb-Oseen + halo L=0.5*D', W(F_lamb_halo), 0.2, lambda m, D, Z, c, cp, o, e: 0.5 * D, False), + ('(c) Lamb-Oseen + halo L=1*D', W(F_lamb_halo), 0.2, lambda m, D, Z, c, cp, o, e: 1.0 * D, False), + ('(c) Lamb-Oseen + halo L=1.5*outer', W(F_lamb_halo), 0.2, lambda m, D, Z, c, cp, o, e: 1.5 * o, False), + ('(c) Lamb-Oseen + halo L=2*outer', W(F_lamb_halo), 0.2, lambda m, D, Z, c, cp, o, e: 2.0 * o, False), + ('(e) conserving halo L=1.5*ext', F_conserve, 0.2, lambda m, D, Z, c, cp, o, e: 1.5 * e, False), + ('(e) conserving halo L=2*ext', F_conserve, 0.2, lambda m, D, Z, c, cp, o, e: 2.0 * e, False), +] +pick = [int(x) for x in sys.argv[3].split(',')] if len(sys.argv) > 3 else range(len(LAWS)) + +def zone_struct(z): + idx = {int(r['idx']): k for k, r in enumerate(z)} + n = len(z) + par = np.array([idx.get(int(r['parentIdx']), -1) for r in z]) + D = np.array([float(r['distance']) for r in z]); m = np.array([float(r['mass']) for r in z]) + body = np.array([r['empty'] == '0' and r['belt'] == '0' for r in z]) + cls = ['Belt' if r['belt'] == '1' else r['type'] for r in z] + fam = np.eye(n, dtype=bool) # fam[i, j]: j is i or a descendant of i + depth = np.zeros(n, int) + for j in range(n): + a = par[j] + while a >= 0: + fam[a, j] = True; depth[j] += 1; a = par[a] + order = list(np.argsort(depth, kind='stable')) + return par, D, m, body, cls, fam, order + +def place(par, D, order): + pos = np.zeros((len(par), 2)); th = rng.uniform(0, 2 * np.pi, len(par)) + for k in order: + if par[k] >= 0: pos[k] = pos[par[k]] + D[k] * np.array([np.cos(th[k]), np.sin(th[k])]) + return pos + +def field(pos_b, c, R, s, tail, x, cf): + """Per-body acceleration contributions at points x: (n, nb, 2), each toward its body.""" + d = pos_b[None, :, :] - x[:, None, :]; r = np.linalg.norm(d, axis=2) + rs = np.maximum(r, 1e-6) + F = np.where(r < 1e-6, 0, tail(c[None, :], R[None, :], rs, s[None, :], cf[None, :])) + return d / rs[:, :, None] * F[:, :, None] + +def disc_acc(x, C, Z): # field of a uniform disc of total strength C over radius Z centred at 0 (inward) + r = np.linalg.norm(x, axis=1, keepdims=True); rs = np.maximum(r, 1e-9) + return -x / rs * np.where(r < Z, C * r / Z ** 2, C / rs) + +def q(v, ps=(10, 50, 90)): + v = np.sort(np.asarray(v)) + return [round(float(v[min(len(v) - 1, int(p / 100 * len(v)))]), 2) for p in ps] if len(v) else [] + +pre = {k: zone_struct(z) for k, z in zones.items()} +sunm = np.sort([float(r['mass']) for r in rows if r['type'] == 'Sun' and r['empty'] == '0']) +SUN_LO, SUN_HI = sunm[int(.4 * len(sunm))], sunm[int(.6 * len(sunm))] # 'median sun': mass between p40 and p60 +for li in pick: + name, tail, p, scale, disc = LAWS[li] + t0 = time.time() + ratio = collections.defaultdict(list); vpar = collections.defaultdict(list); vsp = []; vmed = [] + inward = collections.defaultdict(list); farc = {2: [], 10: []}; nLeR = [0, 0]; weak = collections.Counter() + for key, z in zones.items(): + par, D, m, body, cls, fam, order = pre[key] + Z = float(z[0]['zoneRadius']) + c = np.where(body, strength(m, p), 0.0); R = core_radius(m) + cpar = np.array([c[par[k]] if par[k] >= 0 else 0.0 for k in range(len(par))]) + Dz = np.where(par < 0, Z, D) + outer = np.array([max([D[j] for j in range(len(par)) if par[j] == k], default=Dz[k]) for k in range(len(par))]) + # own family extent: max path length from the body to any descendant; no descendants -> D (Z for a root) + n_ = len(par); plen = np.zeros((n_, n_)) + for j in range(n_): + a, acc = j, 0.0 + while par[a] >= 0: + acc += D[a]; a = par[a]; plen[a, j] = acc + ext = np.where(plen.max(axis=1) > 0, plen.max(axis=1), Dz) + # c_far: the body's linear mass share of its parent family's compressed strength (root: the whole zone) + mb = np.where(body, m, 0.0) + MF = np.array([mb[fam[par[k]]].sum() if par[k] >= 0 else mb.sum() for k in range(n_)]) + cf = np.where(body, strength(np.maximum(MF, 1e-9), p) * mb / np.maximum(MF, 1e-9), 0.0) + s = scale(m, Dz, Z, c, cpar, outer, ext) + nLeR[0] += int(np.sum(body & (s <= R))); nLeR[1] += int(np.sum(body)) + bi = np.where(body)[0]; bpos = {j: t for t, j in enumerate(bi)} + kids = [k for k in range(len(par)) if par[k] >= 0 and body[par[k]]] + a_own = {k: float(tail(c[par[k]], R[par[k]], D[k], s[par[k]], cf[par[k]])) for k in kids} + for k in kids: + v = np.sqrt(max(a_own[k], 0) * D[k]); vpar[cls[k]].append(v) + if a_own[k] < 0.5 * c[par[k]] / D[k]: weak[cls[k]] += 1 + if cls[k] == 'Planet' and z[par[k]]['type'] == 'Sun': + vsp.append(v) + if SUN_LO <= m[par[k]] <= SUN_HI: vmed.append(v) + acc_rel = collections.defaultdict(float); mids = []; farc[2].append(0.0); farc[10].append(0.0) + fb = fam[:, bi][:, :, None] + for dr in range(DRAWS): + pos = place(par, D, order) + A = field(pos[bi], c[bi], R[bi], s[bi], tail, pos, cf[bi]) + Aout = np.where(fb, 0, A).sum(axis=1) # at node k, from bodies outside k's family + if disc: Aout = Aout - disc_acc(pos, c.sum(), Z) + for k in kids: + pk = par[k]; u = pos[pk] - pos[k]; u /= np.linalg.norm(u) + acc_rel[k] += np.dot(Aout[k] - Aout[pk], u) / DRAWS # Aout[k] includes the parent's own pull + th = np.arange(8) * np.pi / 4; X = 0.5 * Z * np.stack([np.cos(th), np.sin(th)], 1) + AX = field(pos[bi], c[bi], R[bi], s[bi], tail, X, cf[bi]).sum(axis=1) + for fz in (2, 10): # far-field strength: mean inward pull x radius at 2Z and 10Z, over the merged body's c + XF = fz * Z * np.stack([np.cos(th), np.sin(th)], 1) + AF = field(pos[bi], c[bi], R[bi], s[bi], tail, XF, cf[bi]).sum(axis=1) + farc[fz][-1] += np.mean(-(AF * XF).sum(1)) / DRAWS / strength(mb.sum(), p) + if disc: AX = AX - disc_acc(X, c.sum(), Z) + mids.append(np.mean(-(AX * X).sum(1) / np.linalg.norm(X, axis=1))) + inward[Z].append(np.mean(mids)) + for k in kids: + if a_own[k] > 0: ratio[cls[k]].append(np.sqrt(max(acc_rel[k], 0) / a_own[k])) + vsp = np.array(vsp) + print(f'== [{li}] {name} ({time.time() - t0:.0f}s, {DRAWS} draws)') + for k in ['Sun', 'GasGiant', 'Planet', 'Planetoid', 'Belt']: + print(f' {k:9s} ratio p10/50/90/99 {q(ratio[k], (10, 50, 90, 99))} parent-only v p10/50/90 {q(vpar[k])} n {len(ratio[k])}') + print(f' planets round suns: parent-only v p10/50/90 {q(vsp)}; share in 12-18 u/s {np.mean((vsp >= 12) & (vsp <= 18)):.2f}') + vmed = np.array(vmed) + print(f' planets round a median sun (mass p40-p60, {SUN_LO:.0f}-{SUN_HI:.0f}): v p10/50/90 {q(vmed)}; share in 12-18 {np.mean((vmed >= 12) & (vmed <= 18)):.2f}; n {len(vmed)}') + nz = {k: sum(1 for x in vpar[k] if x == 0) for k in vpar} + print(' children with no parent pull at their distance: ' + str({k: v for k, v in nz.items() if v})) + print(f' bodies with scale s <= eps (law {name}): {nLeR[0]} of {nLeR[1]}') + print(' children whose parent pull < 0.5 c/D (v < 0.71 of the log rail): ' + str(dict(weak))) + print(' zone far-field c over merged-body c(M_zone), p10/50/90: at 2Z ' + str(q(farc[2])) + ', at 10Z ' + str(q(farc[10]))) + print(' mid-zone inward pull by Z: ' + ', '.join(f'{Z:.0f}: {np.median(v):.2f}' for Z, v in sorted(inward.items()))) + sys.stdout.flush() diff --git a/docs/research/gravity-smooth-wells/seams-output.txt b/docs/research/gravity-smooth-wells/seams-output.txt new file mode 100644 index 000000000..a6cc10c6b --- /dev/null +++ b/docs/research/gravity-smooth-wells/seams-output.txt @@ -0,0 +1,20 @@ +A&S K0 at x = 2 (phi = -c K0(r/lam); jumps in K0^(n), multiply by -c/lam^n): + d^0: jump -7.58e-9 1/lam^0 true value 0.113894 |jump|/|true| 6.65e-8 + d^1: jump -2.48e-7 1/lam^1 true value -0.139866 |jump|/|true| 1.77e-6 + d^2: jump 8.43e-7 1/lam^2 true value 0.183827 |jump|/|true| 4.59e-6 + d^3: jump -7.68e-5 1/lam^3 true value -0.266746 |jump|/|true| 0.000288 +A&S K1 at x = 2 (force = (c/lam) K1(r/lam)): + d^0: jump 1.48e-9 1/lam^0 true value 0.139866 |jump|/|true| 1.06e-8 + d^1: jump 1.92e-7 1/lam^1 true value -0.183827 |jump|/|true| 1.04e-6 + d^2: jump -5.14e-6 1/lam^2 true value 0.266746 |jump|/|true| 1.93e-5 + d^3: jump 4.64e-5 1/lam^3 true value -0.44408 |jump|/|true| 0.000105 +A&S E1 at u = 1 (phi = (c/2) Ein(r^2/eps^2), seam at r = eps; jumps in units c/eps^n): + d^0: jump 2.2e-9 1/eps^0 true value 0.3983 |jump|/|true| 5.52e-9 + d^1: jump -2.19e-6 1/eps^1 true value 0.632121 |jump|/|true| 3.46e-6 + d^2: jump -0.000128 1/eps^2 true value 0.103638 |jump|/|true| 0.00124 + d^3: jump -0.00512 1/eps^3 true value -0.943036 |jump|/|true| 0.00543 +Today's core at r = R (phi: r^2/2R^2 inside, ln(r/R)+1/2 beyond; units c/R^n): + d^0: jump 0.0 1/R^0 true value 0.5 |jump|/|true| 0.0 + d^1: jump 0.0 1/R^1 true value 1.0 |jump|/|true| 0.0 + d^2: jump -2.0 1/R^2 true value -1.0 |jump|/|true| 2.0 + d^3: jump 2.0 1/R^3 true value 2.0 |jump|/|true| 1.0 diff --git a/docs/research/gravity-smooth-wells/seams.py b/docs/research/gravity-smooth-wells/seams.py new file mode 100644 index 000000000..f7f648272 --- /dev/null +++ b/docs/research/gravity-smooth-wells/seams.py @@ -0,0 +1,35 @@ +# Seam jumps of the piecewise shader forms (A&S K0/K1 at x = 2; A&S E1 at u = 1) and of today's core at R. +# Each branch is evaluated as its own analytic function in mpmath (50 digits) and differentiated at the seam, +# so the jump is the branch mismatch itself, not finite-difference noise. +import mpmath as mp +mp.mp.dps = 50 +def i0(x): t = (x / 3.75) ** 2; return 1 + t*(3.5156229 + t*(3.0899424 + t*(1.2067492 + t*(0.2659732 + t*(0.0360768 + t*0.0045813))))) +def i1(x): t = (x / 3.75) ** 2; return x * (0.5 + t*(0.87890594 + t*(0.51498869 + t*(0.15084934 + t*(0.02658733 + t*(0.00301532 + t*0.00032411)))))) +def k0_lo(x): y = (x / 2) ** 2; return -mp.log(x / 2) * i0(x) + (-0.57721566 + y*(0.42278420 + y*(0.23069756 + y*(0.03488590 + y*(0.00262698 + y*(0.00010750 + y*0.00000740)))))) +def k0_hi(x): z = 2 / x; return mp.exp(-x) / mp.sqrt(x) * (1.25331414 + z*(-0.07832358 + z*(0.02189568 + z*(-0.01062446 + z*(0.00587872 + z*(-0.00251540 + z*0.00053208)))))) +def k1_lo(x): y = (x / 2) ** 2; return mp.log(x / 2) * i1(x) + (1 / x) * (1 + y*(0.15443144 + y*(-0.67278579 + y*(-0.18156897 + y*(-0.01919402 + y*(-0.00110404 + y*(-0.00004686))))))) +def k1_hi(x): z = 2 / x; return mp.exp(-x) / mp.sqrt(x) * (1.25331414 + z*(0.23498619 + z*(-0.03655620 + z*(0.01504268 + z*(-0.00780353 + z*(0.00325614 + z*(-0.00068245))))))) +def e1_lo(u): return -mp.log(u) + (-0.57721566 + u*(0.99999193 + u*(-0.24991055 + u*(0.05519968 + u*(-0.00976004 + u*0.00107857))))) +def e1_hi(u): + p = 0.2677737343 + u*(8.6347608925 + u*(18.0590169730 + u*(8.5733287401 + u))) + q = 3.9584969228 + u*(21.0996530827 + u*(25.6329561486 + u*(9.5733223454 + u))) + return mp.exp(-u) / u * p / q + +def report(title, lo, hi, true, x0, unit): + print(title) + for n in range(4): + a, b, t = mp.diff(lo, x0, n), mp.diff(hi, x0, n), mp.diff(true, x0, n) + j = b - a + print(f' d^{n}: jump {mp.nstr(j, 3):>10s} {unit}^{n} true value {mp.nstr(t, 6):>10s} |jump|/|true| {mp.nstr(abs(j / t), 3)}') + +# phi = -c K0(r/lam): phi^(n) jump = -c/lam^n * (K0_hi^(n) - K0_lo^(n)) at x = 2, i.e. at r = 2 lam +report('A&S K0 at x = 2 (phi = -c K0(r/lam); jumps in K0^(n), multiply by -c/lam^n):', k0_lo, k0_hi, lambda x: mp.besselk(0, x), 2, '1/lam') +report('A&S K1 at x = 2 (force = (c/lam) K1(r/lam)):', k1_lo, k1_hi, lambda x: mp.besselk(1, x), 2, '1/lam') +# Lamb-Oseen phi = (c/2) Ein(r^2/eps^2): seam at r = eps. Differentiate in r (eps = 1, c = 1) -> units c/eps^n +ein_lo = lambda r: (e1_lo(r * r) + mp.log(r * r) + mp.euler) / 2 +ein_hi = lambda r: (e1_hi(r * r) + mp.log(r * r) + mp.euler) / 2 +ein_true = lambda r: (mp.e1(r * r) + mp.log(r * r) + mp.euler) / 2 +report('A&S E1 at u = 1 (phi = (c/2) Ein(r^2/eps^2), seam at r = eps; jumps in units c/eps^n):', ein_lo, ein_hi, ein_true, 1, '1/eps') +# today's core at R (c = 1, R = 1): inside r^2/2, outside ln r + 1/2 +report("Today's core at r = R (phi: r^2/2R^2 inside, ln(r/R)+1/2 beyond; units c/R^n):", + lambda r: r * r / 2, lambda r: mp.log(r) + 0.5, lambda r: mp.log(r) + 0.5, 1, '1/R')