Этап 2: модуль Власова (частицы) для 3D-расчётов — проект

Дата: 30.09.2026. Контекст: этап 1 (план §6.12) показал, что скалярное поле как прокси пыли непригодно для режима малых μ: предел Каупа M_max = 0.633/m превышает массу первого горизонта 1D-карты при любой доступной m, а стоимость ∝ m (осцилляции поля задают шаг). Физически правильный источник — бесстолкновительные частицы (уравнение Власова), как в 1D-коде V1. Этот документ фиксирует выбор платформы, уравнения, объём работ и тесты.

1. Платформа #

Выбор: GRTeclyn, ветка AHFinder (или main после слияния). Риск: ветка не доведена — тогда AH для e = 0 считать по сферическим срезам самостоятельно (поиск R = 2M по лучам), для e ≠ 0 — ждать слияния.

2. Уравнения #

Частица p с массой m_p (комовингный вес), положение x^i_p, ковариантный импульс u_i (пространственная часть 4-скорости, нижние индексы), Γ = √(1 + γ^{ij}u_i u_j) (γ^{ij} = χ h̃^{ij}). Уравнения геодезических в 3+1 (Vincent+2012, как в 1D-коде V1, обобщение на 3D):

• dx^i/dt = α γ^{ij} u_j / Γ − β^i
• du_i/dt = −Γ ∂_i α + u_j ∂_i β^j − (α/2Γ) u_j u_k ∂_i γ^{jk}

Источники в правой части CCZ4 (нормальный наблюдатель), депозит на сетку с весом 1/(√γ dV) (√γ = χ^{−3/2}):

• ρ = Σ m_p Γ_p W(x − x_p) / (√γ dV)
• S_i = Σ m_p u_{i,p} W / (√γ dV)
• S_{ij} = Σ m_p u_{i,p} u_{j,p} W / (Γ_p √γ dV), S = γ^{ij} S_ij

Ядро W: CIC/TSC (AMReX ParticleToMesh с пользовательской функцией депозита); на границах уровней — депозит на своём уровне + усреднение вниз (AMReX average_down) и заполнение призрачных ячеек, как для полей.

Класс материи ParticleMatter для CCZ4RHSWithMatter: e-m тензор читается из сеточных переменных (ρ, S_i, S_ij — 10 компонент, «эволюционные» с нулевой правой частью, чтобы интерполяция по уровням и рестарт работали как для полей); частицы — отдельный ParticleContainer с толчком после каждого шага RK4 (или на каждой стадии — как в V1: KDK внутри стадий; для начала — толчок на шаге с полями, взятыми в середине шага, что даёт 2-й порядок по времени для частиц; поля — 4-й).

3. Объём работ (оценка) #

1. Сборка GRTeclyn (AHFinder) + AMReX на сервере, прогон примеров ScalarField/BinaryBH — 1 день.
2. ParticleMatter: депозит ρ, S_i, S_ij; e-m тензор; связи — 2 дня.
3. Толчок частиц: интерполяция α, β, χ, h̃, их производных в точках частиц (расширение ParticleInterpolator), интегратор, перераспределение по боксам/уровням после толчка — 3 дня.
4. Начальные данные: FLRW-пыль (равномерные частицы, u_i = 0), LTB при t₀ (генератор ltb_to_grchombo.py → частицы с весами m'(r)dr и u_i = 0 в синхронном срезе: связи точны, как для скаляра) — 1 день.
5. Тесты (см. §4) — 3 дня. Итого ≈ 2 недели чистой работы.

4. Тесты #

1. FLRW-пыль: ρ ∝ a⁻³, K = −3H, невязка связей от дискретности ∝ 1/√N_p.
2. LTB e = 0, μ = 0.3 и 0.1 (холодная пыль): ρ(0,t), R(t) оболочек, первый AH и кривая M_AH(t) против 1D-кода (таблица в devlog от 30.09) — цель: совпадение времени достижения 0.1 и 0.3 M_H в пределах 5 %.
3. Тёплая пыль σ_ent: та же карта против 1D (σ = 0.016, 0.054).
4. Затем e ≠ 0 (§6.12, этап 1b): e = 0.13, 0.19 при μ = 0.04–0.06 из ICs Зельдовича/PM-3D при t₀ (связи — GRTresna по депонированным ρ, S_i).

5. Что уже есть и переиспользуется #

• Генератор LTB-среза в изотропных координатах (ltb_to_grchombo.py), калибровка K − f(x)K_far (§ devlog 29–30.09), вложенные сферы AMR с порогом плотности, анализ M_AH(t) (pbh_analyse.py).
• 1D-код V1 как эталон (карты холодная/тёплая, runs/v1/ref3d).

6. Реализация: раскладка кода (GRTeclyn, ветка AHFinder fd76177 + AMReX development) #

Каталог pbhgr/grteclyn/PBHVlasov/ (копируется в GRTeclyn/Examples/PBHVlasov/):

Порядок: (1) скелет с нулём частиц компилируется и даёт вакуум/FLRW при ручном ρ; (2) депозит + FLRW-пыль (u_i = 0, равномерные частицы) — тест a(t), K(t); (3) толчок — тест: свободная частица в FLRW (u_i = const, dx/dt = α γ^{ij}u_j/Γ), затем LTB e = 0 при μ = 0.3 против 1D (таблица от 30.09). Депозит только на уровне частицы; на грубые уровни — average_down; сглаживание ядра: CIC на всех уровнях, при необходимости TSC. Число частиц: ≥ 8 на ячейку самого мелкого уровня в ядре (в 1D нужны были 32 оболочки на ячейку для сходимости релаксации пересечений — здесь 3D-аналог проверяется сходимостью по N_p).