Считаю лучи…

Чёрная дыра EN

Чёрная дыра

Это не рисунок. Для каждого пикселя компьютер прямо сейчас решает уравнения Эйнштейна: откуда пришёл свет, огибая вращающуюся дыру, и каким он стал в пути.

Тяните, чтобы облететь. Ниже — статья с формулами.

Как это устроено

Картинка наверху — не рисунок и не видеозапись. Для каждого пикселя браузер решает уравнения движения света в пространстве-времени вращающейся чёрной дыры, а цвет считает по законам излучения. Ниже — какие формулы работают, что упрощено и как проверено, что всё сходится.

Пространство-время Керра

Вращающуюся незаряженную чёрную дыру описывает решение Керра (1963). В единицах, где гравитационная постоянная, скорость света и масса дыры равны единице ($G = c = M = 1$), у неё один параметр — спин $a = J/M$, от 0 до 1. Все расстояния дальше измеряются в гравитационных радиусах $GM/c^2$: для Стрельца A* это 6,3 млн км, для дыры в 10 масс Солнца — 15 км.

Мы считаем в декартовых координатах Керра–Шильда: метрика в них — плоская плюс «поправка» вдоль светоподобного вектора $l_\mu$:

$$g_{\mu\nu} = \eta_{\mu\nu} + f\, l_\mu l_\nu,\qquad f = \frac{2 r^3}{r^4 + a^2 z^2},$$ $$l_\mu = \Big(1,\ \frac{r x + a y}{r^2 + a^2},\ \frac{r y - a x}{r^2 + a^2},\ \frac{z}{r}\Big),$$

где $r$ задан неявно: $r^4 - (x^2 + y^2 + z^2 - a^2)\, r^2 - a^2 z^2 = 0$. Главное достоинство этих координат — в них нет особенности на горизонте. Поэтому одна и та же программа ведёт свет и снаружи, и сквозь горизонт, и внутри, когда камера падает.

Горизонты — там, где $\Delta = r^2 - 2r + a^2$ обращается в ноль: $r_\pm = 1 \pm \sqrt{1 - a^2}$. Внешний $r_+$ — горизонт событий. Между ним и поверхностью $r = 1 + \sqrt{1 - a^2\cos^2\theta}$ лежит эргосфера: там нельзя стоять на месте, само пространство увлекает всё вокруг оси. Поэтому неподвижная камера не подпускается к дыре ближе эргосферы.

Как идёт свет

Луч света — изотропная геодезическая. Удобнее всего записать её как движение в гамильтоновой механике: с импульсом $p_\mu$ и гамильтонианом

$$H = \tfrac12\, g^{\mu\nu} p_\mu p_\nu = 0,\qquad \frac{dx^\mu}{d\lambda} = \frac{\partial H}{\partial p_\mu},\qquad \frac{dp_\mu}{d\lambda} = -\frac{\partial H}{\partial x^\mu}.$$

В форме Керра–Шильда обратная метрика записывается так же просто: $g^{\mu\nu} = \eta^{\mu\nu} - f\, l^\mu l^\nu$. Производные по координатам мы берём аналитически, без численного дифференцирования. Интегрирует видеокарта: метод Рунге–Кутты 4-го порядка с шагом, пропорциональным расстоянию до центра. В среднем лучу хватает 40–60 шагов.

Лучи идут назад во времени — от камеры к источнику. Направление пикселя — это направление в собственной системе отсчёта камеры (ортонормированная тетрада вокруг её 4-скорости). Поэтому аберрация и смещение частоты для движущейся камеры получаются сами. Луч заканчивается в одном из трёх мест:

  • на диске — пересёк экваториальную плоскость между внутренним и внешним краем; точка пересечения уточняется секущими внутри шага;
  • на небе — ушёл дальше 15 000 M; направление его движения и есть точка неба, откуда пришёл свет;
  • в дыре — оказался внутри самой внутренней круговой фотонной орбиты и движется внутрь. Оттуда свет не разворачивается, так что пиксель чёрный.

Вдоль луча сохраняются энергия $E = -p_t$, момент $L = p_\phi$ и константа Картера $Q$. Эталонный интегратор (Дорман–Принс с контролем ошибки, двойная точность) держит $H$ на нуле с точностью $10^{-10}$, а $Q$ — до $10^{-8}$ относительной ошибки.

Тень и фотонное кольцо

У дыры Керра есть сферические фотонные орбиты: свет на них вечно наматывается вокруг дыры. Орбита радиуса $r$ несёт константы (Бардин, 1973; Тео, 2003)

$$\xi = \frac{L}{E} = -\frac{r^3 - 3r^2 + a^2 r + a^2}{a\,(r - 1)},\qquad \eta = \frac{Q}{E^2} = -\frac{r^3\,(r^3 - 6r^2 + 9r - 4a^2)}{a^2\,(r - 1)^2},$$

а $r$ пробегает отрезок между прямой и обратной круговыми фотонными орбитами, $r_{\text{ph}} = 2\big(1 + \cos\big(\tfrac23 \arccos(\mp a)\big)\big)$. Лучи с этими константами, дошедшие до камеры, очерчивают край тени. Пунктир «Аналитическая тень» строится именно так: константы переводятся в импульс в координатах Бойера–Линдквиста, потом в Керра–Шильда, потом в систему камеры. Рендер при этом не используется.

Для невращающейся дыры тень — круг радиусом $\sqrt{27}\,M \approx 5{,}2\,M$ (горизонт — $2M$). При вращении тень сдвигается и со стороны, вращающейся к нам, сплющивается: лучи, летящие по вращению, могут подойти ближе. Тонкая яркая линия по краю тени — фотонное кольцо, наложение изображений диска, обернувшихся вокруг дыры один, два и больше раз.

Диск

Диск — модель Новикова–Торна (1973): тонкий, непрозрачный, газ на круговых кеплеровских орбитах с угловой скоростью $\Omega = 1/(r^{3/2} + a)$. Внутренний край — последняя устойчивая круговая орбита (ISCO, Бардин–Пресс–Тьюкольски, 1972):

$$r_{\text{ISCO}} = 3 + Z_2 - \sqrt{(3 - Z_1)(3 + Z_1 + 2 Z_2)},$$ $$Z_1 = 1 + \sqrt[3]{1 - a^2}\,\big(\sqrt[3]{1 + a} + \sqrt[3]{1 - a}\big),\quad Z_2 = \sqrt{3a^2 + Z_1^2}.$$

Это 6 M без вращения и 1,24 M при $a = 0{,}998$. Поток энергии с единицы площади — формула Пейджа–Торна (1974):

$$F(r) = -\frac{\dot M}{4\pi \sqrt{-g}}\,\frac{\Omega_{,r}}{(E - \Omega L)^2} \int_{r_{\text{ISCO}}}^{r} (E - \Omega L)\, L_{,r}\, dr.$$

Каждый участок диска светит как чёрное тело: $\sigma T^4 = F$. Поэтому у самого внутреннего края температура падает до нуля (там нет трения), максимум приходится немного дальше, а вдали $T \propto r^{-3/4}$. Абсолютная температура зависит от темпа аккреции и массы. Настоящие диски горячее — ультрафиолет и рентген, и в видимом свете они были бы однородно голубовато-белыми. Мы берём пик около 6500 K, чтобы цвета помещались в видимую часть спектра; его можно менять в параметрах.

Цвет и яркость

Свет, пришедший с диска, сдвинут по частоте в $g$ раз:

$$g = \frac{\nu_{\text{набл}}}{\nu_{\text{изл}}} = \frac{-p_\mu u^\mu_{\text{камеры}}}{-p_\mu u^\mu_{\text{газа}}} = \frac{1}{u^t\,(-p_t - \Omega\, p_\phi)},$$

если импульс нормирован на единичную энергию в системе камеры. В $g$ сразу входит всё: движение газа к нам или от нас (Доплер), замедление его часов на орбите и подъём света из гравитационной ямы.

Дальше помогает теорема Лиувилля: $I_\nu/\nu^3$ не меняется вдоль луча. Значит, чёрное тело температуры $T$, увиденное со смещением $g$, — это в точности чёрное тело температуры $gT$. Цвет берётся из таблицы «температура → цвет». Её мы считаем, интегрируя закон Планка с функциями цветового соответствия CIE 1931 (в аппроксимации Уаймана, Слоана и Ширли, 2013). Яркость, растущая как $g^4$, выходит отсюда сама.

Выключатели в параметрах разделяют эффекты. «Без Доплера» — газ заменяется наблюдателем с нулевым моментом импульса (ZAMO): он покоится относительно местного пространства, и вокруг дыры его несёт только увлечение систем отсчёта. «Без гравитационного смещения» — делим на смещение для такого наблюдателя. Если выключить оба, получится $g = 1$, как в «Интерстелларе».

Вдоль луча интегрируется и время $t$, поэтому узор диска берётся в момент излучения, а не в момент приёма. Свет с дальней стороны диска идёт дольше, и его изображение чуть «отстаёт».

Небо

Фон — настоящее небо Земли: Млечный Путь с панорамы NASA «Deep Star Maps 2020» и 9096 звёзд Йельского каталога ярких звёзд. Звёзды рисуются не текстурой, а точками. Каждый пиксель знает, какой кусок неба в него попадает — по соседним лучам, — и если в этот кусок попала звезда, её поток умножается на коэффициент усиления гравитационной линзы

$$\mu = \frac{\Omega_{\text{пикселя}}}{\Omega_{\text{его следа на небе}}}.$$

Поэтому звезда у края тени вспыхивает, раздваивается, а её изображения движутся навстречу друг другу. Кратные изображения никто специально не рисует: разные пиксели просто получают свет от одной и той же звезды. Цвет звезды — чёрное тело её эффективной температуры, тоже со смещением $g$.

Падение

Кнопка «Упасть» отпускает камеру с небольшим боковым толчком: его подбирают, заранее просчитав всю мировую линию, — самый сильный из тех, при которых камера всё же падает в дыру и не проходит сквозь диск. Её 4-скорость подчиняется тем же уравнениям Гамильтона, только $H = -\tfrac12$. Интегрирует Дорман–Принс в двойной точности, шаг — по собственному времени камеры $\tau$. Показ ускорен вдали и замедлен у дыры, чтобы на каждую «октаву» расстояния уходило примерно одинаковое время. На панели видны настоящие $r$ и $\tau$.

Горизонт падающий наблюдатель не замечает: локально там нет ничего особенного. Но после него радиус становится временем — он убывает так же неизбежно, как тикают часы. У невращающейся дыры путь заканчивается в сингулярности. От горизонта до неё — не больше $\pi M$ собственного времени: около минуты для Стрельца A*, 0,15 мс для дыры в 10 масс Солнца. У вращающейся раньше встречается внутренний горизонт $r_-$. Там излучение, падающее следом, неограниченно синеет, и сама классическая картина (решение Керра) перестаёт быть надёжной. На нём мы и останавливаемся.

Приливные силы у горизонта растут как $M/r^3$. Дыра звёздной массы растянет человека задолго до горизонта. Сверхмассивную можно пересечь, ничего не почувствовав.

Что упрощено

  • Диск бесконечно тонкий. Нет короны, джета и горячего толстого потока, которые на самом деле доминируют в M87* и Sgr A*: снимки EHT сделаны в радиодиапазоне, и там светится именно такой поток.
  • Сгустки в диске — иллюстрация: фрактальный шум, который вращается вместе с газом и растягивается дифференциальным вращением. Температурный профиль, скорости и смещения от него не зависят; выключатель «Турбулентность» показывает гладкий диск.
  • Температура диска подобрана так, чтобы цвета были видимыми, а небо подсвечено, чтобы звёзды было видно рядом с диском.
  • Нет излучения из области внутри ISCO и нет поляризации.
  • Цвет панорамы Млечного Пути при сильном смещении пересчитывается приближённо, как у чёрного тела в 6500 K; звёзды — точно.
  • Видеокарта считает в 32-битных числах. Лучи, которые наматываются у горизонта дольше 700 шагов, считаются поглощёнными — эталонный интегратор даёт для них тот же ответ.

Как проверено

Физика написана один раз на JavaScript в двойной точности; шейдер повторяет её строка в строку. Автотесты сравнивают её с известными ответами:

ПроверкаОжидаетсяПолучено
Метрика Керра–Шильда, пересчитанная в координаты Бойера–Линдквистаформа Бойера–Линдквиста±10⁻¹⁰
Аналитический градиент гамильтониана против конечных разностейсовпадение±10⁻⁷
ISCO при a = 0 и a = 0,9986 M; 1,2369 Mточно
Критический прицельный параметр Шварцшильда√27 M = 5,19615 M±2·10⁻⁶
Отклонение света далеко от дыры4M/b±1%
Край тени при a = 0,5; 0,9; 0,998 против кривой Бардинасовпадение±2·10⁻⁴ ширины
Смещение частоты диска Шварцшильда, видимого сверху√(1 − 3M/r)±10⁻⁷
Поток Пейджа–Торна: интеграл против замкнутой формулысовпадение±2·10⁻⁴
Цвет чёрного тела при 3000, 6500, 10 000 K против локуса Планка CIE(x, y)±0,003
Время свободного падения до горизонта против циклоидысовпадение±0,01 M
Быстрый RK4 (как в шейдере) против эталонатот же кадр98,5% пикселей; r ±0,2%; небо ±0,05°
Карта лучей с видеокарты (32 бита, настоящий браузер) против эталона, 4 ракурса и спинатот же кадр≥ 99,9% пикселей; r ±0,4%; небо ±0,1°

Источники

  1. R. P. Kerr. Gravitational field of a spinning mass as an example of algebraically special metrics. Phys. Rev. Lett. 11, 237 (1963).
  2. J. M. Bardeen, W. H. Press, S. A. Teukolsky. Rotating black holes: locally nonrotating frames, energy extraction, and scalar synchrotron radiation. ApJ 178, 347 (1972).
  3. J. M. Bardeen. Timelike and null geodesics in the Kerr metric. In Black Holes (Les Astres Occlus), 1973.
  4. I. D. Novikov, K. S. Thorne. Astrophysics of black holes. In Black Holes (Les Astres Occlus), 1973.
  5. D. N. Page, K. S. Thorne. Disk-accretion onto a black hole. ApJ 191, 499 (1974).
  6. J.-P. Luminet. Image of a spherical black hole with thin accretion disk. A&A 75, 228 (1979).
  7. E. Teo. Spherical photon orbits around a Kerr black hole. Gen. Rel. Grav. 35, 1909 (2003).
  8. O. James, E. von Tunzelmann, P. Franklin, K. S. Thorne. Gravitational lensing by spinning black holes in astrophysics, and in the movie Interstellar. Class. Quantum Grav. 32, 065001 (2015).
  9. C. Wyman, P.-P. Sloan, P. Shirley. Simple analytic approximations to the CIE XYZ color matching functions. JCGT 2(2), 2013.
  10. Event Horizon Telescope Collaboration. First M87 Event Horizon Telescope results (2019); First Sagittarius A* results (2022).
  11. NASA Scientific Visualization Studio. Deep Star Maps 2020. Yale Bright Star Catalogue, 5th ed.