Differentially Rotating Polytropes

• Differentially rotating polytropes are equilibrium configurations characterized by variable angular velocity profiles and a polytropic equation of state that governs their structure.
• They employ various rotation laws, including the KEH and generalized power-law profiles, which influence mass-shedding limits and stability thresholds.
• Numerical methods such as multidomain pseudospectral codes and self-consistent field iterations deliver high-accuracy models essential for simulating neutron star mergers and core-collapse events.

Differentially rotating polytropes are equilibrium configurations of self-gravitating, axisymmetric fluids following a polytropic equation of state (EOS), in which the angular velocity varies with cylindrical radius. These models are fundamental in the theory of neutron stars and related astrophysical compact objects, where differential rotation is naturally generated by processes such as stellar core collapse and binary neutron star merger. The profound structural and stability differences between differentially rotating and uniformly rotating polytropes underlie their distinctive mass-shedding limits, equilibrium morphologies, and dynamical instability thresholds.

1. Polytropic Equation of State and Thermodynamics

The polytropic EOS is typically employed in theoretical and numerical studies of rotating stars due to its analytic tractability and capacity to capture the stiffness characteristic of nuclear matter:

• The canonical form is $P = K\rho^{\Gamma}$, where $P$ is the pressure, $\rho$ the rest-mass (baryonic) density, $K$ the polytropic constant, and $\Gamma=1+1/n$ with $n$ the polytropic index.
• For neutron star studies, $n=1$ ($\Gamma=2$) is commonly adopted, leading to convenient thermodynamic properties such as a linear correspondence between enthalpy and density:

$H(\rho) = 2K\rho \quad (n=1)$

and total (energy) density $\epsilon = \rho + N P$.

• The dimensionless enthalpy variable

$H(p) = \int_0^{p} \frac{dp'}{\epsilon(p') + p'}$

provides a computationally efficient parametrization for equilibrium sequences.

This theoretical scaffolding allows exact solution of the Bernoulli integral and Poisson's equation for $n=1$, facilitating both Newtonian and relativistic treatments (Gondek-Rosinska et al., 2016, Razinkova et al., 11 Jan 2026, Iosif et al., 2020, Szewczyk et al., 2023).

2. Laws of Differential Rotation

The angular velocity profile—how rapidly the local azimuthal velocity falls with cylindrical radius—is critical in determining the star's structure and stability. Multiple rotation laws of increasing generality have been explored:

• Komatsu–Eriguchi–Hachisu (KEH) law: The classic "j-constant" law is

$\Omega(\varpi) = \frac{\Omega_c}{1 + (\varpi/A)^2}$

with $\Omega_c$ the central angular velocity and $A$ the differential rotation lengthscale. The degree of differential rotation is characterized by the dimensionless parameter $\widetilde{A} = r_e/A$, where $r_e$ is the equatorial radius. $\widetilde{A} \to 0$ yields uniform rotation (Gondek-Rosinska et al., 2016, Szewczyk et al., 2023).

• Generalized Power-Law Profiles: The parameterization introduced by Galeazzi et al. extends the possible outer power-law slopes via:

$\Omega(R) = \Omega_c [1 + (R/R_0)^2]^{1/\alpha}$

This form can reproduce "j-constant" ($\alpha = -1$), "v-constant" ($\alpha = -2$), Kepler-like ($\alpha = -4/3$), and merger-remnant-like ($\alpha \approx -4$, i.e., $\Omega \propto R^{-1/2}$) profiles (Galeazzi et al., 2011).

• Uryu+ Law: A four-parameter law,

$$\Omega(F) = \Omega_c \frac{1 + \left( \frac{F}{B^2\Omega_c} \right)^p}{1 + \left( \frac{F}{A^2\Omega_c} \right)^{p+q}}$$

with $F = u^t u_\phi$, generalizing both monotonic and nonmonotonic profiles with two distinct length scales and independent control of core and envelope shear. This formulation enables accurate modeling of post-merger neutron star remnants’ rotation profiles (Iosif et al., 2020).

3. Numerical Construction of Equilibria

The computation of differentially rotating polytrope equilibria involves solving nonlinear elliptic PDEs for the gravitational potential and metric, coupled to the hydrostationary first integral. Key approaches include:

• Multidomain Pseudospectral Codes: The AKM method solves the Einstein–Euler system as four elliptic equations for metric potentials, using Chebyshev expansions and Newton–Raphson iteration with spectral resolutions up to $36\times36$ collocation points, achieving accuracies of $10^{-8} - 10^{-12}$ for global quantities (Gondek-Rosinska et al., 2016).
• Self-Consistent Field (SCF) Methods: Post-Newtonian SCF iterations treat both Newtonian and $1$ PN corrections:
  • Quantities are expanded in two control parameters: rotational strength ($\bar{\upsilon}$) and relativistic compactness ($\bar{\sigma}$). Either rotation law and boundary conditions at the stellar surface (mass-shedding) determine the sequence termination (Fotopoulos et al., 2023).
• Grid-Based Newtonian Approaches: The "ROTAT" code solves the linearized Bernoulli–Poisson system at $n=1$ on a uniform grid, with the Rayleigh stability criterion directly imposed via the rotation law (Razinkova et al., 11 Jan 2026).

Convergence is typically constrained by virial error metrics or loss of regularity at the mass-shedding (Keplerian) boundary.

4. Classification of Equilibrium Sequences and Morphologies

The taxonomy of rotating polytrope equilibrium sequences is driven by both the degree and law of differential rotation:

• General Relativistic Sequences (Ansorg–Gondek-Rosińska–Villain):
  • Four types are identified for fixed maximum enthalpy $H_m$ and degree of differential rotation $\widetilde{A}$:
  • Type A: Spheroidal, connected to the static limit, terminating at mass-shedding.
  • Type C: Toroidal, starting at the static star but developing a central density minimum and never reaching mass-shedding before becoming toroidal.
  • Type B: Highly deformed spheroids, no nonrotating limit; start at mass-shedding, end at toroidal onset, and can reach extreme masses.
  • Type D: Double-shed, connecting two mass-shedding points; appears in a narrow $\widetilde{A}$ window (Gondek-Rosinska et al., 2016, Iosif et al., 2020).
• Newtonian $n=1$ Branches (Razinkova–Yudin–Blinnikov):
  • Type T: Spheroid–to–thick-torus sequence, present for all differentialities.
  • Type P: "Hockey-puck" shapes, at moderate differentiality ($\sigma \sim 2$) in narrow flattening.
  • Type M: "Matryoshka" (spheroid + torus) solutions, with core-torus detachment at high rotational support, again only for $\sigma \sim 2$.
  • For both too weak and too strong differential rotation ($\sigma \lesssim 1.7$ or $\sigma \gtrsim 2.2$), only the classical spheroid–torus sequence exists (Razinkova et al., 11 Jan 2026).

The sequence type determines the possible maximum mass and structure. Table 1 summarizes the main morphological sequences:

Branch Type	Structure	Shear Regime
Type A/T	Spheroids	Weak/Moderate
Type C	Toroidal	Strong
Type B	Deformed Spheroids	Moderate
Type D	Double-shedding	Narrow interval
Type P	Hockey-puck	Moderate ($\sim 2$)
Type M	Spheroid+Torus	Moderate ($\sim 2$)

5. Maximum Masses and Stability Thresholds

Differential rotation significantly extends the mass limit above that of non-rotating (TOV) and rigidly rotating polytropes:

• Maximum Mass Scaling:
  • For uniform rotation, the limit rises to about $1.2 M_{\rm TOV}$.
  • Modest degrees of differential rotation (Type A) reach $\sim 1.6 M_{\rm TOV}$.
  • Type C (toroidal, strong shear) achieves up to $1.5 M_{\rm TOV}$, but at extreme toroidality.
  • Type B sequences, newly explored in (Gondek-Rosinska et al., 2016), reach up to $4 M_{\rm TOV}$ for $\Gamma=2$, at intermediate $\widetilde{A} \sim 0.7-0.8$.
  • For $\widetilde{A}=1$ (typical in merger remnants), the allowed gravitational mass is $\sim 1.9 M_{\rm TOV}$ for quasi-toroidal solutions (Szewczyk et al., 2023).
• Stability Criteria:
  • The onset of secular instability for uniform rotation is at $T/|W| \sim 0.14$; the dynamical bar-mode threshold is $T/|W| \sim 0.27$.
  • Strongly differentially rotating models, particularly of Type M in Newtonian gravity, can reach $\tau \gtrsim 0.36$ (rotational-to-binding energy fraction), well into the regime of expected dynamical fragmentation (Razinkova et al., 11 Jan 2026).
  • The turning-point method along constant rest-mass sequences provides a reliable criterion for the end of stability, corroborated by time-dependent simulations (Szewczyk et al., 2023).

Configurations with extreme differential rotation and high mass generally violate the Kerr bound $J/M^2 < 1$ (producing "supra-Kerrian" models), but secular angular momentum redistribution typically restores this inequality over dynamical timescales (Gondek-Rosinska et al., 2016).

6. Physical and Astrophysical Implications

The study of differentially rotating polytropes directly informs models of neutron star formation and merger remnants:

• Post-merger and post-core-collapse neutron stars are expected to be born highly differentially rotating, at least for $\sim 10-100\,\mathrm{ms}$.
• Subclasses with modest shear (Type B or quasi-toroidal Type C for $\widetilde{A} \sim 1$) can support "hypermassive" remnants above the widely cited $1.5 M_{\rm TOV}$ threshold, with maximum support reaching $\sim 4 M_{\rm TOV}$ (Gondek-Rosinska et al., 2016); for more merger-like rotation profiles, $1.9 M_{\rm TOV}$ is attainable (Szewczyk et al., 2023).
• These configurations have profound consequences for the lifetime of post-merger objects, the delay to black hole collapse, the emitted gravitational wave spectrum (high-frequency, rapidly damped oscillations), and the energetics of associated kilonovae and gamma-ray bursts.

Empirically, equilibrium models built using merger-motivated rotation profiles (Uryu+ law, $\lambda_1 \sim 1.7-2.0$, $\lambda_2 \sim 1.0$) reproduce observed rotational structure in simulations and serve as robust initial data for evolutionary studies (Iosif et al., 2020).

7. Sensitivity to Rotation Law and Numerical Schemes

The choice of differential rotation law has a subtle yet consequential influence:

• The total gravitational mass is remarkably insensitive ($<0.1\%$) to the exact rotation law parameters for given global rotation, although the equatorial radius and detailed density structure vary by up to $10\%$ (Iosif et al., 2020).
• Newtonian and post-Newtonian calculations are reliable for compactness parameter $\bar{\sigma} \lesssim 0.03$ for $n=1$ polytropes, but require full general relativity at higher compactness (Fotopoulos et al., 2023).
• In merger-motivated profiles ($\alpha \approx -4$), the secular bar-mode threshold $T/|W| = 0.14$ is never reached before mass-shedding, implying these states should be secularly stable; however, "low-$T/|W|$" or shear-driven dynamical instabilities may still arise (Galeazzi et al., 2011, Iosif et al., 2020).

The combined advances in rotation law generality, numerical accuracy, and classification of morphologies have greatly expanded the known solution space of differentially rotating polytropes, driving new insights into the structure and fate of compact objects in extreme astrophysical environments.