Velocity distribution functions - progenax — differentiable initial conditions for N-body simulations

Velocity distribution functions — model cards¶

Machine-generated from the provenance registry (ADR-0034): every model’s sources (with public DOI/arXiv/ADS pointers), governing equations, parameter meanings, code entry points, and the validation tests that pin them. The hand-authored theory pages hold the derivations; these cards are the citable single source of truth.

Plummer isotropic DF¶

✅ verified · 3 equations · 3 sources · API: PlummerVelocityDF

The closed-form ergodic equilibrium DF of the Plummer sphere, f(E) ∝ E^(7/2) — the canonical worked example of Eddington inversion and the exact equilibrium behind progenax’s Plummer velocity sampler.

Use it for	Not for
✓ exact detailed equilibrium for Plummer positions (Q_vir = 0.5 with NO rescale)	✗ tidally truncated systems (use the King DF)
✓ pedagogical baseline: the one DF whose entire chain is closed-form	✗ intrinsic radial anisotropy (add anisotropy_radius, or use Michie for King-like systems)

Parameters¶

Name	Meaning	Units	Typical range	Code
`r_h`	3-D half-mass radius of the matched Plummer profile (must equal the profile’s)	pc	0.1–10	`PlummerVelocityDF(r_h=...)`
`anisotropy_radius`	optional Osipkov-Merritt anisotropy radius r_a; None = isotropic	pc	>= 0.75 a (Merritt 1985 non-negativity bound)	`PlummerVelocityDF(r_h=..., anisotropy_radius=...)`

Equations¶

f(\mathcal{E}) \;=\; \frac{1}{\sqrt{8}\,\pi^2} \biggl[\,\int_0^{\mathcal{E}} \frac{\mathrm{d}^2\rho}{\mathrm{d}\Psi^2}\,\frac{\mathrm{d}\Psi}{\sqrt{\mathcal{E} - \Psi}} + \frac{1}{\sqrt{\mathcal{E}}}\,\biggl(\frac{\mathrm{d}\rho}{\mathrm{d}\Psi}\biggr)_{\!\Psi=0}\,\biggr]

(1)

Symbols: $\mathcal{E}$ : binding energy -E = Psi - v^2/2 [pc^2 / Myr^2]; $\Psi$ : relative potential -Phi (positive, -> 0 at infinity) [pc^2 / Myr^2]. Assumes: spherical; isotropic (ergodic f(E)). Derivation: theory page.

f(\mathcal{E}) \;=\; \frac{24\sqrt{2}}{7\pi^3}\,\frac{a^2}{G^5\,M^4}\,\mathcal{E}^{7/2},\qquad \mathcal{E} > 0

(2)

Symbols: $a$ : Plummer scale radius [pc]; $M$ : total mass [Msun]. Assumes: ergodic; isotropic; coefficient independently confirmed by Dejonghe 1987 Eq. 16 (3*2^{7/2} = 24 sqrt 2). Derivation: theory page.

f(v \mid r) \;\propto\; v^2\,\mathcal{E}(r, v)^{7/2},\qquad \mathcal{E}(r, v) = \frac{GM}{\sqrt{r^2+a^2}} - \tfrac{1}{2}v^2

(3)

Symbols: $v$ : speed at radius r, in [0, v_esc(r)] [pc / Myr]. Assumes: marginal speed distribution at fixed r (the inverse-CDF sampling target). Derivation: theory page.

Sources¶

Plummer (1911) — the density model the DF inverts (Eqs. 11-12)
Dejonghe (1987) — closed-form cross-checks: isotropic f(E) coefficient (3/7pi^3)(2E)^{7/2} = 24sqrt(2)/(7pi^3) E^{7/2}; sigma_r^2 = (1/6)(1+r^2)^{-1/2}; projected sigma_p^2 = (3pi/64)/sqrt(1+r_p^2) at q=0 (Eq. 16 (DF); Eq. 17 (dispersion); Eq. 43 with q=0 (projected dispersion), p. 24)
Aarseth et al. (1974) — the N-body sampling recipe (inverse-CDF positions, Beta(3/2,9/2) speed draw) (Eq. 1, p. 183)

Code & validation¶

code: PlummerVelocityDF (src/progenax/kinematics/plummer_df.py)
validation: tests/validation/test_plummer_physics.py::TestPlummerVirialEquilibrium::test_virial_ratio
validation: tests/validation/test_plummer_physics.py::TestPlummerVelocityDispersion::test_radial_dispersion_profile
validation: tests/validation/test_plummer_physics.py::TestPlummerBetaDistribution::test_q_squared_mean
validation: tests/validation/test_dispersion_physics.py::test_projection_isotropic_plummer_los_oracle

King lowered-Maxwellian DF¶

✅ verified · 3 equations · 1 source · API: KingVelocityDF

The velocity side of the King (1966) model: an isothermal Maxwellian lowered by a constant so f vanishes at the tidal-boundary energy — sampled per-particle by inverse-CDF on the marginal speed distribution.

Use it for	Not for
✓ true detailed equilibrium for King positions (unscaled Q_vir ~ 0.5)	✗ radial anisotropy (KingVelocityDF is isotropic — use Michie)
✓ sharply truncated clusters with cold outskirts (sigma -> 0 at r_t)	✗ untruncated smooth systems (use Plummer)

Parameters¶

Name	Meaning	Units	Typical range	Code
`W0`	dimensionless central potential (Phi_t - Phi(0))/sigma_0^2 — must match the profile’s	dimensionless	3–12	`KingVelocityDF(W0=...)`
`r_c`	King core radius (sets sigma_0 self-consistently)	pc	0.1–5	`KingVelocityDF(r_c=...)`

Equations¶

f(E) \;=\; \begin{cases} \rho_1\,(2\pi\sigma_0^2)^{-3/2}\,\Bigl[\,e^{(\Phi_t - E)/\sigma_0^2} - 1\,\Bigr], & E < \Phi_t \\ 0, & E \ge \Phi_t \end{cases}

(4)

Symbols: $\sigma_0$ : central 1-D velocity-scale parameter [pc / Myr]; $\Phi_t$ : potential at the tidal radius (the cutoff energy) [pc^2 / Myr^2]. Assumes: isotropic; single-mass. Derivation: theory page.

f(v \mid W) \;\propto\; v^2\,\Bigl[\,e^{W - v^2 / (2\sigma_0^2)} - 1\,\Bigr], \qquad 0 \le v \le \sigma_0\sqrt{2W}

(5)

Symbols: $W$ : local dimensionless potential [dimensionless]. Assumes: marginal speed distribution at fixed W (the inverse-CDF sampling target). Derivation: theory page.

\sigma^2(W) \;=\; \sigma_0^2\,\frac{\int_0^{\sqrt{2W}} v^4\,(e^{W - v^2/2} - 1)\,\mathrm{d}v}{\int_0^{\sqrt{2W}} v^2\,(e^{W - v^2/2} - 1)\,\mathrm{d}v}

(6)

Symbols: $\sigma$ : mean-square speed of the lowered Maxwellian at W (-> 3 sigma_0^2 as W -> inf) [pc / Myr]. Assumes: decreases outward; vanishes at the tidal boundary. Derivation: theory page.

Sources¶

King (1966) — the lowered-Maxwellian DF and its W-form (Eq. 3 (DF); Eq. 9 (dimensionless W); Eq. 17 (density moment))

Code & validation¶

code: KingVelocityDF (src/progenax/kinematics/king_df.py)
validation: tests/validation/test_king_physics.py::TestKingEquilibriumVelocityDF::test_virial_ratio_is_half_unscaled
validation: tests/validation/test_king_physics.py::TestKingEquilibriumVelocityDF::test_dispersion_profile_matches_king_moment
validation: tests/validation/test_king_physics.py::TestKingVelocityDF::test_velocities_bound_against_king_escape_speed

Michie–King anisotropic DF¶

✅ verified · 2 equations · 2 sources · API: MichieVelocityDF

The self-consistent radially anisotropic King model: Michie’s (1963) angular-momentum Gaussian times the King lowering. Unlike the Osipkov–Merritt overlay, the DF is specified first and Poisson is re-solved, so the density itself changes with r_a.

Use it for	Not for
✓ globular clusters with observed outer radial anisotropy	✗ holding a FIXED density while adding anisotropy (use the OM overlay)
✓ self-consistent anisotropic equilibria (pairs with MichieProfile)	✗ r_a below the radial-orbit bound (construction raises: infinite-mass tail)

Parameters¶

Name	Meaning	Units	Typical range	Code
`W0`	dimensionless central potential (as King)	dimensionless	3–12	`MichieVelocityDF(W0=...)`
`r_a`	anisotropy radius of the J^2 Gaussian; r_a -> inf recovers King	pc	>~ 3-4 r_c at W0=7 (below: radial-orbit pathology, raises)	`MichieVelocityDF(r_a=...)`

Equations¶

f(E, J) \propto \exp\!\left(-\frac{J^2}{2 r_a^2 \sigma^2}\right) \left[\exp\!\left(-\frac{E}{\sigma^2}\right) - 1\right],\qquad E \le 0

(7)

Symbols: $J$ : specific angular momentum r v_t [pc^2 / Myr]; $r_a$ : anisotropy radius [pc]. Assumes: single-mass; Michie 1963 J^2 Gaussian x King 1966 lowering. Derivation: theory page.

\frac{1}{\xi^2}\frac{d}{d\xi}\!\left(\xi^2\frac{d\psi}{d\xi}\right) = -9\,\frac{\hat\rho(\psi,\ \xi/\hat r_a)}{\hat\rho(W_0,0)}, \qquad \psi(0)=W_0,\ \psi^\prime(0)=0

(8)

Symbols: $\xi$ : dimensionless radius r/r_c [dimensionless]; $\psi$ : dimensionless relative potential [dimensionless]. Assumes: King-1966 nondimensionalization of Michie 1963 Eq. 5.8 (anisotropy weight e^{-C z^2 eta^2 (1-mu^2)} = e^{-J^2/2 r_a^2 sigma^2}). Derivation: theory page.

Sources¶

Michie (1963) — the anisotropic DF ansatz and the radius-dependent Poisson equation (Eq. 4.0 (f = A e^{-m/m_0(alpha E + beta J^2)} Q); Eq. 5.0 (assembled DF); Eq. 5.8 (Poisson with anisotropy weight, p. 138 — scanned source, visually verified 2026-07-11))
King (1966) — the lowered-Maxwellian cutoff and the factor-9 nondimensionalization (Eq. 3; Eqs. 15-16)

Code & validation¶

code: MichieVelocityDF (src/progenax/kinematics/michie_df.py)
validation: tests/validation/test_michie_physics.py::TestMichieEquilibrium::test_virial_ratio_half_unscaled
validation: tests/validation/test_michie_physics.py::TestMichieAnisotropyProfile::test_beta_increases_outward
validation: tests/validation/test_michie_physics.py::TestMichieAnisotropyProfile::test_beta_matches_df_oracle

EFF Eddington-inversion DF¶

✅ verified · 1 equation · 2 sources · API: EFFVelocityDF

The isotropic ergodic DF of the (truncated) EFF density, obtained by numerical Eddington inversion — tabulated f(E) with the truncation boundary term, sampled by inverse-CDF.

Use it for	Not for
✓ equilibrium velocities for EFF positions (Q_vir ~ 0.5 at gamma=5, mild truncation)	✗ strong truncation r_t <~ few a at steep gamma (a few percent sub-virial — see the EFF chapter)
✓ young-cluster power-law halos where Plummer is too steep	✗ self-consistent anisotropy (OM overlay available via anisotropy_radius)

Parameters¶

Name	Meaning	Units	Typical range	Code
`a`	EFF core scale radius (must match the profile’s)	pc	0.1–10	`EFFVelocityDF(a=...)`
`gamma`	3-D density slope (must match the profile’s)	dimensionless	3.2–7	`EFFVelocityDF(gamma=...)`
`r_t`	truncation radius (sets the Psi = 0 zero point)	pc	5–100	`EFFVelocityDF(r_t=...)`
`anisotropy_radius`	optional Osipkov-Merritt r_a (None = isotropic)	pc	large enough that f(Q) >= 0 (construction-checked)	`EFFVelocityDF(..., anisotropy_radius=...)`

Equations¶

f(\mathcal{E}) \;=\; \frac{1}{\sqrt{8}\,\pi^2} \left[\int_0^{\mathcal{E}} \frac{\mathrm{d}^2\rho}{\mathrm{d}\Psi^2}\frac{\mathrm{d}\Psi}{\sqrt{\mathcal{E}-\Psi}} + \frac{1}{\sqrt{\mathcal{E}}}\left(\frac{\mathrm{d}\rho}{\mathrm{d}\Psi}\right)_{\!\Psi=0}\right]

(9)

Symbols: $\mathcal{E}$ : binding energy relative to the truncation zero point [pc^2 / Myr^2]. Assumes: isotropic; truncated: the boundary term (drho/dPsi)_0 is NONZERO and load-bearing. Derivation: theory page.

Sources¶

Elson et al. (1987) — the density being inverted (Eq. 1 (projected law; 3-D convention per the deviations note on the EFF profile card))
Merritt (1985) — the Osipkov-Merritt augmented-density construction for the anisotropic option (Eqs. 45-46 (per-paper note, PDF-verified))

Code & validation¶

code: EFFVelocityDF (src/progenax/kinematics/eff_df.py)
validation: tests/validation/test_eff_physics.py::TestEFFVelocityDF::test_eff_eddington_virial_ratio_mild_truncation
validation: tests/validation/test_eff_physics.py::TestEFFVelocityDF::test_eff_eddington_f_is_physical
validation: tests/validation/test_eff_physics.py::TestEFFVelocityDF::test_eff_all_particles_bound

Osipkov–Merritt anisotropy overlay¶

✅ verified · 2 equations · 2 sources · API: PlummerVelocityDF

The one-parameter radial-anisotropy construction shared by the Plummer and EFF DFs: replace f(E) by f(Q) built from the augmented density rho_Q = (1 + r^2/r_a^2) rho, holding the TRUE density and potential fixed. beta(r) = r^2/(r^2 + r_a^2) exactly.

Use it for	Not for
✓ adding outer radial anisotropy to a FIXED density (Plummer/EFF)	✗ King-family anisotropy (KingVelocityDF stays isotropic — use Michie, whose density responds to r_a)
✓ differentiable r_a inference (grad-audited registry cases)	✗ r_a below the non-negativity bound (Plummer: r_a >= 0.75 a; construction raises)

Parameters¶

Name	Meaning	Units	Typical range	Code
`anisotropy_radius`	r_a: isotropic inside, radial outside; r_a -> inf recovers isotropy	pc	Plummer: >= 0.75 a (Merritt 1985 Eq. 46); EFF: f(Q) >= 0 checked at construction	`PlummerVelocityDF/EFFVelocityDF(anisotropy_radius=...)`

Equations¶

Q \;\equiv\; \mathcal{E} - \frac{L^2}{2\,r_a^2}

(10)

Symbols: $L$ : specific angular momentum |r x v| [pc^2 / Myr]. Assumes: f = f(Q): the OM augmented integral of motion. Derivation: theory page.

\beta(r) \;\equiv\; 1 - \frac{\sigma_t^2(r)}{\sigma_r^2(r)} \;=\; \frac{r^2}{r^2 + r_a^2}

(11)

Symbols: $\beta$ : Binney anisotropy (single-component sigma_t convention; the BT two-component form gives the SAME beta) [dimensionless]. Assumes: exact for any OM f(Q). Derivation: theory page.

Sources¶

Merritt (1985) — the OM construction, the closed-form anisotropic Plummer DF, and its non-negativity bound (Eq. 45 (anisotropic Plummer f(Q)); Eq. 46 (r_a >= 3/4 r_0 bound))
Dejonghe (1987) — independent cross-check: his q-family beta = (q/2) r^2/(1+r^2); the q=2 member equals the Merritt r_a=1 model (Eq. 21 + following discussion, p. 18; Eq. 23 (beta(r)))

Code & validation¶

code: PlummerVelocityDF (src/progenax/kinematics/plummer_df.py)
code: EFFVelocityDF (src/progenax/kinematics/eff_df.py)
validation: tests/validation/test_dispersion_physics.py::test_eff_om_beta_identity
validation: tests/validation/test_engine_b_physics.py::test_om_beta_profile_realized

References¶

Plummer, H. C. (1911). On the problem of distribution in globular star clusters. Monthly Notices of the Royal Astronomical Society, 71, 460–470. 10.1093/mnras/71.5.460
Dejonghe, H. (1987). A completely analytical family of anisotropic Plummer models. Monthly Notices of the Royal Astronomical Society, 224, 13–39. 10.1093/mnras/224.1.13
Aarseth, S. J., Henon, M., & Wielen, R. (1974). A comparison of numerical methods for the study of star cluster dynamics. Astronomy and Astrophysics, 37, 183–187.
King, I. R. (1966). The structure of star clusters. III. Some simple dynamical models. The Astronomical Journal, 71, 64–75. 10.1086/109857
Michie, R. W. (1963). On the distribution of high energy stars in spherical stellar systems. Monthly Notices of the Royal Astronomical Society, 125, 127–139. 10.1093/mnras/125.2.127
Elson, R. A. W., Fall, S. M., & Freeman, K. C. (1987). The structure of young star clusters in the Large Magellanic Cloud. The Astrophysical Journal, 323, 54–78. 10.1086/165807
Merritt, D. (1985). Spherical stellar systems with spheroidal velocity distributions. The Astronomical Journal, 90, 1027–1037. 10.1086/113810