How crystallization slows white dwarf cooling

Density gives each shell a different melting temperature

Measure radius $r$ outward from the center and give shell $i$ mass $\Delta m_i$ and fixed density $\rho_i$ . The controlled material is pure carbon with charge number $Z=6$ and mass number $A=12$ . It is placed on a source density geometry; the original carbon–oxygen mixture is not being evolved.

The ion number density $n_i$ determines a characteristic inter-ion spacing $a_i$ . Compare the Coulomb interaction energy at that spacing with thermal energy. Their ratio is the dimensionless coupling parameter $\Gamma_i$ . In Gaussian cgs units, the charge $e$ appears without the SI factor involving vacuum permittivity.

For this one-component-plasma exercise we choose melting at $\Gamma_m=175$ . The resulting melting temperature increases as the cube root of density. A denser central shell therefore freezes first during nearly isothermal cooling. A C/O mixture has its own composition-dependent phase diagram; substituting this fixed threshold would change the physical model.

n_i=\frac{\rho_i}{A m_u},\qquad a_i=\left(\frac3{4\pi n_i}\right)^{1/3},\qquad\Gamma_i=\frac{Z^2e^2}{a_i k_BT}

The ion spacing has length units, so Coulomb energy divided by thermal energy is dimensionless.

T_{m,i}=\frac{Z^2e^2}{175a_i k_B}\propto\rho_i^{1/3}

This threshold belongs to the declared pure-carbon OCP approximation; a C/O separation calculation requires a mixture phase diagram.

Sources: Bauer 2023, Carbon-Oxygen Phase Separation in MESA White Dwarf Models, arXiv2303.10110v2, Eq1–8 and §§2–5 · Bauer 2023, fixed MESA inputs and starting models, Zenodo 7846751, CC BY 4.0

A phase change gives energy without requiring a temperature drop

The total ion count in this shell is $N_i=\Delta m_i/(A m_u)$ . Assign it an approximate latent energy $q_i=0.77N_i k_BT_{m,i}$ . This coefficient is a stated teaching approximation, not a universal exact value for every ionic mixture. Bauer discusses composition-dependent departures in §4.1, footnote4.

While this shell freezes, hold its temperature at $T_{m,i}$ . If the prescribed luminosity depends only on that temperature, it stays constant during the event. The time required is the latent energy divided by that luminosity. Partial freezing advances in proportion to energy emitted.

Use the preceding chapter’s $0.6M_\odot$ pure-carbon ion count and specify a melting temperature of $3\times10^6\,\mathrm K$ . The chosen latent reservoir is $1.91\times10^{46}\,\mathrm{erg}$ . At constant luminosity $10^{-3}L_\odot$ , its release takes about $1.58\times10^8\,\mathrm{yr}$ . This is a specified single-reservoir example, not a measured WD delay.

\Delta t_i=\frac{q_i}{L(T_{m,i})},\qquad f_i(t)=\frac{L(T_{m,i})(t-t_{i,0})}{q_i}

The shell’s solid fraction

f_i

runs from zero to one during its latent interval starting at

t_{i,0}

. No artificial temperature-width heating bump is needed.

Sources: Bauer 2023, Carbon-Oxygen Phase Separation in MESA White Dwarf Models, arXiv2303.10110v2, Eq1–8 and §§2–5

Keep the sensible and latent reservoirs separate

Let $U_{\rm sens}(T)$ contain sensible thermal energy only, with the liquid and solid sensible branches joined continuously at melting. The latent reservoir is still available in the fraction that remains liquid. Its energy is therefore proportional to $1-f_i$ .

This choice makes a simple energy check possible. Integrate the sensible heat capacity over the temperature decline, add the latent energy already released, and compare their sum with the light emitted. If instead an EOS already includes latent heat in its internal-energy transition, adding this explicit latent reservoir again would count it twice.

A front moving inside one constant-density shell follows the volume that has frozen. The radius does not advance linearly with the fraction of the shell’s mass. The cube of radius does, because spherical volume is proportional to that cube.

E(T,\{f_i\})=U_{\rm sens}(T)+\sum_i(1-f_i)q_i,\qquad E_\gamma=E(0)-E(t)

The fixed-volume controlled calculation has no additional contraction, nuclear, neutrino or composition-separation term.

r_f^3=r_{\rm in}^3+f_i(r_{\rm out}^3-r_{\rm in}^3)

The shell’s inner and outer radii bound a fixed-density cell. This relation determines the displayed partial-cell front from the same latent energy.

Sources: Bauer 2023, Carbon-Oxygen Phase Separation in MESA White Dwarf Models, arXiv2303.10110v2, Eq1–8 and §§2–5 · Bauer 2023, fixed MESA inputs and starting models, Zenodo 7846751, CC BY 4.0

A crystal’s heat capacity falls in the cold limit

The ions in a solid vibrate collectively. At high temperature many vibrational modes are excited and the heat capacity approaches the classical harmonic value. At low temperature, only the longest-wavelength modes can be excited efficiently, so the ionic heat capacity falls strongly.

The Debye temperature $\Theta_D$ sets the scale in a simplified vibrational spectrum. With $x=\Theta_D/T$ , the following integral gives the ionic heat capacity. The controlled model takes $\Theta_D=0.45\hbar\omega_p/k_B$ with ion plasma frequency $\omega_p$ , the characteristic angular frequency of collective ion motion in inverse seconds; that spectrum is a declared approximation.

At $T/\Theta_D=0.01$ , numerical integration gives $C_D/(N_i k_B)\simeq2.34\times10^{-4}$ , in agreement with the low-temperature limit. This is the ionic term, not the whole star’s heat capacity. Strongly degenerate electrons still contribute a term proportional to temperature. At sufficiently low temperature their contribution can exceed the ionic term.

For a fixed luminosity, a smaller total heat capacity means temperature falls faster. A real cooling history also changes luminosity, so heat capacity alone does not determine whether a star moves quickly through an observed brightness interval.

C_D=9N_i k_B\left(\frac{T}{\Theta_D}\right)^3\int_0^{\Theta_D/T}\frac{u^4e^{-u}}{(1-e^{-u})^2}\,du

The integration variable

u

is dimensionless. The high-temperature limit is

3N_i k_B

.

C_D\simeq\frac{12\pi^4}{5}N_i k_B\left(\frac{T}{\Theta_D}\right)^3\quad(T\ll\Theta_D),\qquad C_e=\frac{\pi^2N_e k_B^2T}{p_Fv_F}

The electron expression is the leading strongly degenerate result;

N_e

is electron number,

p_F

Fermi momentum and

v_F

Fermi velocity. Their product has energy units. Evaluate each thermal contribution per shell before summing it over the core.

Sources: Bauer 2023, Carbon-Oxygen Phase Separation in MESA White Dwarf Models, arXiv2303.10110v2, Eq1–8 and §§2–5 · Mestel 1952, On the theory of white dwarf stars, MNRAS112583–597, thermal §§2–4

Derive the luminosity law before turning heat into time

For a controlled envelope, assume a thin radiative ideal gas with fixed stellar mass $M$ , negligible surface pressure and opacity $\kappa=\kappa_0\rho T^{-7/2}$ . The opacity is a cross-section per unit mass. Hydrostatic balance and radiative diffusion supply the pressure–temperature relation below.

Use ideal-gas pressure $P=\rho k_BT/(\mu m_u)$ , where $\mu$ is the mean particle mass in atomic units. Integrating gives a pressure squared proportional to $T^{17/2}/L$ . To join to an approximately isothermal degenerate core, make the controlled matching choice that nonrelativistic Fermi energy equals thermal energy. This gives $\rho_b=C_\rho T^{3/2}$ .

Eliminating the boundary pressure then yields $L=\mathcal A T_c^{7/2}$ . The coefficient depends on the opacity assumption and mean particle masses. It is not fitted to a native cooling track. The current model fixes $\kappa_0=10^{24}\,\mathrm{cm^5\,g^{-2}\,K^{7/2}}$ , $\mu=0.6$ and electron molecular weight $\mu_e=2$ .

This derivation supplies the teaching boundary for the core; the atmosphere remains unsolved. With the four frozen density geometries it produces very long total clocks, about 119–454 Gyr over the specified intervals. They remain explicitly uncalibrated. Real age work requires the atmosphere, opacities, composition and structure to be validated together.

\frac{dP}{dr}=-\frac{GM\rho}{r^2},\qquad\frac{dP}{dT}=\frac{64\pi\sigma GM}{3\kappa L}T^3

Radius increases outward; both pressure and temperature decrease. Here sigma is the Stefan–Boltzmann constant, and radiation pressure is omitted.

P^2=\frac{256\pi\sigma GMk_B}{51\kappa_0\mu m_u L}T^{17/2},\qquad C_\rho=\frac{\mu_e m_u}{3\pi^2}\left(\frac{2m_ek_B}{\hbar^2}\right)^{3/2}

The matching electron mass is

m_e

;

\hbar

is the reduced Planck constant. The boundary density is a controlled nonrelativistic matching relation.

L(T_c)=\frac{256\pi\sigma GM\mu m_u}{51\kappa_0C_\rho^2 k_B}T_c^{7/2}\equiv\mathcal A T_c^{7/2}

The coefficient has the units needed to turn core temperature to the power three and a half into luminosity.

Sources: Mestel 1952, On the theory of white dwarf stars, MNRAS112583–597, thermal §§2–4 · MIST white dwarf models, arXiv2509.21717v1, cooling clocks, surface boundary and model comparisons · MESA 22.11.1 fixed const_def.f90, matching original source units

A lighter outer shell can cause a longer latent delay

Take two hypothetical shells with inner ion count twice the outer count. Give them melting temperatures $T_{m,1}=3\times10^6\,\mathrm K$ and $T_{m,2}=2\times10^6\,\mathrm K$ . Use the same approximate latent coefficient in each and the derived luminosity power law.

The outer shell contains fewer ions and releases less total latent energy. But it freezes at a lower temperature, where the prescribed star is fainter. Divide each reservoir by the luminosity at that shell’s freezing temperature.

The resulting outer latent interval is about $1.38$ times the inner interval. Less energy can take longer to leave when the outward power is much smaller. These are the latent parts of the time only; sensible cooling between the two events must still be added.

\frac{\Delta t_2}{\Delta t_1}=\frac{N_2}{N_1}\left(\frac{T_{m,2}}{T_{m,1}}\right)^{-5/2}=\frac12\left(\frac23\right)^{-5/2}\simeq1.38

Latent energy is proportional to ion count times melting temperature, while luminosity follows temperature to the power three and a half. The ratio is independently computed in the worked-example producer.

Sources: Mestel 1952, On the theory of white dwarf stars, MNRAS112583–597, thermal §§2–4 · Bauer 2023, Carbon-Oxygen Phase Separation in MESA White Dwarf Models, arXiv2303.10110v2, Eq1–8 and §§2–5

Use an analytic limit to check the cooling program

Remove latent heat and hold the sensible heat capacity constant. With $L=\mathcal A T^{7/2}$ , the cooling equation can be integrated exactly. It gives the familiar Mestel-type temperature and luminosity powers. A numerical solver should approach this limit when given the same assumptions.

The independent reference also integrates the Debye heat capacity and compares it with the internal-energy difference, tests the high and low temperature limits, and checks exact single-cell latent intervals. Passing those tests establishes the controlled calculation’s numerical consistency. Its age calibration is a separate physical question, already exposed by the very long clocks.

T(t)=T_0(1+t/\tau)^{-2/5},\qquad L(t)=L_0(1+t/\tau)^{-7/5},\qquad\tau=\frac{2C}{5\mathcal A T_0^{5/2}}

These expressions apply only to constant heat capacity, the chosen power-law boundary and no additional energy sources.