Plummer Model

The Plummer model (or Plummer sphere) is an analytical, spherically symmetric gravitational potential–density pair introduced by British astronomer Henry Crozier Plummer in 1911 to describe the radial star-count distribution of globular clusters. Formulated as an astrophysical realization of Arthur Schuster's 1883 solution to the Lane–Emden equation for a self-gravitating polytropic gas sphere of index n = 5, the Plummer model expresses both the Newtonian gravitational potential Φ(r) and the three-dimensional mass density ρ(r) in simple algebraic closed forms without transcendental functions. Because every dynamical observable—including the enclosed mass profile M(r), sky-projected surface density Σ(R), isotropic velocity dispersion σ²(r) from the Jeans equations, and phase-space distribution function f(ℰ) via the Eddington inversion formula—can be written analytically in elementary powers, the Plummer sphere serves as a foundational benchmark in theoretical stellar dynamics and as the standard gravitational softening kernel in N-body simulations.

A Plummer sphere of total mass M and characteristic Plummer scale radius a (also denoted r_P or b) is defined by the softened point-mass gravitational potential:

Φ(r) = −GM / (r² + a²)^1/2

where G is the gravitational constant and r is the three-dimensional radial distance from the cluster center. At large distances (r ≫ a), Φ(r) converges asymptotically to the Keplerian point-mass potential −GM/r, whereas near the origin (r ≪ a) the potential expands quadratically as a harmonic oscillator well, Φ(r) ≈ −(GM/a)[1 − r²/(2a²)].

Applying the spherically symmetric Poisson equation, ∇²Φ(r) = (1/r²) d/dr(r² dΦ/dr) = 4πGρ(r), yields the three-dimensional Plummer mass density profile:

ρ(r) = (3M / 4πa³) (1 + r²/a²)^−5/2

Near the center (r → 0), the density approaches a constant harmonic core value ρ₀ = 3M/(4πa³) with zero radial gradient (dρ/dr = 0), contrasting with the divergent inner cusps of the Hernquist profile (ρ ∝ r⁻¹), Jaffe profile (ρ ∝ r⁻²), and Dehnen profile family. At large radii (r ≫ a), the Plummer density falls off steeply as ρ(r) ∝ r⁻⁵, ensuring that the cumulative enclosed mass M(r) converges rapidly to the finite total mass M without requiring an outer tidal truncation:

M(r) = ∫₀^r 4πr′² ρ(r′) dr′ = Mr³ / (r² + a²)^3/2

Setting M(r_h) = ½M determines the three-dimensional half-mass radius r_h = (2^2/3 − 1)^−1/2a ≈ 1.305*a*, while enclosing 90% of the total mass requires r₉₀ ≈ 3.70*a*.

When a spherical Plummer cluster is observed on the celestial sphere at projected cylindrical radius R = (x² + y²)^1/2 from the center, integrating ρ(r) along the line of sight z = (r² − R²)^1/2 via the Abel projection integral gives a simple rational law for the projected surface mass density Σ(R):

Σ(R) = 2 ∫_R^∞ ρ(r) r dr / (r² − R²)^1/2 = (M / πa²) (1 + R²/a²)⁻²

For a stellar population with spatially uniform mass-to-light ratio Υ = M/L, the surface brightness profile I(R) = Σ(R)/Υ has central surface brightness I(0) = L/(πa²). Integrating Σ(R) over circular apertures on the sky yields the cumulative projected mass M_proj(R) = MR²/(R² + a²). Evaluating M_proj(R_e) = ½M shows that the sky-projected effective (half-light) radius R_e is identical to the Plummer scale parameter (R_e = a), while the photometric core radius where Σ(R_c) = ½Σ(0) is R_c = (2^1/2 − 1)^1/2a ≈ 0.6436*a*.

In 1916, Arthur Stanley Eddington showed that the Plummer model is the unique self-consistent spherical stellar system whose ergodic phase-space distribution function f(ℰ)—governing the density of stars in six-dimensional phase space obeying the collisionless Boltzmann equation—is a pure power law of the relative binding energy ℰ = Ψ(r) − ½v², where Ψ(r) = −Φ(r) = GM/(r² + a²)^1/2 > 0. Because ρ(r) ∝ [(r² + a²)^−1/2]⁵, the three-dimensional space density is a fifth-power function of the relative potential:

ρ(Ψ) = (3M / 4πa³) (a / GM)⁵ Ψ⁵

In the theory of stellar polytropes (ρ ∝ Ψⁿ), this relation identifies the Plummer sphere as the index n = 5 polytrope, which is the limiting polytropic index (n ≤ 5) that possesses a finite total mass while extending to infinite radius. Applying the Eddington inversion formula to ρ(Ψ) ∝ Ψ⁵ yields the closed-form isotropic phase-space distribution function:

f(ℰ) = [24(2)^1/2 / 7π³] (a² / G⁵M⁴) ℰ^7/2 (for ℰ > 0)

Solving the isotropic spherical Jeans equations, d(ρσ_r²)/dr = −ρ dΦ/dr, gives the local one-dimensional radial velocity dispersion σ_r²(r) = Ψ(r)/6 = GM/[6(r² + a²)^1/2]. Integrating the translational kinetic energy of stellar motion over all space yields K = (3π/64)(GM²/a) and total gravitational potential energy W = −2K = −(3π/32)(GM²/a), satisfying the scalar virial theorem 2K + W = 0.

Although Plummer (1911) fitted the inner regions of several Milky Way globular clusters with his n = 5 profile, real star clusters orbiting inside a host galaxy are tidally stripped at a finite Jacobi boundary r_t, favoring the tidally truncated King profile that captures mass segregation and stellar evaporation. By contrast, dwarf spheroidal galaxies and ultra-faint Milky Way satellites are routinely modeled with two-dimensional Plummer surface-density profiles to measure their half-light radii R_e = a.

In computational astrophysics, collisional star-cluster codes initialize test ensembles in standard Hénon N-body units (G = M = 1, a = 3π/16) by sampling positions and velocities analytically from the Plummer model, while flattening the potential along the z-axis generates the Miyamoto–Nagai potential for galactic disks. Furthermore, collisionless galaxy and dark matter simulations replace the singular Newtonian pairwise potential −Gm_im_j/r_ij with the Plummer-softened potential −Gm_im_j/(r_ij² + ε²)^1/2 to suppress artificial two-body deflections below the softening length ε during violent relaxation and phase mixing.

Why does the Plummer model correspond to an n = 5 polytrope?

A spherical polytrope of index n has a three-dimensional mass density proportional to the n-th power of the relative gravitational potential, ρ(Ψ) ∝ Ψⁿ. Because the Plummer potential scales as Ψ(r) ∝ (r² + a²)^−1/2 and its density scales as ρ(r) ∝ (r² + a²)^−5/2, substituting Ψ(r) gives ρ(Ψ) ∝ Ψ⁵, satisfying the Lane–Emden equation for n = 5.

How does the Plummer scale radius relate to the half-light and half-mass radii?

When projected onto the two-dimensional plane of the sky, the cylindrical radius enclosing half of the total cluster luminosity (R_e) equals the Plummer scale parameter (R_e = a). In three-dimensional space, the spherical radius enclosing half of the total mass is r_h = (2^2/3 − 1)^−1/2a ≈ 1.305*a*.

How do you sample particle positions and velocities for an N-body Plummer sphere?

Draw a uniform random number X₁ ∈ (0, 1) for the cumulative mass fraction M(r)/M and invert the mass profile to obtain r = a(X₁^−2/3 − 1)^−1/2. Then draw the normalized speed q = v/v_esc(r) ∈ [0, 1) from g(q) ∝ q²(1 − q²)^7/2.

1. Plummer, H. C. (1911). "On the problem of distribution in globular star clusters." Monthly Notices of the Royal Astronomical Society, 71(5), 460–470. https://doi.org/10.1093/mnras/71.5.460
2. Eddington, A. S. (1916). "The distribution of stars in globular clusters." Monthly Notices of the Royal Astronomical Society, 76(7), 572–585. https://doi.org/10.1093/mnras/76.7.572
3. Dejonghe, H. (1987). "Anisotropic Plummer models." Monthly Notices of the Royal Astronomical Society, 224(1), 13–39. https://doi.org/10.1093/mnras/224.1.13
4. Hénon, M. (1971). "The Monte Carlo Method." Astrophysics and Space Science, 14(1), 151–167. https://doi.org/10.1007/BF00649201
5. Binney, J., & Tremaine, S. (2008). Galactic Dynamics (2nd ed.). Princeton University Press. https://doi.org/10.1515/9781400828722