Setting the velocity to zero in the Euler equation: the scale height, the polytropic atmosphere, the Lane–Emden equation, and a tabulated solar model in which one polytrope is confirmed, one is refuted, and the boundary between them disagrees with helioseismology.
1. The problem: a star that should fall in half an hour
Module 2 derived the Euler equation, its equation (7.3), and ended with a table which showed that viscosity may be dropped in the interiors of stars. This module applies that equation in the simplest case it has: a gas that does not move.
The simplest case is also the most consequential one, and one number shows why. Take the Sun's mass $M_\odot$ and radius $R_\odot$, remove its pressure, and let it fall. A pressureless uniform sphere of mean density $\bar\rho$ collapses to its centre in the free-fall time
derived in Problem D3. Here $G$ is the gravitational constant and $\bar\rho = 3M/(4\pi R^3)$. For the Sun $\bar\rho = 1.4081\ \mathrm{g\,cm^{-3}}$, so $t_{\rm ff} = 1770.3\ \mathrm{s} = 29.5$ minutes.
The Sun does not collapse in 29.5 minutes. Its age of $4.6$ Gyr is $8.2\times10^{13}$ free-fall times. Something holds it up, and that something is the pressure gradient. The question this module answers is how precisely the pressure gradient must balance gravity, and what structure that balance forces on a star.
A second timescale sharpens the question. The Kelvin–Helmholtz time
$$ t_{\rm KH} \;=\; \frac{GM^2}{RL} \tag{1.2}$$
is the time in which a star of luminosity $L$ would radiate its gravitational binding energy, to order of magnitude. For the Sun it is $9.9089\times10^{14}$ s, or $31.4$ Myr. If the Sun had no nuclear source it would contract on this timescale. It is $5.5973\times10^{11}$ times longer than the free-fall time. A star can change, but it changes on a timescale so much longer than the dynamical one that at each instant it is in balance.
The module follows the three-part shape of this book. §3 and §4 derive the terrestrial result: the scale height of the Earth's atmosphere and its lapse rate as a polytrope. §5 to §8 change the regime to a self-gravitating sphere and derive the Lane–Emden equation. §9 and §10 check three predictions against a tabulated standard solar model and against a helioseismic measurement. One prediction is confirmed, one is refuted, and the third finds a 29.4σ disagreement between the model and the measured Sun.
in which $\rho$ is the mass density in $\mathrm{g\,cm^{-3}}$, $\mathbf{u}$ is the bulk velocity, $P$ is the scalar pressure in $\mathrm{dyn\,cm^{-2}}$ (Module 2 wrote it $p$), and $\mathbf{g}$ is the gravitational acceleration. This module writes the pressure as a capital $P$ because lower-case $p$ is not needed and $P$ is the stellar-structure convention.
Definition 1 (hydrostatic equilibrium). A fluid is in hydrostatic equilibrium when $\mathbf{u}=\mathbf{0}$ everywhere and at all times considered. Then (2.1) reduces to
$$ \nabla P \;=\; \rho\,\mathbf{g}. \tag{2.2}$$
Equation (2.2) says that the net force on each fluid element is zero: the pressure pushes up with the same force per unit volume with which gravity pulls down. For a spherically symmetric body, the gravitational acceleration at radius $r$ depends only on the mass inside that radius.
Proposition 1 (spherical hydrostatic equilibrium). In a spherically symmetric, static, self-gravitating fluid with density $\rho(r)$, let
$$ m(r) \;\equiv\; \int_0^r 4\pi r'^2\rho(r')\,dr' \tag{2.3}$$
be the mass inside radius $r$. Then
$$ \frac{dP}{dr} \;=\; -\,\frac{G\,m(r)\,\rho(r)}{r^2}, \qquad \frac{dm}{dr} \;=\; 4\pi r^2\rho. \tag{2.4}$$
The second equation of (2.4) is the derivative of the definition (2.3). For the first, apply Gauss's law for gravity, $\nabla\cdot\mathbf{g} = -4\pi G\rho$, to a sphere of radius $r$. Spherical symmetry makes $\mathbf{g} = g_r(r)\,\hat{\mathbf{r}}$, so the divergence theorem gives $4\pi r^2 g_r = -4\pi G\,m(r)$, that is $g_r = -Gm(r)/r^2$. The mass outside $r$ contributes nothing, by the same theorem. Spherical symmetry also makes $\nabla P = (dP/dr)\,\hat{\mathbf{r}}$. Substituting both into (2.2) and taking the radial component gives the first equation of (2.4). ∎
Equation (2.4) is two equations in three unknown functions, $P$, $\rho$ and $m$. It is not closed. Some relation between $P$ and $\rho$ must be supplied, and the choice of that relation is the whole subject of §5. First we ask whether setting $\mathbf{u}=\mathbf{0}$ is justified at all.
2.1 How precisely the balance holds
A star that contracts is not static. The question is how large an imbalance between the two terms on the right of (2.1) a contraction requires. The answer is a ratio of timescales, and it is small enough to justify (2.2) by a very wide margin.
Proposition 2 (the imbalance required by a slow contraction). Let a star of mass $M$ and radius $R$ contract so that its surface accelerates at $|\ddot R| = R/t^2$, which defines the contraction timescale $t$. Then the net force per unit mass at the surface, as a fraction of the surface gravity $GM/R^2$, is
$$ \frac{|\ddot R|}{GM/R^2} \;=\; \frac{8}{\pi^2}\left(\frac{t_{\rm ff}}{t}\right)^{2}, \tag{2.5}$$
with $t_{\rm ff}$ given by (1.1).
The radial component of (2.1) for the surface element, which moves with the surface, is $\ddot R = -(1/\rho)\,dP/dr - GM/R^2$. The sum of the two force terms on the right is therefore $\ddot R$, and its size relative to the gravity term is $|\ddot R|/(GM/R^2) = (R/t^2)(R^2/GM) = R^3/(GMt^2)$. Now express $R^3/(GM)$ through the free-fall time. With $\bar\rho = 3M/(4\pi R^3)$, equation (1.1) gives
$$ t_{\rm ff}^2 \;=\; \frac{3\pi}{32\,G}\cdot\frac{4\pi R^3}{3M} \;=\; \frac{\pi^2}{8}\,\frac{R^3}{GM}, $$
so $R^3/(GM) = (8/\pi^2)\,t_{\rm ff}^2$. Substituting gives (2.5). ∎
Apply (2.5) to the fastest contraction the Sun could sustain from its own luminosity, $t = t_{\rm KH}$. The coefficient-free factor is $(t_{\rm ff}/t_{\rm KH})^2 = 3.192\times10^{-24}$; the exact coefficient $8/\pi^2 = 0.8106$ multiplies it.
Result. If the Sun contracted on its Kelvin–Helmholtz time of $31.4$ Myr, gravity and the pressure gradient would still have to balance to $2.587\times10^{-24}$ of the surface gravity, by (2.5). On the nuclear timescale of $1.036\times10^{10}$ yr, on which the Sun actually evolves, the fraction is $2.379\times10^{-29}$.
Note the direction of this argument. The small imbalance is not assumed in order to simplify the equations. It is forced by the observed luminosity. By (2.5), an imbalance larger than about $10^{-24}$ means contraction on a timescale $t$ shorter than $t_{\rm KH}$. Such a contraction releases gravitational energy of order $GM^2/R$ in time $t$, a power $GM^2/(Rt)$ larger than $L$ by the factor $t_{\rm KH}/t$ from (1.2). By the virial argument of §8, half of that energy must leave the star, and the observed luminosity does not carry it. Setting $\mathbf{u}=\mathbf{0}$ in (2.1) discards a term of that relative size.
The nuclear timescale uses $0.007\times0.1\,Mc^2/L$: $0.7$ per cent of rest mass released by hydrogen burning, applied to the inner $10$ per cent of the mass. Both factors are conventional order-of-magnitude inputs, not results of this book. The solar $GM$, $R$ and $L$ are the IAU 2015 nominal values. The mean density uses the radius printed in the tabulated solar model of §9, $6.9598\times10^{10}$ cm, which differs from the IAU nominal $6.957\times10^{10}$ cm by $0.04$ per cent; the Kelvin–Helmholtz time uses the IAU value.
A sound wave carries the information that restores balance. The isothermal sound speed at the Sun's centre, $\sqrt{P_c/\rho_c}$, is $394.1\ \mathrm{km\,s^{-1}}$, and the sound crossing time $R/c_s$ is $1765.1$ s, or $29.4$ minutes. It differs from $t_{\rm ff}$ by $0.3$ per cent. Both are order-of-magnitude estimates, so the closeness to $0.3$ per cent is partly accidental. That the two are of the same order is not: the virial theorem, Proposition 12 of §8, gives $3\int P\,dV = -\Omega$, and with $\Omega\sim -GM^2/R$ and $\int P\,dV\sim (P/\rho)M$ this makes $P/\rho\sim GM/R$, so $R/c_s\sim\sqrt{R^3/GM}\sim t_{\rm ff}$. Any departure from balance is corrected in about half an hour. Module 4 studies those sound waves directly.
3. The plane-parallel atmosphere and the scale height
Before the whole star, take a thin layer at its surface, or the atmosphere of a planet. If the layer is thin compared with the radius, the curvature of the surface and the change of $g$ with height can both be ignored.
Definition 2 (plane-parallel layer). A plane-parallel layer is a region in which all quantities depend only on a height coordinate $z$, measured upward, and the gravitational acceleration is uniform, $\mathbf{g} = -g\,\hat{\mathbf{z}}$ with $g\gt 0$ constant. Equation (2.2) then becomes
$$ \frac{dP}{dz} \;=\; -\rho\, g. \tag{3.1}$$
Equation (3.1) needs a second relation between $P$ and $\rho$. The ideal gas law supplies one when the temperature is known:
$$ P \;=\; \frac{\rho\, k_B T}{\mu\, m_u}, \tag{3.2}$$
in which $k_B$ is Boltzmann's constant, $T$ the temperature, $m_u$ the atomic mass unit, and $\mu$ the mean molecular weight, the mean mass per free particle in units of $m_u$. For dry air $\mu = 28.9647$; for a fully ionised gas of hydrogen mass fraction $X$, helium mass fraction $Y$ and heavy-element fraction $Z$, $\mu = 1/(2X + \tfrac{3}{4}Y + \tfrac{1}{2}Z)$.
Proposition 3 (isothermal scale height). In a plane-parallel layer of ideal gas at uniform temperature $T$ and uniform $\mu$,
$$ P(z) \;=\; P(0)\,e^{-z/H}, \qquad H \;\equiv\; \frac{k_B T}{\mu\, m_u\, g}. \tag{3.3}$$
The density falls with the same exponential, $\rho(z) = \rho(0)\,e^{-z/H}$. The length $H$ is the pressure scale height.
Solve (3.2) for the density, $\rho = \mu m_u P/(k_BT)$, and substitute into (3.1):
$$ \frac{dP}{dz} \;=\; -\frac{\mu\, m_u\, g}{k_B T}\,P \;=\; -\frac{P}{H}. $$
Because $T$, $\mu$ and $g$ are constant, $H$ is constant, and this linear equation has the solution $P(0)e^{-z/H}$. The density is $\mu m_u P/(k_BT)$ with a constant prefactor, so it has the same dependence on $z$. ∎
The scale height is the height over which the pressure falls by a factor $e$. It is the height through which an energy $k_BT$, the thermal energy per particle to within a factor of order one, would lift one particle against its weight $\mu m_u g$. A hot gas or a light gas or a weak gravity gives a thick atmosphere.
Proposition 3 was proved under the assumption that $g$ is uniform. For a layer on a body of radius $r$ and mass $M$, $g = GM/r^2$, and it changes by a fractional amount of about $2H/r$ across one scale height. The layer's area also changes by the same fraction. So the plane-parallel form is valid when
$$ H \;\ll\; r. \tag{3.4}$$
The table evaluates $H$ in four regimes, with $g = GM/r^2$ in the last two.
Read the last column. On the Earth, $H$ is $1/755$ of the planetary radius. At the solar photosphere it is $1/5167$ of the solar radius, so the solar atmosphere is thinner, relative to its star, than the Earth's atmosphere is relative to the Earth, and plane-parallel models of it are safe. The intracluster medium is the opposite case. Its scale height is $0.64$ of the radius at which it is evaluated, so $g$ changes by a large factor within one scale height and (3.3) does not describe it. Such a gas must be treated with the spherical equation (2.4), in which the scale height is not a constant.
The photospheric temperature is the IAU 2015 nominal solar effective temperature, $5772$ K. The mean molecular weight $\mu=1.3$ is that of a neutral hydrogen–helium mixture and is an input of moderate trust: the photosphere is weakly ionised, and a different $\mu$ changes $H$ in inverse proportion. The molecular-cloud and cluster masses and radii are representative round values chosen to show the regime, not measurements of any named object; $\mu = 2.33$ is the usual value for molecular hydrogen with helium. The Earth values are the sea-level constants of the International Standard Atmosphere, applied here as if the whole atmosphere were at $288.15$ K, which §4 corrects.
4. The polytropic atmosphere: the terrestrial half
The Earth's lower atmosphere is not isothermal. Its temperature falls with height at a nearly constant rate, the lapse rate. This section shows that a constant lapse rate is exactly what a power-law relation between pressure and density produces, and so it introduces the closure that the rest of the module applies to stars.
Definition 3 (polytrope). A fluid is polytropic with polytropic index $n$ and polytropic constant $K$ if, throughout the region considered,
$$ P \;=\; K\,\rho^{1+1/n}, \tag{4.1}$$
with $K$ and $n$ constant, $0\le n\le\infty$. The case $n\to\infty$ is the isothermal gas; $n=0$ is a gas of uniform density.
Definition 3 is a statement about the structure of a body, a relation between $P$ and $\rho$ at different places. It need not be a statement about how one parcel responds when it is compressed, although in §5 it will be one.
Proposition 4 (lapse rate of a polytropic atmosphere). In a plane-parallel layer of ideal gas (3.2) with uniform $\mu$ and $g$, obeying the polytropic relation (4.1), the temperature falls linearly with height at the rate
$$ \left|\frac{dT}{dz}\right| \;=\; \frac{\mu\, m_u\, g}{k_B\,(n+1)}, \tag{4.2}$$
and the pressure and temperature are related by
$$ \frac{P}{P(0)} \;=\; \left(\frac{T}{T(0)}\right)^{n+1}. \tag{4.3}$$
From (4.1), $\rho = (P/K)^{n/(n+1)}$. Substituting into the ideal gas law (3.2) and solving for $T$,
$$ T \;=\; \frac{\mu m_u}{k_B}\,\frac{P}{\rho} \;=\; \frac{\mu m_u}{k_B}\,K^{n/(n+1)}\,P^{1/(n+1)}. $$
The prefactor is constant, so $T\propto P^{1/(n+1)}$, which is (4.3). Differentiate the logarithm of this relation with respect to $z$:
$$ \frac{1}{T}\frac{dT}{dz} \;=\; \frac{1}{n+1}\,\frac{1}{P}\frac{dP}{dz} \;=\; -\frac{1}{n+1}\,\frac{\rho\, g}{P}, $$
using (3.1) in the second step. By (3.2), $\rho/P = \mu m_u/(k_BT)$, so the right-hand side is $-\mu m_u g/(k_B T(n+1))$. Multiplying by $T$ gives $dT/dz = -\mu m_u g/(k_B(n+1))$, a constant, which is (4.2). ∎
Inverting (4.2) turns a measured lapse rate $\Gamma \equiv |dT/dz|$ into a polytropic index, $n = \mu m_u g/(k_B\Gamma) - 1$. Two lapse rates of the Earth's atmosphere are worth converting.
The first is the lapse rate of the International Standard Atmosphere (ISA), $6.5\ \mathrm{K\,km^{-1}}$ up to the tropopause at $11$ km. It gives $n = 4.2558$.
The second is the dry adiabatic lapse rate, the rate at which a parcel of dry air cools when it rises without exchanging heat. A parcel of ideal gas compressed or expanded adiabatically obeys $P\propto\rho^{\gamma}$, with $\gamma$ the ratio of specific heats. If the atmosphere is well mixed, so that every level lies on the same adiabat, the structure obeys $P=K\rho^{\gamma}$ with one $K$, and comparison with (4.1) gives $1+1/n = \gamma$, or
$$ n \;=\; \frac{1}{\gamma-1}. \tag{4.4}$$
For the diatomic molecules of dry air, $\gamma = 7/5$ and $n = 2.5$. Substituting into (4.2) gives $\Gamma = g(\gamma-1)\mu m_u/(\gamma k_B) = g/c_p$, where $c_p = \gamma k_B/((\gamma-1)\mu m_u) = 1004.7\ \mathrm{J\,kg^{-1}\,K^{-1}}$ is the specific heat at constant pressure. Numerically $g/c_p = 9.761\ \mathrm{K\,km^{-1}}$, and inverting (4.2) returns $n = 2.5000$.
Three hydrostatic atmospheres with the same sea-level state, $T(0)=288.15$ K and $P(0)=1013.25$ hPa. Left: temperature against altitude. Right: pressure, from (3.3) for the isothermal case and (4.3) for the two polytropes. The standard atmosphere, $n=4.26$, lies between the isothermal limit and the dry adiabat, $n=2.50$; at $11$ km the dry adiabat is $36$ K colder. The circle marks the tabulated tropopause pressure of $226$ hPa, which the $n=4.26$ curve reaches by construction.
Fig. 1 compares the three. The real troposphere cools more slowly than the dry adiabat, by $35.9$ K over $11$ km. The difference is the latent heat released when water vapour condenses in rising air, which warms the parcel as it rises. The index $4.26$ is therefore not an adiabatic index of any gas. It is an effective index that summarises the combined effect of adiabatic cooling and condensation heating on the structure.
Equation (4.3) also gives the pressure at the tropopause. With $n = 4.2558$, $T(0) = 288.15$ K and a tropopause temperature of $288.15 - 6.5\times11 = 216.65$ K, it gives $226.32$ hPa, in ratio $1.00001$ to the tabulated $226.32$ hPa.
This agreement is a check of the algebra, not a test of nature. The ISA is defined by its sea-level temperature and pressure, its lapse rate and the hydrostatic equation, and its tabulated tropopause values are computed from those constants. Agreement to $1.00001$ confirms that (4.2) and (4.3) are derived and evaluated correctly. It says nothing about whether the real atmosphere is polytropic. The tests against nature in this module are in §10.
The ISA constants are from ISO 2533:1975, identical to the US Standard Atmosphere 1976 below 11 km: $T(0)=288.15$ K, $P(0)=101325$ Pa, lapse rate $6.5\ \mathrm{K\,km^{-1}}$, tropopause at $11$ km geopotential altitude with $216.65$ K and $22632$ Pa. They are definitional and are trusted exactly. The attribution of the gap between $6.5$ and $9.761\ \mathrm{K\,km^{-1}}$ to latent heat is standard meteorology and is not computed here.
5. Why a polytrope: two indices fixed in advance
The spherical equations (2.4) need a closure, and Definition 3 is the simplest one. A closure is only useful if its index can be fixed before the star is examined. If $n$ were chosen to fit a model, agreement with that model would prove nothing. This section derives two indices from physics alone.
Proposition 5 (adiabatic convection gives $n = 3/2$). Suppose a region of an ideal monatomic gas has uniform specific entropy, so that every part of it lies on one adiabat $P = K\rho^{\gamma}$ with $\gamma = 5/3$. Then the region is a polytrope of index $n = 3/2$.
For an ideal gas whose particles have only the three translational degrees of freedom, the internal energy per particle is $\tfrac{3}{2}k_BT$ (Module 2, Proposition 1), so the specific heats per particle are $c_V = \tfrac{3}{2}k_B$ and $c_P = c_V + k_B = \tfrac{5}{2}k_B$, and $\gamma = c_P/c_V = 5/3$. Along an adiabat, $P\rho^{-\gamma}$ is constant; uniform specific entropy means the same constant $K$ applies everywhere in the region. So $P = K\rho^{5/3}$ throughout, and matching the exponent to (4.1), $1 + 1/n = 5/3$, gives $n = 3/2$. ∎
The physical input is the hypothesis of uniform entropy. Convection supplies it: rising and sinking parcels mix the gas, and where convection is efficient it carries the stratification to the adiabatic state. That is asserted here and derived in Module 6 through the Schwarzschild criterion. The index follows from the hypothesis with no free parameter. The monatomic assumption requires the gas to be fully ionised or neutral atomic; in a zone where hydrogen or helium is partly ionised, ionisation absorbs energy and $\gamma$ falls below $5/3$.
The second index comes from radiation pressure. A hot star contains a photon gas as well as a particle gas, and both push.
Definition 4 (gas-pressure fraction). In a gas of ideal particles and blackbody radiation at temperature $T$, the total pressure is $P = P_{\rm gas} + P_{\rm rad}$ with
$$ P_{\rm gas} \;=\; \frac{k_B}{\mu m_u}\,\rho T, \qquad P_{\rm rad} \;=\; \frac{a_{\rm rad} T^4}{3}, \tag{5.1}$$
where $a_{\rm rad} = 7.5657\times10^{-15}\ \mathrm{erg\,cm^{-3}\,K^{-4}}$ is the radiation constant. The gas-pressure fraction is $\beta \equiv P_{\rm gas}/P$, so $1-\beta = P_{\rm rad}/P$.
Proposition 6 (Eddington's standard model gives $n = 3$). If $\beta$ and $\mu$ are the same everywhere in a star, with $0\lt\beta\lt 1$, then
$$ P \;=\; K\rho^{4/3}, \qquad K \;=\; \left[\frac{3}{a_{\rm rad}}\,\frac{1-\beta}{\beta^4}\right]^{1/3}\left(\frac{k_B}{\mu m_u}\right)^{4/3}, \tag{5.2}$$
so the star is a polytrope of index $n = 3$.
Eliminate $T$ between the two parts of (5.1). From $P_{\rm rad} = (1-\beta)P$, the temperature is $T = [3(1-\beta)P/a_{\rm rad}]^{1/4}$. From $P_{\rm gas} = \beta P$,
$$ \beta P \;=\; \frac{k_B}{\mu m_u}\,\rho\,\left[\frac{3(1-\beta)P}{a_{\rm rad}}\right]^{1/4}. $$
Divide both sides by $\beta P^{1/4}$:
$$ P^{3/4} \;=\; \frac{k_B}{\mu m_u}\,\frac{1}{\beta}\left[\frac{3(1-\beta)}{a_{\rm rad}}\right]^{1/4}\rho. $$
Raise both sides to the power $4/3$:
$$ P \;=\; \left(\frac{k_B}{\mu m_u}\right)^{4/3}\beta^{-4/3}\left[\frac{3(1-\beta)}{a_{\rm rad}}\right]^{1/3}\rho^{4/3}, $$
which is (5.2) after collecting $\beta^{-4/3} = (\beta^{-4})^{1/3}$ inside the bracket. If $\beta$ and $\mu$ are constant, $K$ is constant and $P\propto\rho^{4/3}$; by (4.1), $1+1/n = 4/3$, so $n=3$. ∎
Note what the proof needs. It needs $\beta$ to be constant, and it needs $1-\beta$ to be non-zero, since otherwise $K$ vanishes and the relation says nothing. It does not need $1-\beta$ to be large. The hypothesis of constant $\beta$ has a direct test. Since $1-\beta = a_{\rm rad}T^4/(3P)$, constant $\beta$ means $T^4\propto P$, that is $d\ln T/d\ln P = 1/4$ at every depth. §10 measures that slope through the Sun, and measures the size of $1-\beta$ at its centre.
Both indices are now fixed. The number $3/2$ follows from uniform entropy in a monatomic gas, and the number $3$ follows from a uniform radiation share. Neither was read from a solar model, and so each can be tested against one.
6. The Lane–Emden equation
Substituting the polytropic closure (4.1) into the spherical equations (2.4) produces one second-order equation for the density. Scaling out the dimensional constants leaves an equation that depends on $n$ alone.
Proposition 7 (the Lane–Emden equation). Let a spherical star in hydrostatic equilibrium (2.4) obey the polytropic relation (4.1) with index $n\gt 0$, and let $\rho_c$ be its central density. Define the dimensionless density variable $\theta$ and radius $\xi$ by
$$ \rho \;=\; \rho_c\,\theta^{n}, \qquad r \;=\; a\,\xi, \qquad a^2 \;\equiv\; \frac{(n+1)\,K\,\rho_c^{1/n-1}}{4\pi G}. \tag{6.1}$$
Then $\theta(\xi)$ satisfies
$$ \frac{1}{\xi^2}\,\frac{d}{d\xi}\left(\xi^2\,\frac{d\theta}{d\xi}\right) \;=\; -\,\theta^{n}, \tag{6.2}$$
with the boundary conditions
$$ \theta(0) \;=\; 1, \qquad \theta'(0) \;=\; 0. \tag{6.3}$$
The surface of the star is the first zero of $\theta$, at $\xi = \xi_1$, where $R = a\,\xi_1$.
First combine the two equations of (2.4) into one. Multiply the first by $r^2/\rho$, so that $r^2\rho^{-1}\,dP/dr = -Gm$, and differentiate with respect to $r$. The second equation replaces $dm/dr$:
$$ \frac{1}{r^2}\,\frac{d}{dr}\left(\frac{r^2}{\rho}\,\frac{dP}{dr}\right) \;=\; -4\pi G\rho. $$
Next write the pressure in terms of $\theta$. With $\rho = \rho_c\theta^n$, the closure (4.1) gives $P = K\rho_c^{1+1/n}\theta^{n+1}$, so
$$ \frac{dP}{dr} \;=\; (n+1)\,K\rho_c^{1+1/n}\,\theta^{n}\,\frac{d\theta}{dr}, \qquad \frac{1}{\rho}\frac{dP}{dr} \;=\; (n+1)\,K\rho_c^{1/n}\,\frac{d\theta}{dr}. $$
Substituting and dividing by $4\pi G\rho_c$,
$$ \frac{(n+1)\,K\rho_c^{1/n-1}}{4\pi G}\;\frac{1}{r^2}\frac{d}{dr}\left(r^2\frac{d\theta}{dr}\right) \;=\; -\theta^n. $$
The prefactor is $a^2$ by (6.1), and it has dimensions of length squared. Change variable to $\xi = r/a$. Then $d/dr = a^{-1}\,d/d\xi$ and $r^2 = a^2\xi^2$, so $r^{-2}\,d(r^2\,d\theta/dr)/dr = a^{-2}\,\xi^{-2}\,d(\xi^2\,d\theta/d\xi)/d\xi$. The factor $a^{-2}$ cancels the prefactor $a^2$, which gives (6.2).
For the boundary conditions: $\theta(0)=1$ because $\rho(0)=\rho_c$ by definition. For the second, note that $m(r)\approx\tfrac{4}{3}\pi r^3\rho_c$ near the centre, so $dP/dr = -Gm\rho/r^2\approx -\tfrac{4}{3}\pi G\rho_c^2\, r \to 0$ as $r\to 0$. Since $dP/dr$ is proportional to $\theta^n\,d\theta/dr$ and $\theta(0)=1$, this requires $\theta'(0)=0$. The pressure, and therefore the density, first vanish where $\theta$ first vanishes, which is the surface. ∎
The scaling has removed $K$, $\rho_c$ and $G$. Every polytrope of index $n$, whatever its mass and radius, is the same function $\theta_n(\xi)$, stretched in radius by $a$ and in density by $\rho_c$. So one integration for each $n$ describes all stars of that index.
The point $\xi=0$ is a singular point of (6.2), because the term $2\theta'/\xi$ has a zero denominator. A numerical integration cannot start there, and it starts instead from a short series.
Proposition 8 (series near the centre). The solution of (6.2) with (6.3) has the expansion
$$ \theta(\xi) \;=\; 1 \;-\; \frac{\xi^2}{6} \;+\; \frac{n\,\xi^4}{120} \;-\; \cdots. \tag{6.4}$$
Expanding (6.2) gives $\theta'' + 2\theta'/\xi = -\theta^n$. The equation is unchanged under $\xi\to-\xi$, and the conditions (6.3) are even, so the regular solution is an even function. Write $\theta = 1 + b_2\xi^2 + b_4\xi^4 + O(\xi^6)$. Then $\theta' = 2b_2\xi + 4b_4\xi^3$ and $\theta'' = 2b_2 + 12b_4\xi^2$, so
$$ \theta'' + \frac{2\theta'}{\xi} \;=\; 6b_2 + 20b_4\,\xi^2 + O(\xi^4). $$
On the right, by the binomial expansion, $\theta^n = 1 + n\,b_2\,\xi^2 + O(\xi^4)$. Matching the constant terms gives $6b_2 = -1$, so $b_2 = -1/6$. Matching the $\xi^2$ terms gives $20b_4 = -n b_2 = n/6$, so $b_4 = n/120$. ∎
The case $n=0$ needs one remark before it can be used. With $n=0$, the substitution $\rho=\rho_c\theta^n$ makes $\rho$ uniform and says nothing about $\theta$, and $\rho_c^{1/n}$ in (6.1) is undefined. Since $P_c = K\rho_c^{1+1/n}$, the length in (6.1) can be written $a^2 = (n+1)P_c/(4\pi G\rho_c^2)$, and the closure gives $P = P_c\,\theta^{n+1}$. Both forms remain meaningful at $n=0$, and we take them as the definitions of $a$ and $\theta$ there. For uniform density, (2.4) gives $dP/dr = -\tfrac{4}{3}\pi G\rho_c^2 r$, so $P = P_c - \tfrac{2}{3}\pi G\rho_c^2 r^2$; with $\theta = P/P_c$ and $r = a\xi$, $a^2 = P_c/(4\pi G\rho_c^2)$, this is $\theta = 1-\xi^2/6$, which is the solution of (6.2) with $\theta^0 = 1$. So (6.2) and (6.3) hold at $n=0$ with this reading, and (7.2) below applies to it.
\n\n
Three values of $n$ have solutions in closed form. They are the test cases for the numerical integrator, and each shows a different kind of star.
Proposition 9 (exact solutions for $n=0$, $1$, $5$). The solutions of (6.2) and (6.3) for these three indices are
$$ \theta_0 \;=\; 1 - \frac{\xi^2}{6}, \qquad \theta_1 \;=\; \frac{\sin\xi}{\xi}, \qquad \theta_5 \;=\; \left(1+\frac{\xi^2}{3}\right)^{-1/2}. \tag{6.5}$$
Their surfaces are at $\xi_1 = \sqrt6$ with $\theta_0'(\xi_1) = -\sqrt6/3$, and at $\xi_1 = \pi$ with $\theta_1'(\xi_1) = -1/\pi$. The $n=5$ solution has no zero, so its radius is infinite, but $-\xi^2\theta_5'\to\sqrt3$ as $\xi\to\infty$, so its mass is finite.
Case $n=0$. Then $\theta^0=1$ and (6.2) is $(\xi^2\theta')' = -\xi^2$. Integrate once: $\xi^2\theta' = -\xi^3/3 + C$. At $\xi=0$ the left side is zero, so $C=0$ and $\theta' = -\xi/3$. Integrate again with $\theta(0)=1$ to get $\theta_0 = 1-\xi^2/6$. It vanishes at $\xi_1=\sqrt6$, where $\theta' = -\sqrt6/3$. The density $\rho_c\theta^0$ is uniform, but the pressure $\propto\theta^{n+1} = \theta$ is not.
Case $n=1$. Substitute $\theta = \chi/\xi$. Then $\xi^2\theta' = \xi\chi' - \chi$, and its derivative is $\xi\chi''$, so (6.2) becomes $\chi''/\xi = -\chi/\xi$, that is $\chi'' = -\chi$. The general solution is $\chi = A\sin\xi + B\cos\xi$. The term $B\cos\xi/\xi$ diverges at the origin, so $B=0$; and $\sin\xi/\xi\to1$ as $\xi\to0$, so $A=1$. The first zero of $\sin\xi$ after the origin is $\xi_1=\pi$. There $\theta_1' = (\xi\cos\xi - \sin\xi)/\xi^2 = -\pi/\pi^2 = -1/\pi$.
Case $n=5$. Let $s = 1+\xi^2/3$. Differentiating $\theta_5 = s^{-1/2}$ gives $\theta_5' = -(\xi/3)\,s^{-3/2}$, so $\xi^2\theta_5' = -(\xi^3/3)\,s^{-3/2}$. Differentiating again,
$$ (\xi^2\theta_5')' \;=\; -\xi^2 s^{-3/2} + \frac{\xi^3}{3}\cdot\frac{3}{2}\cdot\frac{2\xi}{3}\,s^{-5/2} \;=\; -\xi^2 s^{-5/2}\left(s - \frac{\xi^2}{3}\right) \;=\; -\xi^2 s^{-5/2}. $$
Dividing by $\xi^2$ gives $-s^{-5/2} = -\theta_5^5$, which is (6.2). The conditions (6.3) hold by inspection. $\theta_5\gt0$ for every $\xi$, so there is no surface. For large $\xi$, $s^{-3/2}\to(\xi^2/3)^{-3/2} = 3^{3/2}\xi^{-3}$, so $-\xi^2\theta_5' \to (\xi^3/3)\cdot3^{3/2}\xi^{-3} = \sqrt3$. ∎
The integrator in figs/m03_numbers.py is a fourth-order Runge–Kutta scheme started from (6.4) at $\xi=10^{-5}$. On the three exact cases it returns $\xi_1 = 2.449490$ for $n=0$ and $3.141593$ for $n=1$, both in ratio $1.000000$ to the exact values. It returns $\theta_1'(\xi_1) = -0.318310$ for $n=1$, in ratio $1.000000$. For $n=5$ it returns $\theta(200) = 8.659914\times10^{-3}$, in ratio $1.000000$ to (6.5), and $-\xi^2\theta' = 1.731856$ at $\xi=200$, in ratio $0.999888$ to $\sqrt3$. The largest departure is on $\theta_0'(\xi_1) = -0.816392$, in ratio $0.999871$ to the exact $-\sqrt6/3$: the $n=0$ solution meets the surface at a finite slope, and the last step lands just past it. The integrator is therefore trusted to about $1.3\times10^{-4}$ on the cases $n=3/2$ and $n=3$, which have no closed form.
Upper panel: solutions of the Lane–Emden equation (6.2) for $n = 0$, $1$, $3/2$, $3$, $4$ and $5$. The dot on each curve is its surface $\xi_1$; the $n=5$ solution has none. A higher index gives a more extended profile with a larger $\xi_1$. Lower panel: the density profiles $\rho/\rho_c = \theta^n$ for $n=3/2$ and $n=3$, each scaled to its own radius, against the tabulated solar model of §9 on a logarithmic axis. Neither polytrope follows the Sun over the whole radius; the Sun's central condensation $\rho_c/\bar\rho$ is $106.9$, against $54.2$ for $n=3$.
7. Properties of polytropes
The dimensionless solution fixes every structural ratio of a polytrope. This section extracts the ones that can be compared with a star of known mass and radius.
Proposition 10 (mass, radius, central condensation and central pressure). A polytrope of index $n\lt 5$ with first zero $\xi_1$ and surface slope $\theta_1' \equiv \theta'(\xi_1)$ has
$$ R \;=\; a\,\xi_1, \qquad M \;=\; 4\pi a^3\rho_c\left(-\xi_1^2\theta_1'\right), \tag{7.1}$$
$$ D_n \;\equiv\; \frac{\rho_c}{\bar\rho} \;=\; \frac{\xi_1}{3\,|\theta_1'|}, \qquad W_n \;\equiv\; \frac{P_c\,R^4}{G M^2} \;=\; \frac{1}{4\pi(n+1)\,\theta_1'^{\,2}}, \tag{7.2}$$
and, for $n\neq3$, the mass–radius relation
$$ R \;\propto\; K^{\,n/(3-n)}\,M^{(1-n)/(3-n)}. \tag{7.3}$$
The radius follows from $r=a\xi$ at the surface. For the mass, integrate (2.3) with $\rho=\rho_c\theta^n$ and $r=a\xi$:
$$ M \;=\; 4\pi a^3\rho_c\int_0^{\xi_1}\xi^2\theta^n\,d\xi \;=\; 4\pi a^3\rho_c\int_0^{\xi_1}\left(-\frac{d}{d\xi}\big(\xi^2\theta'\big)\right)d\xi \;=\; 4\pi a^3\rho_c\left(-\xi_1^2\theta_1'\right), $$
using (6.2) in the middle step and $\xi^2\theta'=0$ at $\xi=0$.
For $D_n$: $\bar\rho = 3M/(4\pi R^3) = 3\rho_c(-\xi_1^2\theta_1')/\xi_1^3 = 3\rho_c|\theta_1'|/\xi_1$, and dividing $\rho_c$ by this gives the first of (7.2).
For $W_n$: the central pressure is $P_c = K\rho_c^{1+1/n}$. By (6.1), $K\rho_c^{1/n} = 4\pi G a^2\rho_c/(n+1)$, so $P_c = 4\pi G a^2\rho_c^2/(n+1)$. From (7.1), $\rho_c = M/(4\pi a^3\mu_1)$ with $\mu_1\equiv-\xi_1^2\theta_1'$. Substituting,
$$ P_c \;=\; \frac{4\pi G a^2}{n+1}\cdot\frac{M^2}{16\pi^2a^6\mu_1^2} \;=\; \frac{GM^2}{4\pi(n+1)\,a^4\mu_1^2}. $$
With $a = R/\xi_1$, $a^4\mu_1^2 = R^4\xi_1^{-4}\cdot\xi_1^4\theta_1'^2 = R^4\theta_1'^2$, which gives the second of (7.2).
For (7.3): from (6.1), $a^2\propto K\rho_c^{(1-n)/n}$, and from (7.1), $\rho_c\propto M/a^3$. Substituting, $a^2\propto K M^{(1-n)/n}a^{-3(1-n)/n}$, so $a^{(3-n)/n}\propto K M^{(1-n)/n}$. Raising to the power $n/(3-n)$ and using $R\propto a$ gives (7.3). ∎
The two ratios $D_n$ and $W_n$ are the useful ones. Given only the mass and radius of a star, and an index $n$, they predict its central density and central pressure with no further input. §10 uses them in exactly that way. The integrator gives the following values.
$n$
$\xi_1$
$-\xi_1^2\theta_1'$
$D_n=\rho_c/\bar\rho$
$W_n=P_cR^4/GM^2$
$0$
$2.44949$
$4.89835$
$1.0001$
$0.11940$
$1$
$3.14159$
$3.14159$
$3.2899$
$0.39270$
$3/2$
$3.65375$
$2.71406$
$5.9907$
$0.77014$
$3$
$6.89685$
$2.01824$
$54.1825$
$11.05068$
$4$
$14.97155$
$1.79723$
$622.4079$
$247.55944$
The central condensation $D_n$ rises steeply with $n$, from $1$ for the uniform sphere to $54$ for $n=3$ and $622$ for $n=4$. A high index describes a star with a dense core and an extended, tenuous envelope. The $n=0$ entry $1.0001$ departs from the exact value $1$ by the $1.3\times10^{-4}$ error in $\theta_0'(\xi_1)$ noted in §6.
Check against a published value. Chandrasekhar (1939), An Introduction to the Study of Stellar Structure, chapter IV, Table 4, gives for $n=3$: $\xi_1 = 6.89685$, $-\xi_1^2\theta_1' = 2.01824$ and $\rho_c/\bar\rho = 54.1825$. The integrator gives the same six significant figures for each; the ratios of the unrounded values are $1.000000$, $0.999998$ and $1.000000$.
The mass–radius relation (7.3) has two special cases. For $n=1$ the exponent of $M$ is zero, so all $n=1$ polytropes with the same $K$ have the same radius. For $n=3$ the exponent is undefined, and (7.1) shows instead that $M\propto K^{3/2}$ with no dependence on $\rho_c$: the mass is fixed by $K$ alone. Problem D1 uses both facts, and Problem K3 uses the second.
8. What balance gives with no closure
Everything in §6 and §7 depends on the closure (4.1). Before testing a closure, it is worth asking what hydrostatic equilibrium implies with no closure at all. The answer is two exact results. One is an inequality, and one is an identity.
Proposition 11 (Chandrasekhar's lower bound on the central pressure). In a spherical star of mass $M$ and radius $R$ in hydrostatic equilibrium (2.4), with finite central density and non-negative surface pressure,
$$ P_c \;\ge\; \frac{G M^2}{8\pi R^4}. \tag{8.1}$$
Define $F(r) \equiv P(r) + Gm(r)^2/(8\pi r^4)$. Differentiate, using both equations of (2.4):
$$ \frac{dF}{dr} \;=\; -\frac{Gm\rho}{r^2} \;+\; \frac{2Gm}{8\pi r^4}\cdot4\pi r^2\rho \;-\; \frac{4Gm^2}{8\pi r^5} \;=\; -\frac{Gm^2}{2\pi r^5} \;\le\; 0. $$
The first two terms cancel exactly. So $F$ does not increase outward, and $F(0)\ge F(R)$. At the centre, $m\approx\tfrac{4}{3}\pi r^3\rho_c$ makes $m^2/r^4\propto r^2\to0$, so $F(0)=P_c$. At the surface, $F(R) = P(R) + GM^2/(8\pi R^4) \ge GM^2/(8\pi R^4)$. Combining gives (8.1). ∎
For the Sun, with the tabulated radius of §9, the bound is $4.4750\times10^{14}\ \mathrm{dyn\,cm^{-2}}$. The tabulated central pressure is $2.3380\times10^{17}\ \mathrm{dyn\,cm^{-2}}$, which exceeds the bound by a factor $522.5$. The bound is true and it excludes almost nothing. The factor of $522$ is the information that a closure must supply. In the notation of (7.2), the bound is $W\ge1/(8\pi) = 0.03979$; the tabulated Sun has $W = 20.7881$, and $n=3$ gives $W_3 = 11.0507$.
Proposition 12 (the virial theorem for a static star). Under the hypotheses of Proposition 11, with zero surface pressure $P(R)=0$,
$$ 3\int_0^R P\,dV \;=\; -\,\Omega, \qquad \Omega \;\equiv\; -\int_0^M \frac{G\,m}{r}\,dm, \tag{8.2}$$
where $dV = 4\pi r^2 dr$ and $\Omega$ is the gravitational potential energy of the star.
Multiply the first equation of (2.4) by $4\pi r^3$ and integrate from $0$ to $R$. On the left, integrate by parts:
$$ \int_0^R 4\pi r^3\frac{dP}{dr}\,dr \;=\; \Big[4\pi r^3 P\Big]_0^R - \int_0^R 12\pi r^2 P\,dr \;=\; -3\int_0^R P\,dV, $$
since the boundary term vanishes at $r=0$ and, because $P(R)=0$, at $r=R$. On the right,
$$ -\int_0^R 4\pi r^3\,\frac{Gm\rho}{r^2}\,dr \;=\; -\int_0^R \frac{Gm}{r}\,4\pi r^2\rho\,dr \;=\; -\int_0^M\frac{Gm}{r}\,dm \;=\; \Omega, $$
using $dm = 4\pi r^2\rho\,dr$. Equating the two sides gives $-3\int P\,dV = \Omega$, which is (8.2). ∎
For an ideal monatomic gas, $\tfrac{3}{2}P$ is the internal energy per unit volume (Module 2, Proposition 1), so (8.2) says that the thermal energy $U = \int\tfrac{3}{2}P\,dV$ is $-\Omega/2$. The total energy is then $E = U + \Omega = \Omega/2$. A star that radiates loses energy, so $E$ and $\Omega$ decrease; a more negative $\Omega$ means a more compact star, and $U = -\Omega/2$ increases. The star contracts and its thermal energy rises, and half the released gravitational energy is radiated. That is why $t_{\rm KH}$ in (1.2), the binding energy divided by the luminosity, is the timescale of contraction.
Integrating the tabulated solar model of §9 gives $3\int P\,dV = 6.1594\times10^{48}$ erg and $-\Omega = 6.1564\times10^{48}$ erg, in ratio $1.00049$.
This agreement does not confirm hydrostatic equilibrium. The solar model was built by solving (2.4) numerically. Any correct reading of its table must satisfy (8.2), so a ratio of $1.00049$ confirms that the table was read and integrated correctly, and nothing more. The residual $4.9\times10^{-4}$ has not been decomposed; two candidate sources are the trapezoid rule and the table ending at $0.983\,R$ rather than at $P=0$. Module 2's solar-wind check did not have this circularity, because Venzmer & Bothmer fitted the density and speed exponents with no continuity constraint imposed. Here the tests must be of a closure, given the equation, and §10 is built that way.
9. The tabulated solar model
The comparison object for the checks is a standard solar model: a numerical solution of the equations of stellar structure, evolved from a uniform initial composition to the present solar age and adjusted to reproduce the present solar radius and luminosity. The one used here is BS2005-AGS,OP of Bahcall, Serenelli & Basu (2005), ApJ621, L85–L88. Its table is stored as data/bs05_agsop.dat.
The table has $1284$ rows, from $r/R = 0.00161$ to $0.98308$. Each row gives the enclosed mass fraction, the radius fraction, $T$, $\rho$, $P$, the enclosed luminosity fraction, and six mass fractions. Its first row gives the central values used below: $T_c = 1.5480\times10^{7}$ K, $\rho_c = 150.50\ \mathrm{g\,cm^{-3}}$ and $P_c = 2.3380\times10^{17}\ \mathrm{dyn\,cm^{-2}}$. The surface helium mass fraction in the last row is $0.22905$, against $0.229$ printed in Table 1 of the paper, which confirms that the file is the model the paper describes.
Two quantities computed from the table are needed in §10. The first is an effective mean molecular weight. Inverting the ideal gas law (3.2) at the centre, $\mu = \rho k_B T/(P m_u)$, gives $\mu = 0.82851$. The composition in the same row, $X = 0.36462$, $Y=0.62029$, $Z = 0.01509$, gives $\mu = 0.83195$ for a fully ionised gas. The ratio is $0.99587$. The tabulated $\mu$ is an effective value which absorbs every contribution to $P$ beyond the ideal particle gas; radiation, the next quantity, accounts for $6\times10^{-4}$ of the $4.1\times10^{-3}$ difference, and the rest is not attributed here.
The second quantity is the radiation share. At the centre, $P_{\rm rad} = a_{\rm rad}T_c^4/3 = 1.4482\times10^{14}\ \mathrm{dyn\,cm^{-2}}$, so
Radiation supplies $6.19\times10^{-4}$ of the pressure at the centre of the Sun. Hold that number. It is the premise of Proposition 6 evaluated in the Sun.
The model table is trusted as a faithful copy of the published BS2005-AGS,OP model: its surface helium matches the paper, and the virial identity (8.2) holds to $4.9\times10^{-4}$. It is not trusted as a description of the real Sun; §10 shows where it fails. Its radii and masses are fractions, converted here with the table's own constants $R_\odot = 6.9598\times10^{10}$ cm and $L_\odot = 3.8418\times10^{33}\ \mathrm{erg\,s^{-1}}$, and with $M_\odot = GM_\odot/G$ from the IAU 2015 nominal $GM_\odot$. Its "central" values are those of the first row, at $r = 0.00161\,R$.
10. The checks
The circularity noted in §8 decides what can be tested. The solar model satisfies (2.4) by construction, so the table cannot test hydrostatic equilibrium. It can test a closure: given (2.4), is the relation between $P$ and $\rho$ in the model the one that §5 derived? Each of the two indices of §5 was fixed before the table was read. This section tests each, then takes the boundary between them to a helioseismic measurement.
10.1 Check 1: Eddington's $n=3$ for the whole Sun
Proposition 10 predicts the central density and pressure of an $n=3$ star from its mass and radius alone: $\rho_c = D_3\,\bar\rho$ and $P_c = W_3\,GM^2/R^4$. With the tabulated radius, $\bar\rho = 1.4081\ \mathrm{g\,cm^{-3}}$.
Quantity
$n=3$ prediction
Tabulated model
Ratio
$\rho_c$
$76.293\ \mathrm{g\,cm^{-3}}$
$150.500\ \mathrm{g\,cm^{-3}}$
$0.5069$
$P_c$
$1.2429\times10^{17}\ \mathrm{dyn\,cm^{-2}}$
$2.3380\times10^{17}\ \mathrm{dyn\,cm^{-2}}$
$0.5316$
$T_c$, from (3.2) with $\mu=0.82851$
$1.6233\times10^{7}$ K
$1.5480\times10^{7}$ K
$1.0486$
Check 1: refuted. Eddington's $n=3$ predicts a central condensation $\rho_c/\bar\rho = 54.18$. The tabulated model has $106.88$. The ratio is $0.5069$, and the central pressure ratio is $0.5316$. Both are wrong by about a factor of two. The central temperature is wrong by only $4.9$ per cent, in ratio $1.0486$.
The temperature survives because, by (3.2) at fixed $\mu$, $T\propto P/\rho$. The temperature ratio is therefore the pressure ratio divided by the density ratio, $0.5316/0.5069 = 1.0486$. Both errors have the same sign and nearly the same size, so their quotient is close to one. A central temperature from $n=3$ is close to the model's value even though the structure that produced it is not.
The reason for the failure is the premise of Proposition 6. That proposition requires a constant radiation share $1-\beta$, which §5 showed is the same as $d\ln T/d\ln P = 1/4$ at every depth. In the rows tabulated in §10.2 this slope runs from $0.1937$ at $0.4494\,R$ to $0.3963$ at $0.7500\,R$, so $\beta$ is not constant and the premise fails. The size of the radiation share explains why nothing holds $\beta$ fixed in the Sun. By (9.1) radiation carries $6.19\times10^{-4}$ of the central pressure, so removing radiation pressure would change the central pressure by that fraction, and the structure is set by the gas pressure together with the transport of energy. Problem K3 solves Eddington's own relation for $\beta$ and finds $1-\beta$ of $0.14$ per cent at $1\,M_\odot$, $9.6$ per cent at $10\,M_\odot$ and $42$ per cent at $50\,M_\odot$. Whether $\beta$ is constant in those stars is not tested in this module.
10.2 Check 2: the convection zone as $n=3/2$
To test a polytropic index locally, define it locally. The first step of the proof of Proposition 4 used only the ideal gas law (3.2) and the closure (4.1), and gave $T\propto P^{1/(n+1)}$ when $\mu$ is constant. The logarithmic slope of $T$ against $P$ is therefore constant and fixes $n$.
Definition 5 (temperature gradient and effective index). The dimensionless temperature gradient of a stratified star is
$$ \nabla \;\equiv\; \frac{d\ln T}{d\ln P}, \tag{10.1}$$
evaluated along the radius. The effective polytropic index at each radius is
$$ n_{\rm eff} \;\equiv\; \frac{1}{\nabla} - 1. \tag{10.2}$$
For a polytrope of index $n$ with uniform $\mu$, $\nabla = 1/(n+1)$ and $n_{\rm eff}=n$ everywhere. For the adiabat of a monatomic ideal gas, $\nabla = \nabla_{\rm ad} = (\gamma-1)/\gamma = 0.4$ and $n_{\rm eff} = 3/2$. Where $\mu$ varies with radius, $T\propto\mu P/\rho$ adds a term $d\ln\mu/d\ln P$ to $\nabla$, and $n_{\rm eff}$ is then not the exponent in $P\propto\rho^{1+1/n}$.
The table quotes $P$, $\rho$ and $T$ to four significant figures, so a two-point difference of $\ln T$ against $\ln P$ is dominated by rounding. The script computes $\nabla$ as a least-squares slope over $\pm25$ rows. In the band $0.80$ to $0.92\,R$ the slope is flat, and its median is the plateau value.
Check 2: confirmed. Proposition 5 predicts $n = 3/2$ for an adiabatically stratified monatomic gas. On the plateau from $0.80$ to $0.92\,R$ of BS2005-AGS,OP, $\nabla = 0.3957$, so $n_{\rm eff} = 1.5273$. The ratio to the prediction is $1.0182$: agreement to $1.8$ per cent.
The index was fixed by Proposition 5 before the table was read. The only choice made on the table side is the averaging band. Changing the fitting window from $\pm10$ to $\pm50$ rows moves the plateau $\nabla$ from $0.3957$ to $0.3954$. The band stops at $0.92\,R$ to exclude the outermost layers, where $n_{\rm eff}$ turns upward in Fig. 3. The $1.8$ per cent residual is not attributed to a mechanism here.
The same diagnostic applied through the whole star shows why Check 1 failed. Fig. 3 plots $n_{\rm eff}(r)$ from the centre to $0.97\,R$.
$r/R$
$\nabla$
$n_{\rm eff}$
$L(r)/L$
$0.0500$
$0.3251$
$2.076$
$0.0919$
$0.0999$
$0.3220$
$2.105$
$0.4430$
$0.1496$
$0.3032$
$2.299$
$0.7759$
$0.1996$
$0.2629$
$2.804$
$0.9366$
$0.3016$
$0.2116$
$3.726$
$0.9977$
$0.4494$
$0.1937$
$4.162$
$0.9990$
$0.6491$
$0.2244$
$3.457$
$0.9989$
$0.7001$
$0.2997$
$2.336$
$0.9989$
$0.7500$
$0.3963$
$1.523$
$0.9989$
$0.8493$
$0.3960$
$1.525$
$0.9989$
$0.9203$
$0.3933$
$1.543$
$0.9989$
Inside $0.15\,R$, where the enclosed luminosity is still growing and reaches $0.78$ of the total only at $0.15\,R$, $n_{\rm eff}$ is about $2.1$. Outward it rises through $3$ at $0.217\,R$ to a maximum of $4.16$ at $0.453\,R$. It falls through $3$ again at $0.674\,R$, and then drops to the plateau near $1.52$. The radiative interior passes through Eddington's value twice and remains at it nowhere. No single polytrope describes the Sun: inside $0.70\,R$ the index ranges from $2.05$ to $4.16$. One caveat applies to the innermost values. Hydrogen burning has changed the composition of the core, so $\mu$ varies there (the central $\mu$ from composition is $0.83195$), and by Definition 5 $n_{\rm eff}$ in the core is a statement about $\nabla$, not about the exponent relating $P$ to $\rho$. The convection zone is chemically mixed, and Check 2 is not affected.
10.3 Check 3: the base of the convection zone, against a measurement
Checks 1 and 2 compare a closure with a model. Neither compares anything with the Sun. The profile of Fig. 3 supplies a quantity that can be measured: the radius at which the stratification changes from radiative to adiabatic, the base of the convection zone.
The script locates it as the innermost radius, between $0.60$ and $0.95\,R$, at which $\nabla$ reaches $0.99$ of its plateau value $0.3957$. Applied to BS2005-AGS,OP, this gives $r_b = 0.7269\,R$. Table 1 of Bahcall, Serenelli & Basu (2005) gives $0.7280\,R$ for the same model, a ratio of $0.99851$. The locator therefore reproduces the model's published boundary to $0.15$ per cent.
Helioseismology measures the same boundary in the real Sun. The frequencies of the Sun's acoustic oscillations depend on the sound-speed profile, and the sound speed changes slope where $\nabla$ changes from its radiative to its adiabatic value. Basu & Antia (2004) obtain $r_b = (0.7133\pm0.0005)\,R$ from GONG p-mode frequencies, and $(0.7132\pm0.0005)\,R$ from MDI frequencies.
Upper panel: the effective polytropic index $n_{\rm eff} = 1/\nabla - 1$ of the tabulated solar model, against the two indices derived in §5. The shaded region is the convection zone, where $n_{\rm eff}$ lies on the $n=3/2$ line to $1.8$ per cent (Check 2). The radiative interior crosses $n=3$ twice and remains on it nowhere (Check 1). Lower panel: $\nabla$ near the base of the convection zone. The orange line is the model's published base, $0.7280\,R$; the violet line and band are the helioseismic measurement, $(0.7133\pm0.0005)\,R$. The gap is $0.0147\,R$, $29$ standard deviations of the measurement (Check 3).
Check 3: the model is refuted by a measurement. The model's published base, $0.7280\,R$, lies $0.0147\,R$ above the measured $0.7133\pm0.0005\,R$. That is a fractional error of $2.06$ per cent and a departure of $29.4\sigma$. The base recovered by this module's locator, $0.7269\,R$, is $27.2\sigma$ from the measurement. Bahcall, Serenelli & Basu quote the measurement rounded to $0.713\pm0.001$; against that uncertainty the published base departs by $15.0\sigma$.
The number of standard deviations would mean little if the locator could move the base by a comparable amount. Across a grid of fitting windows of $\pm10$, $\pm25$ and $\pm50$ rows and thresholds of $0.98$, $0.99$ and $0.995$ of the plateau, the recovered base ranges from $0.7255$ to $0.7308\,R$. The spread, $0.0054\,R$, is $0.37$ of the model–measurement gap. The grid value nearest the measurement, $0.7255\,R$, is still $24.3\sigma$ from it.
State exactly what this check refutes. It does not refute hydrostatic equilibrium, and it does not refute the $n=3/2$ closure of the convection zone. Both are the instrument: they define the plateau and locate the boundary, and Check 2 confirms the closure in the region where it is claimed. What fails is an input of this particular model, the physics that sets where $\nabla$ first reaches $\nabla_{\rm ad}$. In a radiative region, that radius is set by the opacity and by the heavy-element abundance on which the opacity depends. BS2005-AGS,OP was built with the lower heavy-element abundances denoted "AGS" in its name, and Table 1 of Bahcall, Serenelli & Basu (2005) lists its base, $0.7280$, beside the measured $0.713\pm0.001$ that the paper quotes in its equation (1). This disagreement between low-abundance solar models and helioseismology is known as the solar abundance problem.
The measured base is from Basu & Antia (2004), ApJ606, L85, with the GONG value adopted. The model's published base and surface helium are from Table 1 of Bahcall, Serenelli & Basu (2005). The sentence that the base is set by opacity and heavy-element abundance is standard stellar-structure physics; this module does not compute it, and Module 6 derives the criterion that sets the boundary. The inference of $r_b$ from oscillation frequencies is itself a model-dependent inversion, described in the source and not reproduced here.
11. What this model captures, and what it misses
The equations of this module are (2.4) with the closure (4.1). The table names each idealisation and the module that repays it.
Idealisation made here
What it costs
Where it is repaid
Static: $\mathbf{u}=\mathbf{0}$
Excludes oscillations, convective motions and flows. §2.1 showed the error in the mean structure is below $10^{-23}$ of gravity, but the small perturbations about the equilibrium are where helioseismology gets its data
Module 4 linearises Euler and continuity about this equilibrium. The acoustic cutoff frequency $\nu_{\rm ac} = c_s/(4\pi H)$ contains the scale height $H$ of Proposition 3. Module 9's Proposition 1 shows the idealisation cannot hold at all for a corona in a $1/r$ potential, where an isothermal hydrostatic atmosphere keeps a non-zero pressure at infinity and a wind follows. Module 10's Check 3 prices it where it does hold: treating the Perseus cluster core as static costs $3.98$ per cent of its mass, or $7.66$ per cent if the large-scale shear doubles the turbulent pressure
A closure $P=K\rho^{1+1/n}$ asserted, not derived
Check 2 shows it holds with $n=3/2$ in the convection zone. Check 1 and Fig. 3 show that no single index holds through the Sun
Module 6 derives the Schwarzschild criterion, which decides where $\nabla=\nabla_{\rm ad}$ and therefore where $n=3/2$ applies
Uniform entropy in the convection zone (Proposition 5)
Correct where convection is efficient; wrong in the thin superadiabatic layer below the photosphere, and $\gamma\ne5/3$ in ionisation zones
Module 6
No rotation
Excludes centrifugal flattening and the rotational support of discs
Module 11 (accretion discs)
No magnetic field
Omits magnetic pressure and tension, which can support gas against gravity
Module 12
Spherical symmetry
Excludes binaries, tidal distortion and any non-spherical structure
Not repaid in this module
Ideal gas (3.2), with radiation only in Proposition 6
Wrong where electrons are degenerate, as in white dwarfs, where the closure becomes the degenerate equation of state; the tabulated $\mu$ already differs from the ideal-gas value by $0.41$ per cent at the Sun's centre
Module 13 (radiation); degeneracy is used only in Problem D1
One limitation belongs to the checks rather than to the equations. Every comparison in §10 is with a model built from (2.4). Only Check 3 reaches the real Sun, and it reaches it through a helioseismic inversion that itself uses a solar model as a reference. The agreement of Check 2 is a statement about BS2005-AGS,OP; it becomes a statement about the Sun only to the extent that the model's convection zone matches the Sun's, which helioseismology tests separately.
11.1 Reproducing the numbers
Every number quoted in this module is produced by figs/m03_numbers.py, every problem answer by figs/m03_problems_check.py, and every figure by figs/m03_build_figs.py. Checks 2 and 3 need only the following, run from afd/figs.
"""Check 2 and Check 3 from the tabulated model BS2005-AGS,OP.
grad = d ln T / d ln P, a least-squares slope over +/- 25 rows;
n_eff = 1/grad - 1. Run from afd/figs.
"""
import m03_numbers as M
ssm = M.load_ssm()
grad, neff = M.effective_polytropic_index(ssm, half=25)
rb, plateau, n_plat = M.convection_zone_base(ssm, grad)
print(f"plateau grad = {plateau:.4f} n_eff = {n_plat:.4f}"
f" ratio to 3/2 = {n_plat/1.5:.4f}")
sigma = (rb - M.RCZ_MEAS)/M.RCZ_MEAS_ERR
print(f"recovered base = {rb:.4f} R departure = {sigma:.1f} sigma")
assert abs(n_plat/1.5 - 1.0) < 0.02
Running it prints a plateau $\nabla = 0.3957$, $n_{\rm eff} = 1.5273$, a ratio of $1.0182$ to $3/2$, a recovered base of $0.7269\,R$, and a departure of $27.2\sigma$ from the measurement. The assertion fails if the convection zone departs from $n=3/2$ by more than two per cent.
12. Problems
Solutions are given for every problem. Each numerical answer is verified by figs/m03_problems_check.py.
C1. The tabulated solar model satisfies the virial theorem (8.2) to a ratio of $1.00049$. A student claims this confirms that the Sun is in hydrostatic equilibrium to five parts in ten thousand. Explain what the ratio does and does not show, and why Module 2's solar-wind check was not open to the same objection.
Solution
Proposition 12 derives (8.2) from (2.4) alone. Any structure that satisfies (2.4) with $P(R)=0$ satisfies (8.2). The solar model was computed by solving (2.4), so it satisfies (8.2) whether or not the real Sun does. The ratio $3\int P\,dV/(-\Omega) = 6.1594\times10^{48}/6.1564\times10^{48} = 1.00049$ shows that the table was read and integrated correctly. It is a test of the pipeline, not of nature. It would test hydrostatic equilibrium only if $P$, $\rho$ and $m$ had been measured independently of (2.4). In Module 2 the density exponent $\alpha$ and the speed exponent $\beta$ were fitted separately to Helios data, with no continuity constraint imposed in the fit, so $\alpha-\beta = 2$ was not built into the numbers compared.
C2. Check 1 found that $n=3$ predicts $\rho_c$ and $P_c$ in ratios $0.5069$ and $0.5316$ to the tabulated model, but $T_c$ in ratio $1.0486$. Explain why the temperature is so much better, and state what this implies about using a central temperature to validate a stellar model.
Solution
At fixed $\mu$, the ideal gas law (3.2) gives $T = \mu m_u P/(k_B\rho)$. So the predicted-to-tabulated temperature ratio is the pressure ratio divided by the density ratio: $0.5316/0.5069 = 1.0486$. The two structural errors have the same sign and nearly the same size, so their quotient is close to one. The implication is that $T_c$ is a weak test of a structure. A model can be wrong by a factor of two in both $\rho_c$ and $P_c$ and still give $T_c$ to five per cent, so agreement in $T_c$ alone does not show that the density and pressure are right.
C3. The intracluster medium of §3 has $H = 644$ kpc at $r = 1$ Mpc. Explain why the exponential atmosphere (3.3) does not describe it, with a number, and compare the solar photosphere.
Solution
Proposition 3 assumed $g$ uniform. For gravity dominated by a mass $M$ inside $r$, $g\propto r^{-2}$, so across one scale height $g$ falls by the factor $(r/(r+H))^2$. For the cluster, $H = 0.64\,r$, and the factor is $0.370$: gravity falls to $37$ per cent of its value within one e-folding of the pressure, so the constant $H$ used to define that e-folding has no single value. The spherical equation (2.4) must be used instead. At the solar photosphere, $H/R = 1.94\times10^{-4}$ and the factor is $0.999613$, so $g$ changes by $3.9\times10^{-4}$ across one scale height and (3.3) holds to that precision.
D1. Use the mass–radius relation (7.3) to give the exponent of $M$ in $R\propto M^{x}$ for $n = 0$, $1$, $3/2$, $2$, $5/2$ and $4$ at fixed $K$. Interpret $n=0$, $n=1$ and $n=3/2$, and explain what happens at $n=3$.
Solution
By (7.3), $x = (1-n)/(3-n)$. This gives $x = +1/3$ for $n=0$, $0$ for $n=1$, $-1/3$ for $n=3/2$, $-1$ for $n=2$, $-3$ for $n=5/2$, and $+3$ for $n=4$. For $n=0$ the density is uniform and $K$ plays no role, so $M\propto R^3$ at fixed density and $R\propto M^{1/3}$, as for a solid. For $n=1$ the radius does not depend on the mass. For $n=3/2$, $R\propto M^{-1/3}$: a more massive star is smaller. That is the mass–radius relation of a white dwarf supported by non-relativistic degenerate electrons, whose equation of state has the form $P\propto\rho^{5/3}$ with $K$ fixed by fundamental constants; the form of that equation of state is quoted here, not derived. At $n=3$ the exponent is undefined. From (6.1) with $n=3$, $a^2\propto K\rho_c^{-2/3}$, so $a^3\rho_c\propto K^{3/2}$, and (7.1) gives $M\propto K^{3/2}$ with no dependence on $\rho_c$ or $R$. At fixed $K$ an $n=3$ polytrope has one mass only.
D2. A white dwarf has $M = 0.6\,M_\odot$ and $R = 0.013\,R_\odot$. Use Proposition 11 to bound its central pressure from below, and give the ratio of this bound to the Sun's.
Solution
$M = 1.1930\times10^{33}$ g. With $R_\odot = 6.9598\times10^{10}$ cm, $R = 9.0477\times10^{8}$ cm $= 9048$ km. Then
$$ P_c \;\ge\; \frac{GM^2}{8\pi R^4} \;=\; 5.6405\times10^{21}\ \mathrm{dyn\,cm^{-2}}. $$
The Sun's bound is $4.4750\times10^{14}\ \mathrm{dyn\,cm^{-2}}$. Since the bound scales as $M^2/R^4$, the ratio is $(0.6)^2/(0.013)^4 = 1.2605\times10^{7}$, which the direct division reproduces. The bound needs no equation of state. For the Sun the tabulated $P_c$ exceeds its bound by $522.5$, so the bound lies a factor $522.5$ below the actual value and cannot fix the structure.
D3. Derive the free-fall time (1.1). A uniform sphere of mass $M$ and radius $R$ has no pressure and starts at rest. Show that its outer shell reaches the centre in $t_{\rm ff} = \sqrt{3\pi/(32G\bar\rho)}$, and evaluate it for the Sun. By what factor does it exceed the estimate $\sqrt{R^3/GM}$?
Solution
With no pressure, each shell feels only the mass inside it, and shells do not cross, so the outer shell at radius $r(t)$ obeys $\ddot r = -GM/r^2$ with $r(0) = R$ and $\dot r(0)=0$. Energy conservation gives $\dot r^2 = 2GM(1/r - 1/R)$. Substitute $r = R\cos^2\eta$, so $dr = -2R\cos\eta\sin\eta\,d\eta$ and $1/r - 1/R = \tan^2\eta/R$. Then
$$ dt \;=\; -\frac{dr}{|\dot r|} \;=\; \frac{2R\cos\eta\sin\eta\,d\eta}{\sqrt{2GM/R}\,\tan\eta} \;=\; \sqrt{\frac{2R^3}{GM}}\,\cos^2\eta\,d\eta. $$
The shell reaches the centre at $\eta = \pi/2$, and $\int_0^{\pi/2}\cos^2\eta\,d\eta = \pi/4$. So $t_{\rm ff} = (\pi/4)\sqrt{2R^3/(GM)} = \sqrt{\pi^2R^3/(8GM)}$. With $M = \tfrac{4}{3}\pi R^3\bar\rho$ this is $\sqrt{3\pi/(32G\bar\rho)}$. For the Sun, $\bar\rho = 1.4081\ \mathrm{g\,cm^{-3}}$, so $t_{\rm ff} = 1770.3$ s $= 29.5$ min. The estimate $\sqrt{R^3/GM}$ is $1593.8$ s, and the ratio is $1.1107$, which is $\pi/(2\sqrt2)$ exactly.
K1. Compute the isothermal scale height (3.3) at the surfaces of Mars and Venus from $T = 210$ K, $\mu = 43.34$, $g = 371\ \mathrm{cm\,s^{-2}}$ for Mars, and $T = 737$ K, $\mu = 43.45$, $g = 887\ \mathrm{cm\,s^{-2}}$ for Venus. Compare them with the Earth's $8.435$ km, and say which inputs set the Venus-to-Mars ratio.
Solution
From $H = k_BT/(\mu m_u g)$: for Mars, $H = 1.380649\times10^{-16}\times210/(43.34\times1.6605\times10^{-24}\times371) = 10.859$ km. For Venus, $H = 15.900$ km. The ratios to the Earth are $1.2874$ for Mars and $1.8851$ for Venus, and Venus to Mars is $1.4642$. The two atmospheres are both mostly carbon dioxide and have $\mu$ within $0.25$ per cent of each other, so the ratio is set by $T/g$: the temperature ratio $3.5095$ divided by the gravity ratio $2.3908$ is $1.4679$, and the ratio of $\mu$ accounts for the remaining difference from $1.4642$.
The planetary temperatures, molecular weights and surface gravities are representative round values supplied as inputs. They are not checked against a source here, and the problem tests the formula, not the planets.
K2. Model a neutron star of $M = 1.4\,M_\odot$ and $R = 12$ km as an $n=1$ polytrope. Find $K$, $\rho_c$, $\rho_c/\bar\rho$ and $P_c$, compare $\rho_c$ with nuclear saturation density, and check $P_c$ against Proposition 11.
Solution
For $n=1$, (6.1) gives $a^2 = 2K/(4\pi G) = K/(2\pi G)$, independent of $\rho_c$, and $R = a\xi_1 = \pi a$. So $K = 2\pi G R^2/\pi^2 = 2GR^2/\pi = 6.1185\times10^{4}$ in CGS units. From (7.2) with $\xi_1 = \pi$ and $|\theta_1'| = 1/\pi$, $D_1 = \pi^2/3 = 3.2899$, so $\rho_c = (\pi^2/3)\cdot3M/(4\pi R^3) = \pi M/(4R^3)$. With $M = 2.7838\times10^{33}$ g and $R = 1.2\times10^{6}$ cm, $\rho_c = 1.2653\times10^{15}\ \mathrm{g\,cm^{-3}}$ and $\bar\rho = 3.8459\times10^{14}\ \mathrm{g\,cm^{-3}}$. Nuclear saturation density, $0.16$ nucleons per $\mathrm{fm^3}$ times $m_u$, is $2.657\times10^{14}\ \mathrm{g\,cm^{-3}}$, so $\rho_c$ is $4.76$ times it. The central pressure is $P_c = K\rho_c^2 = 9.7951\times10^{34}\ \mathrm{dyn\,cm^{-2}}$. The bound of Proposition 11 is $9.9245\times10^{33}\ \mathrm{dyn\,cm^{-2}}$, and $P_c$ exceeds it by $9.87$. This must equal $8\pi W_1 = 8\pi\times0.39270 = \pi^2 = 9.87$, and it does.
The saturation density $0.16\ \mathrm{fm^{-3}}$ is a standard nuclear-physics value quoted as an input, not checked against a source here. The $n=1$ polytrope is a toy: it ignores general relativity, which is not small in a neutron star.
K3. For an $n=3$ star, (7.1) fixes $K$ from $M$, and (5.2) then fixes $\beta$. Taking $\mu = 0.832$, solve for $1-\beta$ at $1$, $10$ and $50\,M_\odot$, and compare the $1\,M_\odot$ value with the Sun's central $1-\beta$ of (9.1).
Solution
For $n=3$, $a^2 = K\rho_c^{-2/3}/(\pi G)$, so $a^3\rho_c = (K/\pi G)^{3/2}$ and (7.1) gives $M = 4\pi(K/\pi G)^{3/2}\mu_1$ with $\mu_1 = 2.01824$. Inverting, $K = \pi G\,[M/(4\pi\mu_1)]^{2/3}$. Rearranging (5.2),
$$ \frac{1-\beta}{\beta^4} \;=\; \frac{a_{\rm rad}}{3}\,K^3\left(\frac{\mu m_u}{k_B}\right)^{4}, $$
which is Eddington's quartic. Its left side decreases monotonically from $+\infty$ to $0$ as $\beta$ rises from $0$ to $1$, so it has one root, found by bisection. At $1\,M_\odot$, $K = 3.8409\times10^{14}$, the right side is $1.4329\times10^{-3}$, and $1-\beta = 1.4247\times10^{-3}$. At $10\,M_\odot$, $K = 1.7828\times10^{15}$, the right side is $1.4329\times10^{-1}$, and $1-\beta = 9.5784\times10^{-2}$. At $50\,M_\odot$, $K = 5.2130\times10^{15}$, the right side is $3.5821$, and $1-\beta = 4.1617\times10^{-1}$. So radiation carries $0.14$ per cent of the pressure at $1\,M_\odot$, $9.6$ per cent at $10\,M_\odot$ and $42$ per cent at $50\,M_\odot$; the share rises by a factor $67.2$ from $1$ to $10\,M_\odot$. Eddington's $1\,M_\odot$ value is $2.30$ times the tabulated central value $6.194\times10^{-4}$. Both are below one per cent. Proposition 6 does not require a large share, so this does not by itself refute $n=3$; Check 1 refutes it through the measured variation of $\nabla$. What K3 shows is that at one solar mass radiation pressure cannot be the agent that fixes the structure, while at $10$ and $50\,M_\odot$ it carries a share large enough to matter in (2.4).
Notation
Symbol
Meaning
Units (CGS)
First used
$P$, $\rho$, $T$
pressure, mass density, temperature
$\mathrm{dyn\,cm^{-2}}$, $\mathrm{g\,cm^{-3}}$, K
(2.1)
$\mathbf{u}$, $\mathbf{g}$
bulk velocity, gravitational acceleration
$\mathrm{cm\,s^{-1}}$, $\mathrm{cm\,s^{-2}}$
(2.1)
$m(r)$
mass inside radius $r$
g
(2.3)
$M$, $R$, $L$
stellar mass, radius, luminosity
g, cm, $\mathrm{erg\,s^{-1}}$
(1.2)
$\bar\rho$
mean density $3M/(4\pi R^3)$
$\mathrm{g\,cm^{-3}}$
(1.1)
$t_{\rm ff}$, $t_{\rm KH}$
free-fall time, Kelvin–Helmholtz time
s
(1.1), (1.2)
$\mu$, $m_u$
mean molecular weight, atomic mass unit
—, g
(3.2)
$X$, $Y$, $Z$
mass fractions of hydrogen, helium, heavier elements
Every published value this module compares against, with what was taken from it.
The tabulated standard solar model. J. N. Bahcall, A. M. Serenelli and S. Basu, New solar opacities, abundances, helioseismology, and neutrino fluxes, The Astrophysical Journal621, L85–L88 (2005), arXiv:astro-ph/0412440. Quoted: Table 1 for BS05(AGS,OP), convection-zone base $R_{\rm CZ}/R_\odot = 0.7280$ and surface helium $Y = 0.229$, and the helioseismic value quoted in the paper, $0.713\pm0.001$. The table data/bs05_agsop.dat was downloaded from http://www.sns.ias.edu/~jnb/SNdata/Export/BS2005/bs05_agsop.dat.
The measured base of the convection zone. S. Basu and H. M. Antia, Constraining solar abundances using helioseismology, The Astrophysical Journal606, L85 (2004), arXiv:astro-ph/0403485. Quoted: $r_b = (0.7133\pm0.0005)\,R_\odot$ from GONG and $(0.7132\pm0.0005)\,R_\odot$ from MDI p-mode frequencies.
Polytrope constants. S. Chandrasekhar, An Introduction to the Study of Stellar Structure (University of Chicago Press, 1939), chapter IV, Table 4. Quoted for $n=3$: $\xi_1 = 6.89685$, $-\xi_1^2\theta_1' = 2.01824$, $\rho_c/\bar\rho = 54.1825$.
The standard atmosphere. ISO 2533:1975, Standard Atmosphere, identical below 11 km to the U.S. Standard Atmosphere, 1976. Quoted: $T(0) = 288.15$ K, $P(0) = 101325$ Pa, lapse rate $6.5\ \mathrm{K\,km^{-1}}$, tropopause at 11 km with $216.65$ K and $22632$ Pa.
Solar constants. IAU 2015 Resolution B3, nominal conversion constants: $GM_\odot = 1.3271244\times10^{26}\ \mathrm{cm^3\,s^{-2}}$, $R_\odot = 6.957\times10^{10}$ cm, $L_\odot = 3.828\times10^{33}\ \mathrm{erg\,s^{-1}}$, $T_{\rm eff} = 5772$ K.
Earlier modules. The Euler equation (7.3) and the equipartition identity (Proposition 1) of Module 2; Module 2's solar-wind check, quoted in §8 and Problem C1.