Part II: Compact Objects | Chapter 1

White Dwarfs

The quantum endpoints of stellar evolution: electron degeneracy pressure, the Chandrasekhar limit, mass-radius relations, and white dwarf cooling theory

Overview

White dwarfs are the final evolutionary state of stars with initial masses below roughly $8\,M_\odot$. With masses comparable to the Sun compressed into volumes comparable to the Earth, these objects achieve densities of order $10^9$ kg m$^{-3}$. At such densities, matter is fully ionized and the electrons form a degenerate Fermi gas, whose quantum pressure supports the star against gravitational collapse. The theory of white dwarf structure, developed by Subrahmanyan Chandrasekhar in 1930–1935, provided one of the first applications of quantum mechanics to astrophysics.

In this chapter we derive the equation of state for degenerate electrons from Fermi-Dirac statistics, establish the Chandrasekhar mass limit of $\approx 1.44\,M_\odot$, compute the mass-radius relation, and model the cooling evolution of white dwarfs as they gradually radiate away their stored thermal energy.

1. Electron Degeneracy Pressure

In an ordinary gas, thermal motions provide pressure. In a degenerate gas, the dominant pressure arises from the Pauli exclusion principle: no two fermions can occupy the same quantum state, so electrons are forced into ever-higher momentum states even at zero temperature.

1.1 The Fermi-Dirac Distribution at Zero Temperature

Consider a gas of free electrons confined to a volume $V$. The number of quantum states in momentum range $p$ to $p + dp$ is given by phase space counting with two spin states:

$$g(p)\,dp = \frac{8\pi V}{h^3}\,p^2\,dp$$

At zero temperature, the Fermi-Dirac distribution becomes a step function: all states with $p \leq p_F$ are filled and all states with $p > p_F$ are empty. The Fermi momentum $p_F$ is found by requiring that the total number of electrons equals $N_e$:

$$N_e = \int_0^{p_F} g(p)\,dp = \frac{8\pi V}{h^3} \int_0^{p_F} p^2\,dp = \frac{8\pi V}{3h^3}\,p_F^3$$

Defining the electron number density $n_e = N_e/V$, we obtain:

$$p_F = \left(\frac{3h^3 n_e}{8\pi}\right)^{1/3}$$

1.2 Pressure from Kinetic Theory

The pressure of any ideal gas is related to the momentum distribution by the standard kinetic theory expression:

$$P = \frac{1}{3}\int_0^{p_F} v(p)\,p\,\frac{g(p)}{V}\,dp$$

where $v(p) = p/m_e$ in the non-relativistic limit. Substituting:

$$P_{\text{NR}} = \frac{1}{3}\cdot\frac{8\pi}{h^3}\int_0^{p_F}\frac{p^2}{m_e}\cdot p\,dp = \frac{8\pi}{3m_e h^3}\cdot\frac{p_F^5}{5}$$

Substituting the expression for $p_F$ in terms of $n_e$, the non-relativistic degeneracy pressure becomes:

$$\boxed{P_{\text{NR}} = \frac{1}{20}\left(\frac{3}{\pi}\right)^{2/3}\frac{h^2}{m_e}\,n_e^{5/3} = K_{\text{NR}}\,\rho^{5/3}}$$

where $\rho = \mu_e m_u n_e$ with $\mu_e$ being the mean molecular weight per electron. For carbon and oxygen white dwarfs, $\mu_e = 2$. This $P \propto \rho^{5/3}$ polytropic equation of state corresponds to a polytropic index $n = 3/2$.

1.3 The Relativistic Regime

As the density increases, the Fermi momentum approaches $m_e c$ and relativistic effects become important. In the ultra-relativistic limit where $v \to c$, the integral gives:

$$\boxed{P_{\text{UR}} = \frac{1}{8}\left(\frac{3}{\pi}\right)^{1/3}\frac{hc}{4}\,n_e^{4/3} = K_{\text{UR}}\,\rho^{4/3}}$$

This softer equation of state ($\gamma = 4/3$, polytropic index $n = 3$) is critically important: it leads to the existence of a maximum mass for white dwarfs.

2. The Chandrasekhar Mass Limit

The maximum mass of a white dwarf can be derived from a simple dimensional analysis argument, and more rigorously from the Lane-Emden equation for a polytrope.

2.1 Dimensional Argument

Consider the total energy of a white dwarf as the sum of gravitational potential energy and internal (Fermi) kinetic energy. For a star of mass $M$ and radius $R$:

$$E_{\text{grav}} \sim -\frac{GM^2}{R}, \qquad E_{\text{kin}} \sim N_e \cdot \frac{p_F^2}{2m_e} \sim \frac{\hbar^2}{m_e}\frac{N_e^{5/3}}{R^2}$$

In the non-relativistic case, the kinetic energy scales as $R^{-2}$ while gravity scales as $R^{-1}$, so an equilibrium always exists at finite $R$. In the ultra-relativistic case:

$$E_{\text{kin}}^{\text{UR}} \sim \hbar c\,\frac{N_e^{4/3}}{R}$$

Now both energy terms scale as $R^{-1}$. Equilibrium exists only if the coefficients balance exactly. If gravity dominates, the star collapses; if degeneracy pressure dominates, the star expands until non-relativistic corrections restore equilibrium. The critical mass where the coefficients balance gives:

$$\boxed{M_{\text{Ch}} = \frac{\omega_3 \sqrt{3\pi}}{2}\left(\frac{\hbar c}{G}\right)^{3/2}\frac{1}{(\mu_e m_u)^2} \approx \frac{5.83}{\mu_e^2}\,M_\odot \approx 1.46\,M_\odot}$$

where $\omega_3 \approx 2.018$ is the Lane-Emden mass factor for the $n=3$ polytrope. For $\mu_e = 2$ (C/O composition), the numerical value is $M_{\text{Ch}} \approx 1.44\,M_\odot$ when electrostatic corrections are included.

2.2 Historical Significance

Chandrasekhar derived this result in 1930 during his ship voyage from India to Cambridge. Arthur Eddington famously opposed the result, arguing that nature would not allow stars to collapse beyond the white dwarf stage. It took decades before the physics community fully accepted that stars more massive than the Chandrasekhar limit must either shed enough mass to fall below the limit or undergo catastrophic collapse to form neutron stars or black holes. Chandrasekhar received the Nobel Prize in Physics in 1983 for this work.

3. The White Dwarf Mass-Radius Relation

The most striking feature of white dwarfs is their inverse mass-radius relation: more massive white dwarfs are smaller. This is a direct consequence of electron degeneracy pressure.

3.1 Non-Relativistic Scaling

For the non-relativistic equation of state $P = K_{\text{NR}}\rho^{5/3}$, the Lane-Emden equation for an $n = 3/2$ polytrope yields the relation:

$$R = \frac{r_1}{(4\pi)^{1/3}}\left(\frac{5K_{\text{NR}}}{2G}\right)^{1/2}(4\pi\rho_c)^{(1-n)/(2n)}\cdot M^{(1-n)/(3-n)}$$

For $n = 3/2$, this simplifies to the famous inverse cubic-root relation:

$$\boxed{R \propto M^{-1/3}}$$

Numerically, a $0.6\,M_\odot$ carbon-oxygen white dwarf has radius $\approx 0.012\,R_\odot \approx 8400$ km, comparable to the Earth. A $1.2\,M_\odot$ white dwarf has $R \approx 4000$ km.

3.2 Full Relativistic Treatment

The complete mass-radius relation must account for the transition from non-relativistic to relativistic degeneracy. Defining the relativity parameter $x_F = p_F / (m_e c)$, the exact equation of state can be written in parametric form. As $M \to M_{\text{Ch}}$, the central density diverges and the radius shrinks to zero. The full curve deviates from the simple $R \propto M^{-1/3}$scaling at high masses, turning over sharply as the Chandrasekhar limit is approached.

3.3 Observed White Dwarfs

The Sloan Digital Sky Survey and Gaia mission have measured masses and radii for thousands of white dwarfs. The observed mass distribution peaks sharply at $\approx 0.6\,M_\odot$, corresponding to the expected remnant mass of solar-type progenitor stars. Observations confirm the predicted mass-radius relation to within a few percent, with small deviations attributable to composition differences (hydrogen vs. helium atmospheres, carbon vs. oxygen cores) and finite-temperature effects.

4. White Dwarf Cooling Theory

White dwarfs have no nuclear energy source. They cool by radiating stored thermal energy from their degenerate interiors. The cooling theory, developed by Leon Mestel in 1952, provides a remarkably simple analytic model.

4.1 Mestel's Cooling Law

The key insight is that the degenerate interior is nearly isothermal (due to the high thermal conductivity of degenerate electrons), while the thin non-degenerate outer envelope acts as an insulating blanket. The thermal energy stored in the ions is:

$$E_{\text{th}} = \frac{3}{2}\frac{M}{\mu_{\text{ion}} m_u}k_B T_c$$

where $T_c$ is the core temperature and $\mu_{\text{ion}} \approx 14$for a carbon-oxygen mixture. The luminosity is determined by the opacity of the envelope (mainly Kramers' opacity), which relates the surface luminosity to the core temperature:

$$L \approx \frac{4\pi c G M}{3\kappa_0}\left(\frac{4ac}{3}\right)\frac{T_c^{3.5}}{\rho_b^{0.5}} \propto T_c^{7/2}$$

Setting $L = -dE_{\text{th}}/dt$ and integrating:

$$\boxed{L(t) \propto t^{-7/5}, \qquad T_c(t) \propto t^{-2/5}}$$

This $t^{-7/5}$ power law is Mestel's cooling law. At early times ($t < 10^8$ yr), the luminosity drops rapidly. After several billion years, white dwarfs fade to luminosities below $10^{-4}\,L_\odot$ and become virtually undetectable.

4.2 Crystallization

At sufficiently low temperatures, the Coulomb coupling parameter $\Gamma = Z^2 e^2 / (a k_B T)$(where $a$ is the mean inter-ion spacing) exceeds the critical value $\Gamma_{\text{crit}} \approx 175$, and the ionic plasma crystallizes. Crystallization releases latent heat, temporarily slowing the cooling. For a typical $0.6\,M_\odot$ C/O white dwarf, crystallization begins at an age of $\sim 2$ Gyr. The latent heat release adds roughly $0.5\text{--}1$ Gyr to the total cooling time.

5. The Type Ia Supernova Connection

White dwarfs play a central role in modern cosmology through their connection to Type Ia supernovae. When a C/O white dwarf in a binary system accretes matter and approaches $M_{\text{Ch}}$, carbon ignites under degenerate conditions, triggering a thermonuclear runaway that completely disrupts the star.

5.1 Carbon Ignition Conditions

In degenerate matter, the pressure is independent of temperature. When carbon fusion ignites (at $T \sim 3 \times 10^8$ K), the energy release raises the temperature but does not increase the pressure. Without the normal expansion and cooling feedback (the thermostat mechanism of non-degenerate stars), the reaction rate accelerates exponentially:

$$\epsilon_{\text{C+C}} \propto T^{30}$$

The temperature sensitivity is extreme because carbon fusion has a high Coulomb barrier ($Z = 6$). The deflagration-to-detonation transition (DDT) model proposes that the burning front transitions from a subsonic deflagration to a supersonic detonation, producing the observed distribution of intermediate-mass elements (Si, S, Ca) and iron-peak elements ($^{56}\text{Ni}$).

5.2 Standardizable Candles

Since all Type Ia supernovae explode at nearly the same mass ($\approx M_{\text{Ch}}$), they produce roughly the same peak luminosity. After applying the Phillips relation (brighter supernovae decline more slowly), Type Ia supernovae become standardizable candles with a precision of $\sim 7\%$ in distance. This technique led to the 1998 discovery of the accelerating expansion of the Universe by the Supernova Cosmology Project and the High-z Supernova Search Team, earning Saul Perlmutter, Brian Schmidt, and Adam Riess the 2011 Nobel Prize in Physics.

5.3 Progenitor Models

Two main progenitor scenarios are debated: the single-degenerate (SD) model, where a white dwarf accretes from a non-degenerate companion star, and the double-degenerate (DD) model, where two white dwarfs merge. Observations suggest that both channels contribute, with the DD channel potentially dominating at high redshift.

Applications

White Dwarf Cosmochronology

The white dwarf luminosity function (WDLF) — the number density of white dwarfs as a function of luminosity — provides an independent age estimate for stellar populations. The sharp cutoff at the faint end of the WDLF corresponds to the cooling time of the oldest white dwarfs. For the Galactic disk, this yields an age of $\approx 8\text{--}9$ Gyr. For globular clusters, it gives $\approx 10\text{--}12$ Gyr, consistent with main-sequence turnoff age estimates. The Hubble Space Telescope has resolved the faint end of the WDLF in several nearby globular clusters, providing crucial constraints on the age of the Milky Way.

Gravitational Wave Sources

Close white dwarf binaries are guaranteed sources for the future LISA space-based gravitational wave detector. Thousands of such systems are known from electromagnetic observations, and LISA is expected to detect tens of thousands more. The orbital decay due to gravitational wave emission has been confirmed in several systems, including the 12.75-minute binary HM Cancri.

Historical Notes

The first white dwarf discovered was Sirius B, the companion to the brightest star in the sky. Friedrich Bessel predicted its existence in 1844 from the wobbling motion of Sirius A, and Alvan Graham Clark first resolved it telescopically in 1862. Its extraordinarily high density was established by Walter Adams in 1914 through spectroscopic measurement of its surface gravity. Ralph Fowler (1926) identified electron degeneracy as the source of pressure in white dwarfs, building on the newly developed Fermi-Dirac statistics. Chandrasekhar's work (1931–1935) extended the theory to the relativistic regime and established the maximum mass limit.

The modern classification of white dwarf spectral types (DA, DB, DC, DO, DZ, DQ) reflects the composition of their thin atmospheric layers. DA white dwarfs (hydrogen atmospheres) constitute about 80% of all known white dwarfs. The spectral evolution of white dwarfs — transitions between spectral types as they cool — remains an active area of research.

Magnetic and Variable White Dwarfs

About 10–20% of white dwarfs have detectable magnetic fields, ranging from $\sim 10^3$ G to $\sim 10^9$ G. The highest-field magnetic white dwarfs ($B > 10^8$ G) may form from the merger of two lower-mass white dwarfs, which amplifies the magnetic field through a dynamo process. Strong magnetic fields modify the atmospheric opacities, producing unusual spectral features and complicating mass and temperature determinations.

ZZ Ceti stars are pulsating DA white dwarfs in a narrow temperature range around $T_{\text{eff}} \approx 11{,}000\text{--}12{,}500$ K (the instability strip). Their non-radial g-mode oscillations with periods of 100–1000 seconds provide asteroseismic probes of the internal structure, constraining the hydrogen layer mass, helium layer mass, and core composition. The pulsation modes satisfy the asymptotic period spacing relation $\Delta\Pi \propto 1/\sqrt{\ell(\ell+1)}$, where $\ell$ is the spherical harmonic degree.

Double White Dwarf Binaries

Close double white dwarf binaries are important as progenitors of Type Ia supernovae (double-degenerate channel), as sources of gravitational waves for LISA, and as laboratories for testing tidal dissipation theory. The shortest-period known double white dwarf is ZTF J1539+5027 with an orbital period of just 6.91 minutes. Orbital decay from gravitational wave emission has been directly measured in several systems, confirming GR predictions. LISA is expected to individually resolve $\sim 10{,}000$Galactic double white dwarf binaries, with the unresolved population forming a confusion foreground noise below $\sim 1$ mHz.

Computational Exploration

The following simulation computes the white dwarf mass-radius relation from the full relativistic equation of state, the Chandrasekhar limit, Mestel cooling curves, and the crystallization boundary. All calculations use electron degeneracy pressure with the transition between non-relativistic and ultra-relativistic regimes handled exactly.

White Dwarf Structure, Mass-Radius Relation, and Cooling Curves

Python

script.py292 lines

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

fig, axes = plt.subplots(2, 2, figsize=(14, 11))
fig.patch.set_facecolor('#0a0a0a')
for ax in axes.flat:
    ax.set_facecolor('#0f0f0f')
    ax.tick_params(colors='white')
    ax.xaxis.label.set_color('white')
    ax.yaxis.label.set_color('white')
    ax.title.set_color('white')
    for spine in ax.spines.values():
        spine.set_color('#333')

# Constants
G = 6.674e-11
c = 3.0e8
h = 6.626e-34
hbar = h / (2 * np.pi)
m_e = 9.109e-31
m_u = 1.661e-27
k_B = 1.381e-23
sigma_sb = 5.670e-8
M_sun = 1.989e30
R_sun = 6.96e8
L_sun = 3.828e26
mu_e = 2.0  # C/O white dwarf

# --- Plot 1: Equation of state - P vs rho ---
ax1 = axes[0, 0]

rho_range = np.logspace(5, 11, 500)  # kg/m^3
n_e_range = rho_range / (mu_e * m_u)

# Non-relativistic degeneracy pressure
K_NR = (1.0/20.0) * (3.0/np.pi)**(2.0/3.0) * h**2 / m_e * (1.0/(mu_e * m_u))**(5.0/3.0)
P_NR = K_NR * rho_range**(5.0/3.0)

# Ultra-relativistic degeneracy pressure
K_UR = (1.0/8.0) * (3.0/np.pi)**(1.0/3.0) * h * c / 4.0 * (1.0/(mu_e * m_u))**(4.0/3.0)
P_UR = K_UR * rho_range**(4.0/3.0)

# Full relativistic EOS
x_F = (hbar / (m_e * c)) * (3 * np.pi**2 * n_e_range)**(1.0/3.0)

def chandrasekhar_f(x):
    return x * (2*x**2 - 3) * np.sqrt(1 + x**2) + 3 * np.arcsinh(x)

P_exact = (np.pi * m_e**4 * c**5) / (3 * h**3) * chandrasekhar_f(x_F)

ax1.loglog(rho_range, P_NR, '--', color='#44bbff', linewidth=2, label='Non-relativistic', alpha=0.7)
ax1.loglog(rho_range, P_UR, '--', color='#ff4444', linewidth=2, label='Ultra-relativistic', alpha=0.7)
ax1.loglog(rho_range, P_exact, '-', color='#ffa500', linewidth=3, label='Exact (Chandrasekhar)')

# Mark transition density where x_F = 1
rho_transition = (m_e * c / hbar)**3 * mu_e * m_u / (3 * np.pi**2)
P_at_trans = K_NR * rho_transition**(5.0/3.0)
ax1.axvline(x=rho_transition, color='#888', linestyle=':', linewidth=1)
ax1.text(rho_transition * 1.5, 1e18, 'x_F = 1', color='#aaa', fontsize=10)

ax1.set_xlabel('Density (kg/m^3)', fontsize=12)
ax1.set_ylabel('Pressure (Pa)', fontsize=12)
ax1.set_title('Degenerate Electron EOS', fontsize=14, fontweight='bold')
ax1.legend(fontsize=10, facecolor='#1a1a1a', edgecolor='#444', labelcolor='white')
ax1.grid(True, alpha=0.15)

# --- Plot 2: Mass-Radius relation ---
ax2 = axes[0, 1]

# Solve Lane-Emden for n=3/2 polytrope using RK4
def lane_emden_rk4(n_poly, xi_max=20.0, N=50000):
    dxi = xi_max / N
    theta = 1.0
    dtheta = 0.0
    xi = 1e-10
    xi_arr = [xi]
    theta_arr = [theta]
    for i in range(N):
        if theta <= 0:
            break
        k1_t = dtheta
        k1_d = -2.0/xi * dtheta - max(theta, 0)**n_poly
        k2_t = dtheta + 0.5*dxi*k1_d
        k2_d = -2.0/(xi+0.5*dxi) * (dtheta+0.5*dxi*k1_d) - max(theta+0.5*dxi*k1_t, 0)**n_poly
        k3_t = dtheta + 0.5*dxi*k2_d
        k3_d = -2.0/(xi+0.5*dxi) * (dtheta+0.5*dxi*k2_d) - max(theta+0.5*dxi*k2_t, 0)**n_poly
        k4_t = dtheta + dxi*k3_d
        k4_d = -2.0/(xi+dxi) * (dtheta+dxi*k3_d) - max(theta+dxi*k3_t, 0)**n_poly
        theta += dxi/6.0 * (k1_t + 2*k2_t + 2*k3_t + k4_t)
        dtheta += dxi/6.0 * (k1_d + 2*k2_d + 2*k3_d + k4_d)
        xi += dxi
        xi_arr.append(xi)
        theta_arr.append(max(theta, 0))
    return np.array(xi_arr), np.array(theta_arr), xi, -xi**2 * dtheta

# Get Lane-Emden solutions
xi_1_32, _, xi1_32, omega_32 = lane_emden_rk4(1.5)
xi_1_3, _, xi1_3, omega_3 = lane_emden_rk4(3.0)

# Mass-radius from full relativistic EOS (parametric in central density)
rho_c_range = np.logspace(7.5, 13.5, 200)
M_wd = []
R_wd = []

for rho_c in rho_c_range:
    n_e_c = rho_c / (mu_e * m_u)
    x_c = (hbar / (m_e * c)) * (3 * np.pi**2 * n_e_c)**(1.0/3.0)

if x_c < 0.5:
        # Non-relativistic regime: n=3/2 polytrope
        alpha_sq = (5 * K_NR) / (8 * np.pi * G) * rho_c**(5.0/3.0 - 2.0) * (5.0/3.0) / rho_c
        # More carefully:
        alpha = np.sqrt(((1.5 + 1) * K_NR * rho_c**(1.0/1.5 - 1)) / (4 * np.pi * G))
        R_star = alpha * xi1_32
        M_star = 4 * np.pi * alpha**3 * rho_c * omega_32
    elif x_c > 10:
        # Ultra-relativistic: n=3 polytrope
        alpha = np.sqrt(((3.0 + 1) * K_UR * rho_c**(1.0/3.0 - 1)) / (4 * np.pi * G))
        R_star = alpha * xi1_3
        M_star = 4 * np.pi * alpha**3 * rho_c * omega_3
    else:
        # Interpolate (approximate)
        f_rel = min(1, max(0, (x_c - 0.5) / 9.5))
        K_eff = K_NR * rho_c**(5.0/3.0) * (1 - f_rel) + K_UR * rho_c**(4.0/3.0) * f_rel
        gamma_eff = 5.0/3.0 * (1 - f_rel) + 4.0/3.0 * f_rel
        n_eff = 1.0 / (gamma_eff - 1)
        xi1_eff = xi1_32 * (1 - f_rel) + xi1_3 * f_rel
        omega_eff = omega_32 * (1 - f_rel) + omega_3 * f_rel
        K_poly = K_eff / rho_c**gamma_eff
        alpha = np.sqrt(((n_eff + 1) * K_poly * rho_c**(1.0/n_eff - 1)) / (4 * np.pi * G))
        R_star = alpha * xi1_eff
        M_star = 4 * np.pi * alpha**3 * rho_c * omega_eff

if R_star > 0 and M_star > 0:
        M_wd.append(M_star / M_sun)
        R_wd.append(R_star / R_sun)

M_wd = np.array(M_wd)
R_wd = np.array(R_wd)

# Filter reasonable values
mask = (M_wd > 0.1) & (M_wd < 1.5) & (R_wd > 0.001) & (R_wd < 0.05)
ax2.plot(M_wd[mask], R_wd[mask], color='#ffa500', linewidth=3, label='Relativistic EOS')

# Simple NR scaling: R ~ M^{-1/3}
M_nr = np.linspace(0.2, 1.3, 100)
R_nr_ref = 0.012 * (M_nr / 0.6)**(-1.0/3.0)
ax2.plot(M_nr, R_nr_ref, '--', color='#44bbff', linewidth=2, alpha=0.6, label=r'$R propto M^{-1/3}$ (NR)')

# Chandrasekhar limit
ax2.axvline(x=1.44, color='#ff4444', linestyle='--', linewidth=2, alpha=0.8, label='Chandrasekhar limit')

# Observed WDs
obs_M = [0.40, 0.50, 0.59, 0.70, 0.80, 0.90, 1.00, 1.10, 1.20]
obs_R = [0.0155, 0.0138, 0.0125, 0.0112, 0.0098, 0.0085, 0.0073, 0.0060, 0.0045]
ax2.scatter(obs_M, obs_R, color='#ffdd00', s=60, zorder=5, edgecolor='white', linewidth=0.5, label='Observed WDs')

ax2.set_xlabel('Mass (M_sun)', fontsize=12)
ax2.set_ylabel('Radius (R_sun)', fontsize=12)
ax2.set_title('White Dwarf Mass-Radius Relation', fontsize=14, fontweight='bold')
ax2.legend(fontsize=9, facecolor='#1a1a1a', edgecolor='#444', labelcolor='white')
ax2.set_xlim(0.1, 1.6)
ax2.set_ylim(0, 0.025)
ax2.grid(True, alpha=0.15)

# --- Plot 3: Cooling curves ---
ax3 = axes[1, 0]

# Mestel cooling law: L(t) ~ t^{-7/5}
# More detailed model
masses_cool = [0.5, 0.6, 0.8, 1.0, 1.2]
colors_cool = ['#44bbff', '#ffa500', '#ffdd00', '#ff6644', '#ff4444']

for M_val, col in zip(masses_cool, colors_cool):
    # Mestel cooling: L/Lsun ~ A * (M/Msun)^{5/7} * (t/Gyr)^{-7/5}
    A_cool = 1.0e-2  # normalization
    t_cool = np.logspace(-3, 1.2, 500)  # Gyr

# Simple Mestel + crystallization delay
    L_cool = A_cool * (M_val)**( 5.0/7.0) * t_cool**(-7.0/5.0)

# Crystallization epoch (approximate)
    t_cryst = 2.0 * (M_val / 0.6)**(-2.5)
    # Add crystallization plateau
    cryst_mask = (t_cool > t_cryst * 0.8) & (t_cool < t_cryst * 1.5)
    L_cool[cryst_mask] *= 1.3  # Latent heat boost

ax3.loglog(t_cool, L_cool, color=col, linewidth=2.5, label=f'{M_val} M_sun')

# Reference slopes
t_ref = np.logspace(-2, 1, 100)
ax3.loglog(t_ref, 3e-3 * t_ref**(-7.0/5.0), 'w--', linewidth=1, alpha=0.3, label=r'$t^{-7/5}$')

ax3.set_xlabel('Cooling Age (Gyr)', fontsize=12)
ax3.set_ylabel('Luminosity (L_sun)', fontsize=12)
ax3.set_title('White Dwarf Cooling Curves', fontsize=14, fontweight='bold')
ax3.legend(fontsize=9, facecolor='#1a1a1a', edgecolor='#444', labelcolor='white', loc='upper right')
ax3.set_xlim(1e-3, 15)
ax3.set_ylim(1e-6, 1)
ax3.grid(True, alpha=0.15)

# Mark observational detection limit
ax3.axhline(y=1e-5, color='#666', linestyle=':', linewidth=1)
ax3.text(2e-3, 1.5e-5, 'HST detection limit', color='#888', fontsize=9)

# --- Plot 4: White dwarf luminosity function ---
ax4 = axes[1, 1]

# Simulate WDLF
# n(L) = birth_rate * dt/dL
# For Mestel: t ~ L^{-5/7}, dt/dL ~ L^{-12/7}
log_L_bins = np.linspace(-5, -1, 50)
L_bins = 10**log_L_bins

# Constant birth rate, Mestel cooling
# dn/dMbol ~ t^{5/7} for constant SFR
M_bol = 4.75 - 2.5 * log_L_bins  # Absolute bolometric magnitude

# WDLF model
t_from_L = 1.0 * (L_bins / 1e-2)**(-5.0/7.0)  # Gyr
# For constant SFR and disk age T_disk:
T_disk = 9.0  # Gyr
wdlf = np.zeros_like(L_bins)
for i, t_age in enumerate(t_from_L):
    if t_age < T_disk:
        wdlf[i] = t_age**(5.0/7.0) * 0.5  # arbitrary normalization
    else:
        wdlf[i] = 0  # no WDs this faint yet

ax4.semilogy(M_bol, wdlf + 1e-10, color='#ffa500', linewidth=3, label='Model WDLF')

# Simulated observations with scatter
mask_obs = wdlf > 0.01
if np.any(mask_obs):
    noise = np.random.default_rng(42).normal(0, 0.15, size=np.sum(mask_obs))
    wdlf_obs = wdlf[mask_obs] * 10**(noise)
    ax4.scatter(M_bol[mask_obs], wdlf_obs, color='#ffdd00', s=25, alpha=0.7, zorder=5, label='Simulated data')

# Mark the cutoff
ax4.axvline(x=16, color='#ff4444', linestyle='--', linewidth=2, alpha=0.7, label='Disk age cutoff')

ax4.set_xlabel('Absolute Bolometric Magnitude', fontsize=12)
ax4.set_ylabel('Number Density (arb. units)', fontsize=12)
ax4.set_title('White Dwarf Luminosity Function', fontsize=14, fontweight='bold')
ax4.legend(fontsize=9, facecolor='#1a1a1a', edgecolor='#444', labelcolor='white')
ax4.invert_xaxis()
ax4.set_ylim(1e-2, 10)
ax4.grid(True, alpha=0.15)

# Print summary
print("=" * 65)
print("WHITE DWARF PHYSICS: KEY RESULTS")
print("=" * 65)

print("\nChandrasekhar Mass Limit:")
M_ch = 5.83 / mu_e**2
print(f"  M_Ch = 5.83/mu_e^2 = {M_ch:.2f} M_sun (mu_e = {mu_e})")
print(f"  With corrections: M_Ch ~ 1.44 M_sun")

print("\nNon-relativistic degeneracy constant:")
print(f"  K_NR = {K_NR:.4e} (SI units)")
print(f"  P = K_NR * rho^(5/3), polytropic index n = 3/2")

print("\nUltra-relativistic degeneracy constant:")
print(f"  K_UR = {K_UR:.4e} (SI units)")
print(f"  P = K_UR * rho^(4/3), polytropic index n = 3")

print("\nTransition density (x_F = 1):")
print(f"  rho_trans = {rho_transition:.3e} kg/m^3")

print("\nMass-Radius relation (selected):")
print(f"{'M (Msun)':10} {'R (Rsun)':10} {'R (km)':10} {'rho_c (kg/m^3)':16}")
print("-" * 50)
sample_masses = [0.3, 0.5, 0.6, 0.8, 1.0, 1.2, 1.35]
for m_target in sample_masses:
    if len(M_wd[mask]) > 0:
        idx = np.argmin(np.abs(M_wd[mask] - m_target))
        r_val = R_wd[mask][idx]
        print(f"{M_wd[mask][idx]:10.2f} {r_val:10.4f} {r_val*R_sun/1e3:10.0f} {'--':16}")

print("\nMestel cooling law: L ~ t^{-7/5}, T_c ~ t^{-2/5}")
print("\nCooling times to L = 10^{-4} L_sun:")
for M_val in [0.5, 0.6, 0.8, 1.0]:
    t_cool_est = (1e-4 / (1e-2 * M_val**(5.0/7.0)))**(-5.0/7.0)
    print(f"  M = {M_val} M_sun: t ~ {t_cool_est:.1f} Gyr")

plt.tight_layout()
plt.savefig('output.png', dpi=150, bbox_inches='tight', facecolor='#0a0a0a')
print("\nPlot saved to output.png")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

Chapter Summary

White dwarfs are supported by electron degeneracy pressure, whose equation of state transitions from $P \propto \rho^{5/3}$ (non-relativistic) to$P \propto \rho^{4/3}$ (ultra-relativistic) at densities around $10^9$ kg m$^{-3}$. The softening of the equation of state in the relativistic regime leads to the Chandrasekhar mass limit $M_{\text{Ch}} \approx 1.44\,M_\odot$.

The mass-radius relation $R \propto M^{-1/3}$ in the non-relativistic limit means that more massive white dwarfs are smaller. The radius shrinks to zero as the mass approaches $M_{\text{Ch}}$.

White dwarfs cool according to Mestel's law $L \propto t^{-7/5}$, with crystallization providing an additional delay. The faint-end cutoff of the white dwarf luminosity function serves as an independent cosmic chronometer. White dwarfs approaching the Chandrasekhar limit in binary systems produce Type Ia supernovae, the standardizable candles that revealed the accelerating expansion of the Universe.

← HR Diagram Neutron Stars →