---
title: "Термодинамика расширяющейся Вселенной"
topic: thermal-history-eq
order: 4
description: "Равновесные распределения в метрике ФЛРУ и их сохранение при свободном разлёте; плотности n, ρ, s и эффективные числа степеней свободы g_*, g_*s Стандартной модели; сохранение энтропии в сопутствующем объёме; связь времени и температуры; температура реликтовых нейтрино; сегодняшние плотности фотонов, нейтрино, энтропии и η"
---
Ранняя Вселенная была почти идеальным газом в тепловом равновесии:
спектр реликтового излучения — планковский с точностью $10^{-4}$
[У]{.status}. В расширяющемся пространстве равновесие не мешает
эволюции, а задаёт её: при медленном (по сравнению со скоростями
реакций) расширении состояние в каждый момент определяется одной
температурой $T$, а её падение — сохранением энтропии. Отклонения от
равновесия (отрыв частиц, нуклеосинтез, рекомбинация) — тема
следующего уровня; здесь — равновесный фон, на котором они происходят.
Единицы натуральные: $\hbar = c = k_B = 1$.
## Равновесные распределения в расширяющемся пространстве
В равновесии числа заполнения — распределения Бозе–Эйнштейна и
Ферми–Дирака:
$$
f(p) = \frac{1}{e^{(E - \mu)/T} \mp 1}, \qquad E = \sqrt{p^2 + m^2},
$$
и плотности получаются интегрированием по импульсам с весом $g$ —
числом внутренних состояний.
::: {#def-thermo-densities}
Для сорта частиц с $g$ внутренними состояниями
$$
n = g\int\frac{\dd^3p}{(2\pi)^3}f,\quad
\rho = g\int\frac{\dd^3p}{(2\pi)^3}Ef,\quad
p = g\int\frac{\dd^3p}{(2\pi)^3}\frac{p^2}{3E}f,\quad
s = \frac{\rho + p - \mu n}{T}.
$$
:::
**Свободный разлёт сохраняет форму распределения безмассовых частиц.**
После того как частицы перестали взаимодействовать, импульс каждой
убывает как $1/a$ (тема «Метрика Фридмана–Леметра–Робертсона–Уокера»),
а число частиц в сопутствующем объёме сохраняется. Поэтому
$f(p, t) = f_{\rm dec}(p\,a/a_{\rm dec})$. Для $m = 0$ и $\mu = 0$ это
$1/(e^{p a/(a_{\rm dec}T_{\rm dec})} \mp 1)$ — снова равновесное
распределение с температурой
$$
T(t) = T_{\rm dec}\,\frac{a_{\rm dec}}{a(t)} .
$$
Поэтому реликтовые фотоны сохраняют планковский спектр через
13,8 млрд лет после отрыва [У]{.status}. Для нерелятивистских частиц
$E \approx m + p^2/2m$, и форма $e^{-p^2/2mT}$ сохраняется с
$T \propto a^{-2}$: оторвавшийся холодный газ остывает быстрее
излучения.
**Пределы.** Для ультрарелятивистских частиц ($T \gg m$, $\mu = 0$)
интегралы берутся точно. Разложим $1/(e^x \mp 1) = \sum_{k\ge1}(\pm1)^{k+1}e^{-kx}$
и проинтегрируем почленно — SymPy суммирует ряды в дзета-функции:
```{python}
import math
import astropy.constants as const
import astropy.units as u
import matplotlib.pyplot as plt
import numpy as np
import sympy as sp
from astropy.cosmology import Planck18
from scipy.integrate import cumulative_trapezoid
from scipy.special import zeta
from cosmevo.thermo import (BOSON, FERMION, equilibrium, g_star,
neutrino_temperature_ratio)
from cosmevo.units import natural
x = sp.Symbol("x", positive=True)
k = sp.Symbol("k", integer=True, positive=True)
def bose_fermi_integrals(n):
"""∫₀^∞ xⁿ/(eˣ ∓ 1) dx как суммы рядов по k."""
term = sp.integrate(x**n * sp.exp(-k * x), (x, 0, sp.oo)) # n!/k^{n+1}
bose = sp.summation(term, (k, 1, sp.oo))
fermi = sp.summation((-1) ** (k + 1) * term, (k, 1, sp.oo))
return bose, fermi
I_n_b, I_n_f = bose_fermi_integrals(2)
I_rho_b, I_rho_f = bose_fermi_integrals(3)
print("∫x³/(eˣ−1) =", I_rho_b, " ∫x²/(eˣ−1) =", I_n_b)
print("фермионы/бозоны: ρ —", sp.nsimplify(sp.simplify(I_rho_f / I_rho_b)),
", n —", sp.nsimplify(sp.simplify(I_n_f / I_n_b)))
```
С множителем $g/(2\pi^2)$ из угловой части:
$$
\rho = \frac{\pi^2}{30}gT^4,\quad
n = \frac{\zeta(3)}{\pi^2}gT^3,\quad
p = \frac{\rho}{3},\quad
s = \frac{2\pi^2}{45}gT^3
$$
для бозонов и те же формулы с множителями $7/8$ ($\rho$, $p$, $s$) и
$3/4$ ($n$) для фермионов. В нерелятивистском пределе ($T \ll m$)
$n = g(mT/2\pi)^{3/2}e^{-(m-\mu)/T}$: при $\mu = 0$ частицы
экспоненциально исчезают — пары аннигилируют, когда $T$ падает ниже $m$.
Численные интегралы для любой массы — функция `cosmevo.thermo.equilibrium`
(проверки предельных случаев — в тестах пакета):
```{python}
for y in (0.0, 1.0, 3.0, 10.0): # y = m/T
b = equilibrium(y, 1.0, 1, BOSON)
f = equilibrium(y, 1.0, 1, FERMION)
print(f"m/T = {y:4}: ρ_b/(π²T⁴/30) = {b['rho'] / (math.pi**2 / 30):.4f}, "
f"ρ_f/(π²T⁴/30) = {f['rho'] / (math.pi**2 / 30):.4f}, p_f/ρ_f = {f['p'] / f['rho']:.3f}")
```
## Эффективные числа степеней свободы
Полная плотность энергии и энтропии газа частиц разных сортов
записывается через температуру фотонов $T$:
::: {#def-gstar}
$$
\rho = \frac{\pi^2}{30}\,g_*(T)\,T^4, \qquad
s = \frac{2\pi^2}{45}\,g_{*s}(T)\,T^3,
$$
$$
g_* = \sum_{\rm бозоны}g_i\left(\frac{T_i}{T}\right)^4
+ \frac78\sum_{\rm фермионы}g_i\left(\frac{T_i}{T}\right)^4,
\qquad
g_{*s} = \sum_{\rm бозоны}g_i\left(\frac{T_i}{T}\right)^3
+ \frac78\sum_{\rm фермионы}g_i\left(\frac{T_i}{T}\right)^3
$$
в ультрарелятивистском пределе; $T_i$ — температура сорта $i$ (у
оторвавшихся частиц она может отличаться от $T$). Для массивных частиц
вклад вычисляется по точным интегралам.
:::
Выше масс всех частиц Стандартной модели: бозоны — фотон (2), глюоны
($8\times2$), $W^\pm$ ($2\times3$), $Z$ (3), бозон Хиггса (1) — итого 28;
фермионы — кварки ($6\times3\times2\times2 = 72$), заряженные лептоны
($3\times4 = 12$), нейтрино ($3\times2 = 6$) — итого 90. Отсюда
$g_* = 28 + \frac78\cdot90 = 106{,}75$. При понижении температуры частицы
с $m \gg T$ выбывают. Построим $g_*(T)$ модулем `cosmevo.thermo`:
```{python}
#| label: fig-gstar
#| fig-cap: "Эффективные числа степеней свободы Стандартной модели g_* (энергия) и g_*s (энтропия) как функции температуры фотонов. Расчёт идеального газа с точными интегралами для массивных частиц; переход кварк-глюонной плазмы в адронный газ (пионы) условно резкий при 150 МэВ — по расчётам КХД на решётке это плавный кроссовер. Ниже ~1 МэВ g_* и g_*s расходятся: нейтрино оторвались и холоднее фотонов."
T_grid = np.logspace(-5.3, 3.5, 160) # ГэВ
gs = np.array([g_star(T) for T in T_grid])
fig, ax = plt.subplots(figsize=(7.5, 4.2))
ax.loglog(T_grid, gs[:, 0], color="tab:blue", label=r"$g_*$")
ax.loglog(T_grid, gs[:, 1], "--", color="tab:orange", label=r"$g_{*s}$")
for T_mark, label in [(173 / 3, "t, H, Z, W"), (4.2 / 3, "b, c, τ"), (0.15, "КХД-переход"),
(0.106 / 3, "μ, π"), (0.000511 / 3, "e±")]:
ax.axvline(T_mark, color="0.75", lw=0.7)
ax.text(T_mark / 1.25, 2.6, label, rotation=90, fontsize=8, color="0.35", va="bottom")
for val in (106.75, 10.75, 3.91, 3.36):
ax.axhline(val, color="0.85", lw=0.6, ls=":")
ax.set(xlabel="T, ГэВ", ylabel="число степеней свободы", ylim=(2.5, 150),
xlim=(T_grid[0], T_grid[-1]))
ax.invert_xaxis()
ax.legend(loc="upper right", fontsize=9)
ax.grid(alpha=0.2, which="major")
fig.tight_layout()
for T in (1e3, 1.0, 0.01, 1e-5):
g_r, g_e = g_star(T)
print(f"T = {T:7.0e} ГэВ: g_* = {g_r:7.2f}, g_*s = {g_e:7.2f}")
```
Ось температур направлена справа налево — в сторону времени. Ступени
следуют за массами частиц с отставанием: частица выбывает, когда
$T \approx m/3$–$m/6$, а не при $T = m$. Плато $g_* = 10{,}75$
(фотоны, $e^\pm$ и три сорта нейтрино: $2 + \frac78\cdot10$) держится
примерно от 10 до 1 МэВ — в эпоху, предшествующую нуклеосинтезу.
## Сохранение энтропии
Из уравнения непрерывности $\dot\rho + 3H(\rho + p) = 0$ (тема
«Уравнения Фридмана») и термодинамического тождества $\dd p/\dd T = s =
(\rho + p)/T$ (при $\mu = 0$) следует сохранение энтропии в сопутствующем
объёме. Проверим символьно для произвольного уравнения состояния:
```{python}
t = sp.Symbol("t")
a = sp.Function("a")(t)
T = sp.Function("T")(t)
rho_f, p_f = sp.Function("rho"), sp.Function("p")
rho, p = rho_f(T), p_f(T)
S_comoving = a**3 * (rho + p) / T
dS = sp.diff(S_comoving, t)
# уравнение непрерывности: ρ'(T)·Ṫ = −3(ȧ/a)(ρ + p); тождество: p'(T) = (ρ + p)/T
Tdot = sp.solve(sp.Eq(sp.diff(rho, t), -3 * sp.diff(a, t) / a * (rho + p)), sp.diff(T, t))[0]
dp_dT = (rho + p) / T
dS = dS.subs(sp.Subs(sp.Derivative(p_f(sp.Symbol("xi")), sp.Symbol("xi")), sp.Symbol("xi"), T), dp_dT)
dS = dS.subs(sp.diff(p, t), dp_dT * sp.diff(T, t)).subs(sp.diff(T, t), Tdot)
print("d(s a³)/dt =", sp.simplify(dS))
```
$$
\frac{\dd}{\dd t}\left(s a^3\right) = 0
\quad\Longrightarrow\quad
g_{*s}(T)\,T^3a^3 = \text{const}, \qquad
T \propto g_{*s}^{-1/3}\,a^{-1} .
$$
Пока $g_{*s}$ постоянно, $T \propto 1/a$. Когда сорт частиц аннигилирует,
его энтропия переходит к оставшимся частицам плазмы, и температура
падает медленнее $1/a$ [У]{.status}. Энтропия в сопутствующем
объёме — главный «счётчик» ранней Вселенной: отношения $n_i/s$
сохраняются, пока число частиц сорта $i$ не меняется.
## Температура реликтовых нейтрино
Нейтрино взаимодействуют только слабо и отрываются от плазмы при
$T_{\rm dec} \approx 1$–2 МэВ (оценка — в теме «Кинетическое уравнение в
расширяющейся Вселенной»), *до* аннигиляции $e^\pm$ при $T \sim m_e/3$.
После отрыва $T_\nu \propto 1/a$, а энтропия $e^\pm$ достаётся только
фотонам. Сохранение энтропии фотонов и $e^\pm$:
$\left(2 + \frac78\cdot4\right)(Ta)^3_{\rm до} = 2\,(Ta)^3_{\rm после}$, и
так как $T_\nu a$ не меняется,
$$
\frac{T_\nu}{T} = \left(\frac{4}{11}\right)^{1/3} .
$$
```{python}
T0 = Planck18.Tcmb0.value # К
T_nu0 = (4 / 11) ** (1 / 3) * T0
print(f"T_ν/T = (4/11)^(1/3) = {(4 / 11) ** (1 / 3):.4f}; T_ν0 = {T_nu0:.3f} К")
r_model = neutrino_temperature_ratio(1e-8)
print(f"модель cosmevo.thermo (мгновенный отрыв при 2 МэВ, m_e учтена): T_ν/T = {r_model:.4f}")
```
Небольшое отличие модели от $(4/11)^{1/3}$: при 2 МэВ электроны ещё не
вполне безмассовы. Точный расчёт с учётом того, что отрыв не мгновенен и
нейтрино чуть-чуть подогреваются аннигиляцией, даёт эффективное число
сортов нейтрино $N_{\rm eff} = 3{,}044$ вместо 3 (Фрусте, Питру, Вольпе,
2020; Беннетт и соавторы, 2021) [У]{.status}. Измерения по реликтовому
излучению согласуются с этим (Planck 2018: $N_{\rm eff} = 2{,}99 \pm 0{,}17$)
[У]{.status}. Реликтовые нейтрино напрямую не зарегистрированы, но их
вклад в расширение виден в реликтовом излучении и нуклеосинтезе.
## Время и температура
В эпоху излучения $H^2 = \frac{8\pi G}{3}\rho = \frac{8\pi^3}{90}g_*\frac{T^4}{\Mpl^2}$,
т. е. $H = 1{,}66\,g_*^{1/2}T^2/\Mpl$, и при постоянном $g_*$
($a \propto t^{1/2}$, $H = 1/2t$)
$$
t = \frac{1}{2H} = 0{,}301\,g_*^{-1/2}\frac{\Mpl}{T^2}
\approx 2{,}42\,g_*^{-1/2}\left(\frac{T}{1~\text{МэВ}}\right)^{-2}\ \text{с}.
$$
Когда $g_*$ меняется, формула верна лишь приближённо. Точная связь:
$\dd t = \dd a/(aH)$, где $a(T)$ задаётся сохранением энтропии,
$\dd\ln a = -\dd\ln T - \frac13\dd\ln g_{*s}$:
```{python}
M_pl = natural((const.hbar * const.c / const.G) ** 0.5).to(u.GeV).value
sec = natural(1 * u.s).to(1 / u.GeV).value # 1 с в ГэВ⁻¹
c_H = math.sqrt(8 * math.pi**3 / 90)
print(f"H = {c_H:.3f} g*^(1/2) T²/M_Pl; t = {1 / (2 * c_H):.3f} g*^(-1/2) M_Pl/T²"
f" = {1 / (2 * c_H) * M_pl / 1e-6 / sec:.3f} g*^(-1/2) (T/МэВ)⁻² с")
# точное t(T): интеграл по ln T от высокой температуры (там t = 1/(2H))
lnT = np.log(T_grid[::-1]) # от горячего к холодному
g_r, g_e = gs[::-1, 0], gs[::-1, 1]
H_T = c_H * np.sqrt(g_r) * T_grid[::-1] ** 2 / M_pl
dlna_dlnT = -1 - np.gradient(np.log(g_e), lnT) / 3
t_exact = 1 / (2 * H_T[0]) + cumulative_trapezoid(-dlna_dlnT / H_T, -lnT, initial=0)
def t_approx(T):
return 1 / (2 * c_H * math.sqrt(g_star(T)[0])) * M_pl / T**2
for T in (100.0, 1.0, 0.1, 1e-3, 1e-4):
te = np.exp(np.interp(-math.log(T), -lnT, np.log(t_exact)))
print(f"T = {T:7.0e} ГэВ: t = {te / sec:9.3e} с (точно), {t_approx(T) / sec:9.3e} с (формула с g*(T))")
```
Формула с локальным $g_*(T)$ ошибается сильнее всего сразу после
быстрого изменения $g_*$: на 11 % при $T = 0{,}1$ ГэВ (после
КХД-перехода) и на 2–3 % при $T = 0{,}1$ МэВ (во время аннигиляции
$e^\pm$); на плато — около процента. Опорные точки: $T = 1$ МэВ —
0,75 с, $T = 100$ ГэВ (электрослабый масштаб) — $2{,}4\cdot10^{-11}$ с
[У]{.status} (в той мере, в какой Стандартная модель
верна до этих энергий: прямо проверено коллайдерами при $T \lesssim$
ТэВ).
## Сегодняшние плотности
Из $T_0 = 2{,}7255$ К (Planck18):
```{python}
T0_GeV = natural(Planck18.Tcmb0).to(u.GeV).value
cm3 = natural(1 * u.cm).to(1 / u.GeV).value ** 3 # 1 см³ в ГэВ⁻³
n_gamma = 2 * zeta(3) / math.pi**2 * T0_GeV**3 * cm3
g_s0 = 2 + 7 / 8 * 6 * 4 / 11
g_r0 = 2 + 7 / 8 * 6 * (4 / 11) ** (4 / 3)
s0 = 2 * math.pi**2 / 45 * g_s0 * T0_GeV**3 * cm3
n_nu = 3 / 4 * 2 * zeta(3) / math.pi**2 * (4 / 11) * T0_GeV**3 * cm3 # ν + ν̄ одного сорта (g = 2)
rho_c = (3 * Planck18.H0**2 / (8 * math.pi * const.G)).to(u.g / u.cm**3)
n_b = (Planck18.Ob0 * rho_c / const.m_p).to(u.cm**-3).value
print(f"g_*0 = {g_r0:.3f}, g_*s0 = {g_s0:.3f}")
print(f"n_γ = {n_gamma:.1f} см⁻³; n_ν (на сорт, ν + ν̄) = {n_nu:.0f} см⁻³; s₀ = {s0:.0f} см⁻³")
print(f"n_b = {n_b:.3e} см⁻³; η = n_b/n_γ = {n_b / n_gamma:.2e}; s₀/n_b = {s0 / n_b:.2e}")
```
На каждый барион во Вселенной приходится около $1{,}6\cdot10^9$ фотонов
и $1{,}1\cdot10^{10}$ единиц энтропии [У]{.status}. Барион-фотонное
отношение $\eta = 6{,}1\cdot10^{-10}$ — главный параметр первичного
нуклеосинтеза, и его совпадение с измеренным по реликтовому излучению —
одна из самых сильных проверок горячей модели (темы «Первичный
нуклеосинтез» и «Горячая Вселенная — наблюдательные основания»).
Почему $\eta$ мало, но не равно нулю, — открытый вопрос (тема
«Барионная асимметрия Вселенной»). Почти вся энтропия Вселенной в
пределах горизонта — у реликтовых фотонов и нейтрино (если не считать
сверхмассивных чёрных дыр; тема «Что такое космическая эволюция»).
## Физика и наблюдательный статус
- Тепловое равновесие при $T \lesssim$ нескольких МэВ проверено
нуклеосинтезом и реликтовым излучением: спектр, $N_{\rm eff}$,
$\eta$ [У]{.status}.
- $g_*(T)$ при $T \gg$ МэВ — расчёт в рамках Стандартной модели,
проверенной коллайдерами при энергиях до ТэВ [У]{.status}; область
КХД-перехода требует расчётов на решётке (резкий переход в нашем
расчёте — грубое приближение).
- Если существуют неизвестные лёгкие частицы, $g_*$ при высоких $T$
больше; ограничение $\Delta N_{\rm eff}$ по реликтовому излучению
проверяет такие сценарии [Г]{.status} (упражнение 3).
## Итоги
- Свободный разлёт сохраняет равновесную форму спектра безмассовых
частиц с $T \propto 1/a$.
- $\rho = \frac{\pi^2}{30}g_*T^4$, $s = \frac{2\pi^2}{45}g_{*s}T^3$;
$g_* = 106{,}75$ выше масс всех частиц СМ, $10{,}75$ при 1–10 МэВ,
сегодня $g_* = 3{,}36$, $g_{*s} = 3{,}91$.
- $sa^3 = \text{const}$, $T \propto g_{*s}^{-1/3}a^{-1}$;
$T_\nu/T = (4/11)^{1/3}$, $T_{\nu 0} = 1{,}945$ К.
- $t \approx 2{,}42\,g_*^{-1/2}(T/\text{МэВ})^{-2}$ с; при 1 МэВ — 0,75 с.
- Сегодня $n_\gamma = 411$ см⁻³, $s_0 = 2891$ см⁻³,
$\eta = 6{,}1\cdot10^{-10}$.
## Упражнения
1. Посчитайте $g_*$ при $T \approx 1$ ГэВ, считая кварки $u, d, s$, глюоны,
$e$, $\mu$, нейтрино и фотон безмассовыми, а остальные частицы —
выбывшими. Сравните с расчётом `g_star(1.0)` и объясните разницу.
2. Найдите плотность энергии реликтовых нейтрино сегодня в предположении,
что они безмассовы, и долю $\Omega_\nu h^2$. Сравните с $\Omega_\gamma h^2$.
3. Гипотетическая безмассовая скалярная частица ($g = 1$) оторвалась от
плазмы при $T \gg 100$ ГэВ. Какова её температура сегодня в единицах
$T_\nu$ и вклад в $N_{\rm eff}$
($\Delta N_{\rm eff} = \frac87\rho_X/\rho_{\nu,1}$, где $\rho_{\nu,1}$ — плотность
одного сорта нейтрино)?
4. Во сколько раз выросла бы температура фотонов относительно закона
$T \propto 1/a$ при аннигиляции $e^\pm$? При КХД-переходе (по нашей
модели, от $g_{*s}$ сразу выше 150 МэВ к сразу ниже)?
5. Каким было бы $N_{\rm eff}$, если бы нейтрино оторвались не до, а
*после* аннигиляции $e^\pm$ (при $T \ll m_e$)?
::: {.callout-tip collapse="true"}
## Ответы
1. Бозоны: $2 + 16 = 18$; фермионы: $3\cdot12 + 4 + 4 + 6 = 50$;
$g_* = 18 + \frac78\cdot50 = 61{,}75$. Расчёт даёт 76,4: при 1 ГэВ
$c$-кварк ($m = 1{,}27$ ГэВ) и $\tau$-лептон ($1{,}78$ ГэВ) ещё дают
заметный вклад.
2. $\rho_\nu = 3\cdot\frac78\cdot\frac{4}{11}^{4/3}\rho_\gamma = 0{,}681\rho_\gamma$;
$\Omega_\gamma h^2 = 2{,}47\cdot10^{-5}$, $\Omega_\nu h^2 = 1{,}68\cdot10^{-5}$
(для безмассовых нейтрино; массивные сегодня нерелятивистские и дают
$\sum m_\nu/(93{,}1~\text{эВ})$ по стандартной формуле).
3. Энтропия, отданная после отрыва $X$, сжимает его относительно нейтрино:
$T_X/T_\nu = (10{,}75/106{,}75)^{1/3} = 0{,}465$;
$\Delta N_{\rm eff} = \frac47\left(\frac{10{,}75}{106{,}75}\right)^{4/3} = 0{,}027$.
Это на порядок ниже нынешней точности ($\pm0{,}17$); достичь такой
чувствительности — цель будущих наблюдений реликтового излучения.
4. $e^\pm$: $(11/4)^{1/3} = 1{,}401$. КХД-переход: в нашей модели
$g_{*s}$ падает с 61,3 до 16,7, и $(g_{*s}^{+}/g_{*s}^{-})^{1/3} = 1{,}54$
(в реальном кроссовере — плавно, но с тем же итоговым множителем).
5. Тогда нейтрино нагрелись бы вместе с фотонами, $T_\nu = T$, и их
плотность энергии была бы больше в $(11/4)^{4/3}$ раз:
$N_{\rm eff} = 3\,(11/4)^{4/3} = 11{,}6$ — это исключено
реликтовым излучением ($2{,}99 \pm 0{,}17$).
:::
## Литература
- Горбунов Д. С., Рубаков В. А. *Введение в теорию ранней Вселенной:
Теория горячего Большого взрыва.* — М.: ЛКИ, 2008. Гл. 3, 5.
- Kolb E. W., Turner M. S. *The Early Universe.* — Redwood City:
Addison-Wesley, 1990. Гл. 3.
- Husdal L. On effective degrees of freedom in the early Universe.
*Galaxies* **4**, 78 (2016), arXiv:1609.04979.
- Borsanyi S. et al. Calculation of the axion mass based on
high-temperature lattice quantum chromodynamics. *Nature* **539**, 69
(2016).
- Froustey J., Pitrou C., Volpe M. C. Neutrino decoupling including
flavour oscillations and primordial nucleosynthesis. *JCAP* **12**, 015
(2020), arXiv:2008.01074.
- Bennett J. J. et al. Towards a precision calculation of $N_{\rm eff}$
in the Standard Model II: neutrino decoupling in the presence of
flavour oscillations and finite-temperature QED. *JCAP* **04**, 073
(2021), arXiv:2012.02726.
- Planck Collaboration. Planck 2018 results. VI. Cosmological parameters.
*Astron. Astrophys.* **641**, A6 (2020), arXiv:1807.06209.
```{python}
#| echo: false
#| output: false
# Проверка формул и ответов (не показывается на странице)
assert I_rho_b == sp.pi**4 / 15 and I_n_b == 2 * sp.zeta(3)
assert sp.simplify(I_n_f.rewrite(sp.zeta) - sp.Rational(3, 2) * sp.zeta(3)) == 0
assert sp.simplify((I_rho_f / I_rho_b).rewrite(sp.zeta) - sp.Rational(7, 8)) == 0
assert sp.simplify((I_n_f / I_n_b).rewrite(sp.zeta) - sp.Rational(3, 4)) == 0
assert sp.simplify(dS) == 0
assert 28 + 7 / 8 * 90 == 106.75 and 2 + 7 / 8 * 10 == 10.75
assert abs(g_star(1e3)[0] - 106.75) < 0.1 and abs(g_star(0.01)[0] - 10.75) < 0.05
assert round(T_nu0, 3) == 1.945
assert round(c_H, 2) == 1.66 and round(1 / (2 * c_H), 3) == 0.301
assert round(1 / (2 * c_H) * M_pl / 1e-6 / sec, 2) == 2.42
assert round(n_gamma) == 411 and round(s0) == 2891 and round(n_nu) == 112
assert round(g_r0, 2) == 3.36 and round(g_s0, 2) == 3.91
assert round(n_b / n_gamma * 1e10, 1) == 6.1
assert round(n_gamma / n_b / 1e9, 1) == 1.6 and round(s0 / n_b / 1e10, 1) == 1.1
t1MeV = np.exp(np.interp(-math.log(1e-3), -lnT, np.log(t_exact))) / sec
assert round(t1MeV, 2) == 0.75
t100 = np.exp(np.interp(-math.log(100), -lnT, np.log(t_exact))) / sec
assert round(t100 * 1e11, 1) == 2.4
t01 = np.exp(np.interp(-math.log(0.1), -lnT, np.log(t_exact)))
assert round((t01 / t_approx(0.1) - 1) * 100) == 12 or round((t01 / t_approx(0.1) - 1) * 100) == 11
t1e4 = np.exp(np.interp(-math.log(1e-4), -lnT, np.log(t_exact)))
assert 0.02 <= t1e4 / t_approx(1e-4) - 1 <= 0.03
# 1
assert 18 + 7 / 8 * 50 == 61.75 and round(g_star(1.0)[0], 1) == 76.4
# 2
ratio_nu = 3 * 7 / 8 * (4 / 11) ** (4 / 3)
Og_h2 = Planck18.Ogamma0 * Planck18.h**2
assert round(ratio_nu, 3) == 0.681 and round(Og_h2 * 1e5, 2) == 2.47
assert round(ratio_nu * Og_h2 * 1e5, 2) == 1.68
# 3
assert round((10.75 / 106.75) ** (1 / 3), 3) == 0.465
assert round(4 / 7 * (10.75 / 106.75) ** (4 / 3), 3) == 0.027
# 4
assert round((11 / 4) ** (1 / 3), 3) == 1.401
g_above, g_below = g_star(0.1501)[1], g_star(0.1499)[1]
assert round(g_above, 1) == 61.3 and round(g_below, 1) == 16.7
assert round((g_above / g_below) ** (1 / 3), 2) == 1.54
# 5
assert round(3 * (11 / 4) ** (4 / 3), 1) == 11.6
```
