Brief verification and V0 / V1 / V1b implementation plan

Date: 23 September 2026. Brief under review: /root/PBH/preview.md ("PBH collapse: oscillating-scalar EKG → Vlasov limit — minimum decisive experiment", 22 Sep 2026). Existing code: /root/PBH/pbhgr (spherical Einstein–Vlasov, polar-areal gauge, Stage A step 1).

Method. Every claim the brief tags [FACT] was re-read in the source (full arXiv PDF or HTML fetched on 23 Sep 2026; verbatim quotes with section and equation numbers are in agent_reports/*.md). Every [DERIVED] number was recomputed independently (check_scalings.py, ltb_yoo_profile.py, ltb_collapse_times.py; outputs in *.out and inline below). Every [HYPOTHESIS] and [SPECULATION] is listed in Section 3 with the measurement that decides it; none is used as an input to the plan in Sections 5–7. A second independent audit of the same brief (pbh-memo.md, from pbh-memo.zip) was cross-checked on 23 September; see Section 4a.

Status labels used below: VERIFIED (statement is in the source as the brief reads it), CORRECTED (true in substance, wrong in a detail that matters), NOT IN SOURCE (the source does not say it), REFUTED, and RECOMPUTED (my own calculation agrees with the brief unless a number is given).

0. Bottom line #

1. The brief's factual base survives. All Yoo et al. 2026 numbers (profile, entry at k/aH = √6, CMC start at H_i = 5k, LTB 0.32/0.33, ĥ(0.2,0) ≈ 0.45, dust crash for μ ≤ 0.95 at e = 0.2, particle horizon at μ ≥ 0.050 and none at μ ≤ 0.045, 80³, flagged constraint violations) are in the paper. Milligan et al. do not state m/H. The de Jong 2022 footnote exists. The abundance coefficients (0.05556 σ⁵, 3.70 σ^{3/2}, 0.2055 σ^{13/2}, 19×) are as quoted.
2. Six details are wrong or missing and change what the plan must do:
   • Yoo et al. evolve with slicing (2.20), which adds a term (3/2)αH_k exp(−C²α_b²/a_k²) with C = 5, not the (2.18) condition the brief quotes. V0 must use (2.20).
   • Their horizon-formation times are not stated anywhere; the brief's "100–250 t_H" is a plot reading. My own reading of Fig. 7 and an exact LTB calculation give t_AH ≈ (0.735 ± 0.02) × t_C^{LTB}(r_m) for every μ from 0.05 to 0.30 (Section 2.5). This makes the deadline reading quantitative, but the same law predicts a horizon at ≈ 220 t_H for μ = 0.045, inside their 250 t_H window, and none was seen. A single-infall deadline therefore does not explain the 0.045 result; a later (second-passage) collapse can. V0 stays decisive, and the brief's revised even-odds credence is about right.
   • The brief's "exact kinetic seed" φ = 0, Π = √(2ρ) is exactly the initial-data choice of Milligan et al. 2025 (Sec. III.3.2: φ̃_b0 = χ̃_0 = 0, all information in the kinetic term). It is correct but not new, and it should be cited as theirs.
   • de Jong, Aurrekoetxea & Lim 2022 solved a radial ODE for ψ with K = −3H_0, not CTTK; only the 2023 paper uses CTTK. The brief's decaying-mode diagnosis of flat-slice data is right (Section 2.8) and is echoed by de Jong's thesis (Sec. 6.2).
   • Two existing full-GR collisionless results bear directly on hypotheses H4/H7/F6 and are absent from the brief: East 2019 (3D, prolate spheroids of collisionless particles all form black holes, contra Shapiro–Teukolsky) and Ames, Andréasson & Rinne 2023 (axisymmetric Einstein–Vlasov, apparent horizon in every prolate case, hoop ratio ≤ 1.12). Both are asymptotically flat, not cosmological, but they remove the Shapiro–Teukolsky argument the brief leans on.
   • Observational inputs: n_s = 0.9743 ± 0.0034 (P-ACT-LB); the De la Torre Luque erratum is PRD 112, 109904 (Nov 2025); Gottlieb et al. 2026 weaken the Esser limit through capture inefficiency, not through star survival under slow accretion; the poltergeist calculation is not applicable at P_ζ ~ 10⁻³ (linear theory fails at kη_R ≈ 18), so the "TENSION" verdict for sudden reheating rests on no valid calculation either way.
3. The causal bound μ_V,sph ≤ 0.33 is sound. My growing-mode LTB reproduces Yoo's mapping exactly (their eqs. 4.10 and 4.12 are my m(r) and E(r) with m_Yoo = 6m), passes the Hellaby–Lake no-shell-crossing test for all μ ≤ 1, and the domain-of-dependence step has a rigorous precedent in Andréasson, Kunze & Rein 2011 and Andréasson & Rein 2025 (Oppenheimer–Snyder-type collapse of a collisionless gas forms a trapped surface). My independent recomputation puts the boundary between μ = 0.325 and 0.330 (Section 2.4), and shows that the trapped region guaranteed at μ = 0.33 encloses 0.94 M_H(t_k).
4. The novelty claim holds: no spherically symmetric Einstein–Vlasov PBH threshold in an expanding matter-dominated background exists in the literature through September 2026 (arXiv and INSPIRE searches returned zero hits). The direct competitor is Yoo et al. themselves (3D, 80³).
5. Codes: Yoo's COSMOS is public but its released version has no particle module; CosmoGRaPH (public, MIT) has a particle component in full BSSN and is the best public candidate for V0; OllinSphere-BiB (public, no license file) already runs φ² scalar dark matter with a Gaussian perturbation on an expanding background with an apparent-horizon finder and is the natural second code for V1b. Milligan's Misner–Sharp code is not public.
6. Plan changes relative to the brief: V1 must treat the velocity dispersion as a physical axis, not a regulator, and report μ_th as a surface over (deadline, σ_v, profile); V1 must include a run-by-run test of the causal bound at μ = 0.36 and 0.40; V0 depends on obtaining a 3D GR particle code, and the cheapest route is to ask Yoo's group to extend their own runs while an independent capability is built on CosmoGRaPH.

1. Literature ledger #

1.1 Yoo, Escrivà, Harada & Kohri 2026 (arXiv:2609.14218v1, 13 Sep 2026) #

1.2 Scalar-field papers #

1.3 Abundance, Newtonian and correspondence papers #

1.4 Einstein–Vlasov numerics, rigorous results, LTB theory #

1.5 Observational and cosmological inputs #

1.6 Novelty and codes #

2. Derived claims recomputed #

2.1 Scalings (Sections 2, 4, 9, 11, 12 of the brief) #

All numbers below are from check_scalings.py.

2.2 The growing-mode LTB solution equals Yoo's #

From γ_ij = a²e^{2ζ}δ_ij with ζ = μ e^{−k²r²/6} and pure growing mode (t_B = 0): m(r) = (H_i²/2)(e^{ζ}r)³, E(r) = rζ′(1 + rζ′/2). Yoo's (4.10) is m_Yoo = 4Ψ⁶/(3t_i²) = 6m and their (4.12) is k̃r² = −2E, term by term. The central singularity time is t_C(0)/t_k = 0.589 e^{3μ} μ^{−3/2} (8.36 t_k at μ = 0.33, matching Yoo Fig. 1's ≈ 2.05×10³/k = 8.4 t_H). The linear-theory version 0.589 μ^{−3/2} (δ_c = 1.686 applied to the central δ = 2.4μ) omits the e^{3μ} factor; it is a 15% effect at μ = 0.05 and a factor 2.6 at μ = 0.32.

2.3 No shell crossing before the central singularity (Hellaby–Lake) #

For the growing-mode LTB of the Yoo profile, E < 0 everywhere (all shells bound), M′ > 0, t_B′ = 0, and d ln t_C/dr > 0 for all r > 0 at every μ from 0.05 to 1.0 (ltb_yoo_profile.out, block 1). Hellaby–Lake's elliptic conditions are therefore satisfied: the first singularity is the central shell-focusing event p, as the causal-bound argument assumes. Yoo et al. do not report this check. For the monochromatic (sinc) mean profile the shells beyond r = 4.49/k₀ are unbound (E changes sign there), but inside that radius, which contains r_m = 2.74/k₀ and the whole overdense core, d ln t_C/dr > 0 for μ = 0.1, 0.3 and 0.5, so the first singularity is again central (block 5 of the same output and the inner-region recheck). The nakedness boundary and the causal-bound radius r₁ for P2 and P3 still have to be computed before V1 runs them; the script takes a profile argument.

2.4 Independent recomputation of the 0.32/0.33 nakedness boundary #

Method (ltb_yoo_profile.py, block 3): the limiting outgoing radial null ray from the regular centre launched at t → t_C(0)⁻ (this is the Cauchy horizon of the central singularity, the "most past null geodesic" of Yoo Sec. IV-C); the ray is integrated with dt/dr = R′/√(1 + 2E) through the exact growing-mode LTB, and classified TRAPPED if it enters a region with 2m/R > 1 and Ṙ < 0 (a future trapped sphere, as opposed to the cosmological anti-trapped region where Ṙ > 0) before reaching the underdense region, or ESCAPED otherwise. Results, with the launch offset δ = 10⁻³, 10⁻⁵, 10⁻⁷ giving the same verdict in every case:

The boundary lies between μ = 0.325 and 0.330, which is Yoo's "globally naked at 0.32, PBH at 0.33" reproduced independently from the profile alone (no input from their paper beyond the profile definition). Two further facts come out of the same table. First, the mass whose trapping is guaranteed by the causal bound is not small near the threshold: at μ = 0.33 the outer trapped shells enclose 0.94 M_H(t_k), because in LTB a shell is trapped only when it has fallen to its own Schwarzschild radius, which encloses all the interior mass. This is the quantity H6 should be compared with in the spherical case; de Jong's 10⁻² M_H formation mass comes from a massless-shell trigger and is not the same situation. Second, the table fixes the pass criterion of the V1 causal-bound runs: at μ = 0.36 the Vlasov apparent horizon must exist at r ≥ 0.51 r_m by t = 10.1 t_k, and at μ = 0.40 at r ≥ 0.37 r_m by t = 8.6 t_k, up to discretization error; a horizon appearing later, or only inside r₁, fails the cold-limit consistency check and triggers inspection of the regularization, initial geometry, boundaries and the domain-of-dependence step before any physical reading (the brief's F3); a finite-dispersion family is not bound by this argument.

2.5 Deadline test: Yoo's particle horizon times against the exact LTB #

Horizon-formation times are not stated in Yoo et al.; I read the midpoint of the α₀ drop in their Fig. 7 (uncertainty ≈ ±10%). The exact LTB collapse time of the shell at r_m (the shell enclosing the mean density the brief's deadline formula uses) is t_C(r_m) = 2π m(r_m)/(−2E(r_m))^{3/2} (ltb_collapse_times.py).

Reading: over a factor 6 in μ and a factor 10 in time, the e = 0.2 particle horizon appears at (0.735 ± 0.02) × t_C^{LTB}(r_m), with no sign of a floor. That is the quantitative form of the brief's "deadline" concern. But the same law puts the μ = 0.045 horizon at ≈ 220 t_H, inside the 250 t_H window, where the α₀ dip at ≈ 200 t_H shows the bulk infall arriving on schedule and then not trapping. So the 0.045 result is not explained by a single-infall deadline. It can still be delayed collapse on a second passage (the outgoing central streams re-collapsing after ≈ 1.5–2 t_C), a resolution artefact, or a real threshold at 80³. V0 to 1000 t_H separates these; the predicted horizon times under the "0.735 law" are 219, 260, 314 and 392 t_H for μ = 0.045, 0.04, 0.035 and 0.03.

2.6 Constraint-violation onset equals the LTB central caustic #

Yoo's Fig. 10 (max-norm Hamiltonian violation) rises from 10⁻² to O(1) at ≈ 25–30, 40, 50, 65–70 and 75 t_H for μ = 0.10, 0.075, 0.06, 0.05 and 0.045. The exact LTB central singularity times are 25.1, 35.9, 48.0, 61.2 and 70.6 t_H. The violation is therefore created at the first central shell crossing, as their footnote 1 says, and grows thereafter; the mean-norm curves keep rising to ≈ 10⁻³. This gives V0 a clean diagnostic: the max-norm onset must stay at t_C(0) while its amplitude must fall with resolution and particle number.

2.7 The causal bound μ_V,sph ≤ 0.33: assessment #

The argument has three steps. (i) At μ ≥ 0.33 the earliest outgoing null ray from the central singularity p is trapped before reaching the underdense region (Yoo Sec. IV-C, verified; my recomputation in 2.4). (ii) Every shell beyond the trapping radius r₁ is already inside its own apparent horizon when that ray arrives: along the ray inside the trapped region R decreases while 2m(r) increases outward, so R < 2m(r₁) < 2m(r) for r > r₁. Hence the apparent horizon at r > r₁ lies outside J⁺(p). (iii) With no shell crossing before p (2.3), cold Vlasov data evolve exactly as dust outside J⁺(p), so the same trapped spheres exist in the Vlasov spacetime, and by Hayward's theorems they persist. The mechanism "outer matter is trapped before anything from the centre can reach it" is the one used rigorously by Andréasson, Kunze & Rein 2011, and Andréasson & Rein 2025 prove that collisionless approximations of Oppenheimer–Snyder data form trapped surfaces. I therefore treat the bound as established for the L_reg → 0 limit, with two residual assumptions to be tested in V1: the Vlasov regularization must not act before p (true as L_reg → 0), and the numerically found horizon must appear at r ≥ r₁ no later than t_AH^{LTB}(r₁). The bound also yields a guaranteed minimum horizon mass, m(r₁), at the moment the ray arrives (value in 2.4), which is a direct check of H6.

2.8 Flat-slice data is the pure decaying mode (Section 4 of the brief) #

With ρ = ρ̄(1 + δ) and local expansion θ = 3H(1 + δ_H), linear continuity gives δ̇ = −3Hδ_H. The matter-era growing mode δ ∝ a has δ_H = −δ/3; the decaying mode δ ∝ a^{−3/2} has δ_H = +δ/2. A spatially flat physical slice γ_ij = δ_ij with K fixed locally by the Hamiltonian constraint (K² = 24πρ_loc) has exactly δ_H = +δ/2, i.e. zero growing-mode content; the same decomposition holds sub-horizon. A conformally flat metric whose factor ψ = e^{ζ/2} carries the curvature perturbation is the growing mode, so the diagnosis concerns the flat physical slice, not CTTK as such. The diagnosis is correct and is consistent with de Jong's thesis remark that perturbing the background field "can have large effects on the local expansion rate".

2.9 The kinetic seed and its generalization #

φ = 0, Π = √(2ρ) on the comoving dust slice satisfies both constraints exactly (E = Π²/2 = ρ, J_i = −Π∂_iφ = 0). This is Milligan et al.'s initial condition. The generalization used in the V1b plan keeps a homogeneous φ̄ and sets Π(x) = √(2[E(x) − V(φ̄)]) with E from the Hamiltonian constraint of the separate-universe metric (γ_ij = a²e^{2ζ}δ_ij, K = −3H_i, Ã_ij = 0): both constraints hold exactly for any oscillation phase, which allows the phase-variation test of gate 5. Existence requires E ≥ V(φ̄) everywhere, which holds for the amplitudes in question; Baumgarte, Clough & Giblin 2025 show that CTTK-type constructions can fail to have solutions, so the algebraic construction here is preferable.

2.10 Is the cold limit (L_reg → 0) regular? #

The brief assumes convergence in L_reg → 0. Newtonian self-similar infall gives ρ ∝ r⁻² inside the virialized region for flat-topped initial profiles (δ̄ ∝ M^{−ε} with ε ≤ 2/3; Fillmore–Goldreich 1984, confirmed in Sikivie–Tkachev–Wang 1997 and Nusser 2001), so 2m/R tends to a constant of order the virial compaction rather than diverging; only steeper profiles (ε > 2/3) give r^{−9ε/(1+3ε)} cusps with 2m/R growing inward. The Yoo profile and the peak-theory mean profiles are flat-topped, so a regular cold limit is plausible, and a horizon then forms when the plateau value of 2m/R reaches 1. Since the plateau scales with the compaction at turnaround, C_ta ≈ 0.88 μ (from shell energy conservation, my estimate), a spherical collisionless threshold of order 0.2–0.4 is plausible but not established. Angular momentum does not create a constant-density core in these solutions (Sikivie et al.), so σ_v is a genuine physical axis. V1 measures the plateau and its L_reg and σ_v dependence directly; the estimate here is not an input.

3. Working hypotheses of the brief, and what decides each #

None of the following is established; each is retained as a hypothesis with its discriminating measurement.

4. Corrections and additions to carry back into the brief #

1. Section 6 and 13: the Yoo gauge is eq. (2.20), not (2.18). Their horizon times are plot readings; replace "100–250 t_H" by the table in 2.5.
2. Section 4: the kinetic seed is Milligan et al.'s construction; cite it. de Jong 2022 used a radial ODE for ψ, not CTTK.
3. Section 3 and 7: add East 2019 and Ames–Andréasson–Rinne 2023 (collisionless prolate collapse forms horizons in GR); drop Shapiro–Teukolsky 1991 as evidence, since East revisited exactly those configurations with better slicing and found black holes.
4. Section 7: add Andréasson–Kunze–Rein 2011 and Andréasson–Rein 2025 as the rigorous basis of the causal bound; note that Yoo do not check shell crossing and that this memo does (2.3).
5. Section 5: Olabarrieta–Choptuik used maximal-areal coordinates; Akbarian–Choptuik is arXiv:1409.5176 and never runs L = 0; Andréasson–Rein 2006 is gr-qc/0601112; Günther et al. 2021 is arXiv:2009.08163.
6. Section 9: at mΔt = 0.1 the RK4 amplitude loss over the 2×10⁵ steps of a q = 100 run to 300 t_H is 1.4×10⁻³ in amplitude, which fails the 10⁻³ energy target; use mΔt ≤ 0.07 or a symplectic integrator for the scalar sector.
7. Section 11: Harada et al. 2023's velocity-dispersion result is a necessary condition for formation, not a floor; Ebrahimian et al. are Newtonian and say nothing about the inner slope.
8. Section 12: n_s = 0.9743 ± 0.0034; De la Torre Luque erratum Nov 2025; Gottlieb's mechanism is capture inefficiency; the poltergeist verdict is undecidable with current (linear) calculations; LISA reaches 10⁻¹³.
9. Section 14: Seidel–Suen's number could not be checked; the 0.607 m_Pl²/m value is Alcubierre et al. 2003.
10. Codes: COSMOS public without particles; CosmoGRaPH public with particles; OllinSphere-BiB public with cosmological φ² runs; GRChombo has no particle module; Milligan's code private.
11. Novelty statement can be strengthened: zero hits for spherical Einstein–Vlasov PBH thresholds through Sep 2026; Padilla et al. 2026 is spherical, quartic, critical-collapse only.

4a. Cross-check against the independent memo of 23 September (pbh-memo.md) #

A second, independent audit of the same brief (pbh-memo.md, 608 lines, same date) reaches the same decision: GO for V0 and V1, V1b as a gated validation track, and no validated 10⁻⁵ EKG→Vlasov claim. Its LTB check agrees with mine to 0.3%: at μ = 0.33 it finds the first protected trapped sphere at k r = 2.037 (here 2.033, r₁/r_m = 0.830) at 16.17 t_k (here 16.12), and no trapping at μ = 0.32. Its RK4 table (0.28% energy loss at q = 100, mΔt = 0.1; 4.3% at q = 1600) matches the per-step numbers of Section 2.1. Twelve of its points are adopted, with the corresponding changes made in this document:

1. Cold regulator (V1, 6.3): a uniform angular momentum for every shell is a centrifugal barrier that can manufacture a threshold, since a regular distribution has L ∝ R² p_t² near the origin. Regularize instead with a small isotropic dispersion (widths 10⁻³, 3×10⁻⁴, 10⁻⁴ c, each converged at three resolutions before shrinking, plus a second distribution shape with matched low moments) and with L = 0 exactly, reflecting the signed radial momentum through the centre. The uniform-ℓ ladder is removed.
2. Outcome rule (6.7): a horizon counts when θ₊ = 0 and θ₋ < 0 are resolution-stable over a short resolved continuation or validated excision; survival for 0.5 t_AH is recorded, not required. RESOLVED NON-COLLAPSE is a finite-deadline statement and never a proof of permanent support.
3. F3 (2.4): a spherical Vlasov result above 0.33 triggers inspection of the regularization, initial geometry, boundaries and the domain-of-dependence step; it is not automatically a code bug, and a finite-dispersion family is not bound by the cold-limit argument.
4. Flat-slice diagnosis (2.8): it concerns the spatially flat physical slice γ_ij = δ_ij with K solved locally; CTTK on a conformally flat metric whose ψ = e^{ζ/2} carries the curvature is the growing mode.
5. V1b gate 5 (7.3): vary the release slice within one background solution (q_i = 0.5, 1, 2 on the same homogeneous scalar trajectory); a fresh release from rest at each q_i changes the cosmology and tests nothing about slice dependence. Phase shifts at fixed macroscopic data can be physical at finite q.
6. Time step per q (7.4): with the leading RK4 energy-loss exponent mT x⁵/72 and mT ≈ 199 q at D = 300, a 10⁻³ target needs mΔt ≤ 0.08 (q = 100), 0.06 (q = 400), 0.047 (q = 1600), or a symplectic/exponential integrator; 0.07 alone is insufficient at q = 1600.
7. V0 domain (5.4): the comoving light-travel distance from entry to D is 2 r_m (D^{1/3} − 1), i.e. 26/k at D = 250 and 44/k at D = 1000 against the 10/k box. A compensated profile weakens the tidal effect of periodic images, but domain convergence must be shown (L = 20/k and 40/k at coarse resolution), not assumed from L = 25/k.
8. The first-horizon mass, not the accreted mass, sets the wave parameter q γ_AH at formation. V1b records γ_AH at first detection; de Jong's 2×10⁻² would put q = 100–400 runs at q γ_AH ≈ 2–8, the soliton regime H2/F8 worry about.
9. Fallback logic: using EKG at q = 100 as the Vlasov solver when GR-PIC fails is circular. The fallback is an independently validated kinetic method (V1's particle code plus a finite-volume phase-space check on selected endpoints); otherwise the 3D comparison is reported as unresolved.
10. Abundance mapping: a lower threshold at one ellipticity does not license β × 9⁵ (the brief's "Yoo particle slope" column), and a falling deadline boundary does not by itself establish β ∝ σ⁵; both need the shape integral. "Any modulus gives q ≥ 3×10⁹" holds only for Planck-suppressed decay with a stated coefficient. Schrödinger–Poisson methods do encode multi-streaming through interference; their limit near a horizon is the weak-field approximation, not single-stream dust.
11. Cost gating (8): allocations follow measured cost per resolved endpoint from pilot runs. The core-hour figures here are estimates from pbhgr's measured 3.8 μs per particle-step in numpy and an assumed 10–100× compiled speed-up. An Eulerian phase-space Vlasov check (the memo's V1.1) needs gigabytes per array at 2048×512×64 and is reserved for selected endpoints.
12. Peak profiles (6.4): the memo's explicit prescription (lognormal spectra with Δ = 0.1 and 0.5, Gaussian smoothing R_s = 1/k₀, peak height ν = 4, conditional mean over the peak-weighted curvature distribution with the negative-definite Hessian condition) replaces the monochromatic sinc and Δ = 1 placeholders for P2 and P3; the sinc profile stays as the Δ → 0 limit check.

Points where the memo and this document differ, kept as stated here:

• The memo calls the causal bound "conditional"; here it is treated as established for the cold limit on the protected exterior (Hellaby–Lake verified, rigorous analogues cited), with the residual assumptions explicit and tested. Both readings prescribe the same V1 test.
• The memo cites Widrow 1996 (relativistic scalar-field representation of Vlasov matter in the short-de-Broglie limit) and Tanaka–Sasaki 2007 (O(ε²) gradient expansion for a single scalar field in uniform-Hubble slicing) as machinery for V1b initial data; both are verified at abstract level and added to the ledger (1.4). Uniform-Hubble slicing is the CMC slice used here, so Tanaka–Sasaki's O(ε²) solution is the natural upgrade of the leading-order algebraic construction if ε_i = 0.2 is kept.
• COSMOS-S (github.com/cmyoo/cosmos-s, BSD-3) is a public spherically symmetric NR code from Yoo's group (perfect fluid with linear equation of state and massless scalar, scale-up coordinates, fixed mesh refinement, OpenMP). Adding a mass term makes it a third EKG cross-check candidate with Yoo's own gauge philosophy, beside OllinSphere-BiB.
• Esser et al. has a v2 (25 June 2025) whose abstract keeps the 3.7σ exclusion; the memo's statement that uniform IMF priors weaken the limit substantially refers to its Section 5.2 and appendices and was not re-read here.

5. V0 — extend the Yoo et al. particle runs #

5.1 Question and decision tree #

V0 decides what the 0.045/0.050 transition at e = 0.2, p = 0 is. Section 2.5 makes the alternatives quantitative:

5.2 Code routes #

No public code reproduces Yoo's setup out of the box. Three routes, in the order I recommend:

1. Ask Yoo, Escrivà, Harada and Kohri to extend their own runs. The extension is trivial for them (their runs "will continue to run unless manually stopped") and the specific asks are short: μ ∈ {0.030, 0.035, 0.040, 0.045} to 1000 t_H at 80³; μ ∈ {0.045, 0.050, 0.060} at 64³, 80³, 112³ to 500 t_H; particles per cell 1 → 8 at 80³ for μ = 0.045 and 0.050; and the horizon-formation times of Fig. 7 as numbers. This answers the deadline question within weeks and costs nothing here. Its weakness is that it is not an independent check of their method.
2. Build the capability on CosmoGRaPH (public, MIT, BSSN with particles and phase_space_sheet components, used in Giblin et al. 2019 for collisionless matter in an expanding box). Work items: (a) Yoo's long-wavelength CMC initial data, eqs. (3.4)–(3.11) with profile (5.1), and particle masses from eq. (6.10); (b) Yoo's slicing (2.20) with C = 5 and initial lapse (2.21), and the shift (2.22)–(2.23); (c) an apparent-horizon finder (a flow finder on a spherical-harmonic surface, or as a first proxy the 2M/R = 1 crossing along the three axes); (d) the App. B spherical test (μ = 1.3, 40³) against LTB, and the e = 0, μ = 0.3 horizon at ≈ 12 t_H; (e) a uniform 128³–160³ grid in place of Yoo's scale-up coordinates (2.24), or that coordinate map if the code allows it. Estimated effort: 2–3 months for one developer who knows BSSN codes. This is the same capability the brief's V2 needs, so the effort is not wasted if V0 is answered by route 1 first.
3. phantomNR (Einstein Toolkit + Phantom SPH dust): full GR and public, but the SPH deposition assumes a single velocity field, so it does not represent multi-streaming after shell crossing. Not suitable for V0's post-crossing question; listed only for completeness.

I recommend doing route 1 immediately and a two-week feasibility spike on route 2 before committing developer time.

5.3 Setup to replicate (from the paper, not from the brief) #

BSSN; box 0 ≤ X^i ≤ L with L = 10/k; scale-up coordinates X = x − (η/(1+η))(L/π) sin(πx/L) with η = 10; H_i = 5k, a_i = 1 (ε_i = 0.2); CMC initial data eqs. (3.4)–(3.11); profile (5.1) with e = 0.2, p = 0; slicing (2.20) with C = 5 and initial lapse (2.21); shift (2.22)–(2.23); particles regularly aligned (one per cell is implied), masses m_p = (E/Γ)√γ ΔxΔyΔz (6.10), top-hat deposition (6.9); time unit t_H = 100√6/k. Note the two "horizon entry" conventions in the paper: a_k in (2.19) uses k = aH; t_H uses k/aH = √6.

5.4 Diagnostics and acceptance #

• Central lapse α₀(t); apparent-horizon search every 0.5 t_H from t = 0.8 t_C(0) onward; M_AH(t) after formation.
• Hamiltonian constraint max-norm and mean; the max-norm onset must sit at t_C(0) = 0.589 e^{3μ} μ^{−3/2} t_H (Section 2.6) and its amplitude must decrease with grid and particle number (this is the brief's Section 6 acceptance test made concrete).
• Central density, velocity dispersion and virial ratio of the halo at μ = 0.045; the time of the second infall.
• Domain test: the comoving light-travel distance from entry to D is 2 r_m (D^{1/3} − 1), i.e. 26/k at D = 250 and 44/k at D = 1000 against the 10/k box, so run one pair (0.045, 0.050) at L = 20/k and 40/k with coarse resolution and show domain convergence rather than assume it.
• Deposition test: one pair with a cloud-in-cell kernel instead of the top-hat.

Cost. The number of time steps to a fixed cosmic time in a comoving box grows only as D^{1/3}, so a 1000 t_H run costs ≈ 1.6× a 250 t_H run; at 80³–112³ this is 10²–10³ core-hours per run. The brief's 10³–10⁴ core-hours for V0 is generous; the binding cost is code, not compute.

6. V1 — spherical GR-Vlasov #

6.1 Objectives #

1. μ_th(D; σ_v; profile) for D ∈ {100, 300, 1000} t_H, a velocity-dispersion scan, and three profiles (6.4), with the brief's Section 9 error budget.
2. The causal-bound gate: μ_th(σ_v = 0) ≤ 0.33 for the Yoo profile, and a run-by-run check at μ = 0.36 and 0.40 that the horizon appears at r ≥ r₁ no later than t_AH^{LTB}(r₁) (Section 2.4 values).
3. The cold-limit behaviour: convergence in L_reg → 0 and in particle number, and the late-time 2m/R profile (plateau value and inner slope, Section 2.10).
4. M_AH(t) from formation to t_dead at every μ (H6), compared with the guaranteed floor m(r₁).
5. The deadline derivative d ln μ_th/d ln D against the top-hat curve (H8/F12).
6. A one-dimensional replica of Yoo's e = 0 particle runs (μ = 0.3: their 3D horizon at ≈ 12 t_H, LTB central singularity at 8.8 t_H; App. B μ = 1.3).

6.2 Formulation: areal radius with constant-mean-curvature slicing #

pbhgr's polar-areal gauge cannot represent an expanding background (polar slicing needs K^θ_θ = 0, FLRW has K^θ_θ = −H) and stops at 2m/R = 1. The replacement keeps the areal radius R and adds a shift and a slicing condition:

ds² = −α²dt² + A²(dR + β dt)² + R² dΩ².

FLRW is α = 1, A = 1, β = −HR. The gauge is horizon-penetrating (g^{RR} = 1/A² − β²/α² = 1 − 2m/R changes sign with A finite) and, with K = −3H(t) imposed on every slice, coincides with cosmic time in the FLRW exterior, which is where the deadline is defined. With normal-frame matter variables E (energy density), J (radial momentum density), p (radial pressure) and q (tangential pressure), the scheme is fully constrained: four radial solves per time step and no evolved metric variable (Appendix B gives the equations and the FLRW check). The apparent horizon is the surface where 1 − 2m/R = 0 with K^θ_θ > 0 (collapsing), which is θ₊ = 0 in this gauge; the cosmological anti-trapped surface has K^θ_θ < 0 and is excluded by the sign.

A second gauge, geodesic slicing (α = 1, Painlevé–Gullstrand-like) with the same areal coordinate, costs one fewer solve and is used as the gauge cross-check for t_AH and M_AH. The existing polar-areal code stays for the isolated benchmarks (6.6).

Grid: comoving coordinate x = R/a(t) on a stretched mesh (uniform in ln x inside 0.1 r_m, uniform outside), so that the domain [0, x_out] with x_out = 8–12 r_m is fixed while R_out = a x_out expands; ζ(x_out) is below 10⁻¹⁶ of μ, so the exterior is FLRW to machine precision and the boundary data are α = 1, K^θ_θ = −H, m = H²R³/2 (the Yoo profile is compensated: total δM = 0). Particles beyond x_out are pure Hubble flow.

Excision: once the horizon's expansion signs are resolution-stable over a short continuation, particles inside 0.5 R_AH are removed and their mass carried in m_in, exactly as pbhgr's r_in/m_in already does; the CMC lapse collapses inside anyway, so excision is an optimization, not a requirement.

6.3 Initial data #

Use Yoo's long-wavelength CMC data restricted to spherical symmetry (eqs. (3.4)–(3.11) at H_i = 5k), transform to the areal gauge (R = aψ²r√(γ̃_θθ/r²), A, K^θ_θ = Ã^θ_θ − H_i), and solve the two constraints for E and J (in this gauge they are the first-order radial equations of Appendix B read for E and J). Sample N shells by inverse-CDF in the enclosed mass; cold data give every shell the single-stream radial momentum p_R = A·J/E and L = 0 (or L_reg). Validate against the exact LTB (Yoo eqs. (4.5), (4.10), (4.12); my LTB class) until t_C(0). This is a faithful spherical version of Yoo's 3D initial data, so V1 and V0 differ only in dimensionality and method.

Velocity dispersion: in each shell's rest frame add an isotropic Maxwellian with one-dimensional dispersion σ_v; the tangential components give the L distribution. Specify σ_v at horizon entry of the r_m scale (free streaming scales it as 1/a before entry): σ_ent ∈ {0, 10⁻³, 3×10⁻³, 10⁻², 3×10⁻², 10⁻¹} in units of c. Report against two physical maps, both as outputs rather than inputs: Harada et al. 2023 eq. (42), σ_v² ≃ h δ_ent(k̃) for small-scale modes that are nonlinear at entry; and Ebrahimian et al. eqs. (37)–(39), σ_v² ~ (GM/R)(k_0/k)² with the mode amplitude from ζ_k = √(k³P).

Cold-limit regularization at σ_v = 0: no uniform angular momentum. A common L for every shell is a centrifugal barrier that can manufacture a threshold, whereas a regular distribution has L ∝ R² p_t² near the origin. Regularize instead with a small isotropic dispersion of one-dimensional width 10⁻³, 3×10⁻⁴ and 10⁻⁴ (units of c), each resolved at three numerical resolutions before the width is reduced, plus one further positive distribution shape with matched low moments; and with L = 0 exactly, reflecting the signed radial momentum through the centre as pbhgr does now, validated first on a non-self-gravitating test. If the regulators do not approach a common answer, the cold boundary is reported as unresolved. The mass-per-shell scale 2M/N must stay well below the smallest pericentre the regulator implies (with N = 10⁵ it is 10⁻⁵ M_H).

6.4 Profiles #

• P1: Yoo Gaussian, ζ = μ e^{−k²r²/6}, r_m = √6/k. It is the peak-theory typical profile for the spectrum P ∝ k̃³ e^{−3k̃²/2k²} (Yoo 2026 eq. 5.2), so it is already a peak-theory mean profile.
• P2: monochromatic mean profile, ζ = μ sin(k_0 r)/(k_0 r) (Yoo, Harada, Garriga & Kohri 2018 eq. 75, sign chosen so that the peak is an overdensity), r_m = 2.74/k_0. Shells beyond 4.49/k₀ are unbound, the bound core passes Hellaby–Lake (Section 2.3), and its nakedness boundary and r₁ must be computed with ltb_yoo_profile.py before use.
• P3: lognormal spectrum of width Δ ln k = 1, mean profile from Yoo et al. 2018 eqs. (17)–(20) at k_* = k_c (γ ≈ 0.6–0.7).

For P2 and P3 recompute r_m, the entry time (aH r_m = 1), t_C(0), and the LTB nakedness boundary before running.

6.5 Run design #

• Every run goes to t_dead = 1000 t_H; the D = 100 and 300 verdicts are read from the same run via t_AH(μ). μ_th(D) is obtained by interpolating t_AH(μ) = D on the bisection ladder, which is cheaper than separate bisections and gives the full deadline curve.
• Bisection per (profile, σ_v) from [0.03, 1.0] to a 1% relative bracket (7–8 runs), at three resolutions; UNRESOLVED endpoints give one-sided bounds, as in the brief's Section 5 rule. Outcome classification is automated with logged manual overrides.
• Resolutions: N_R ∈ {2000, 4000, 8000} cells with Δx_min ∈ {10⁻³, 5×10⁻⁴, 2.5×10⁻⁴} r_m; particles N ∈ {10⁵, 3×10⁵, 10⁶}; regulator ladder (dispersion widths and L = 0) at the middle resolution; Courant 0.4 and 0.2; ε_i = 0.2 and 0.1; x_out = 8 and 12 r_m.
• Causal-bound runs: μ = 0.36 and 0.40 at σ_v = 0 with the regulator ladder, all three resolutions, with r_AH(t) compared to the LTB apparent-horizon curve outside J⁺(p). Pass criteria from Section 2.4: a horizon at r ≥ 0.51 r_m by 10.1 t_k (μ = 0.36) and at r ≥ 0.37 r_m by 8.6 t_k (μ = 0.40).
• Yoo replica: e = 0, μ = 0.30 (σ_v = 0) and μ = 1.3.

Cost per run (compiled pusher, ~30 ns per particle per RHS): the CFL step in cosmic time is 0.4 A a Δx_min/α, so the step count to D = 1000 is ≈ 2 (r_m/Δx_min) D^{1/3}/0.4 ≈ 5×10⁴ at Δx_min = 10⁻³ r_m; with 10⁵ particles that is ≈ 15 ms per step and ≈ 12 minutes per run, ≈ 2 hours at 10⁶ particles. The full matrix (3 profiles × 6 σ_v × 3 resolutions × ~10 runs, plus ℓ, ε_i, domain and causal-bound runs) is ≈ 600 runs and ≈ 500–1000 core-hours, in line with the brief's ≲ 10³. Memory is below 1 GB per run.

6.6 Validation ladder (each with a pass criterion) #

1. Homogeneous FLRW in the new gauge to 1000 t_H: a(t) within 10⁻⁶ of t^{2/3}; α = 1, A = 1, K^θ_θ = −H to round-off; particle Hubble flow preserved.
2. LTB pre-crossing (P1, μ = 0.3 and 1.3): R(t,r), m(R) and 2m/R agree with the exact solution with second-order convergence in Δx and N^{−1/2} statistical scatter until t_C(0).
3. Dust apparent-horizon curve (μ = 0.5): the found AH tracks R = 2m of the LTB at r > r₁ (outside J⁺(p)); first detection time converges toward t_C(0) = 7.5 t_H as Δx_min → 0 (Yoo's 3D dust needed 11–12 t_H at 80³).
4. Isolated benchmarks in the new gauge: Olabarrieta–Choptuik 2002 family (a) (maximal-areal, M₀* ≈ 1.3, critical max 2M/r ≈ 0.76; use K = 0 in place of K = −3H) within 3%; Rein–Rendall–Schaeffer 1998 family (A* ≈ 0.70, M_BH → 0.11) in the existing polar-areal code within 5%; pbhgr's static polytrope remains static at the 1% level in both gauges.
5. Gauge cross-check: t_AH and M_AH from CMC-areal and PG-areal agree within the resolution error; the polar-areal "BH_approach" time at 2m/R = 8/9 precedes both.
6. Constraint monitors: the unused evolution equation (∂_t A) residual and the mass at R_out relative to FLRW, both below 10⁻³ outside the horizon and converging at second order; the momentum-constraint check pbhgr already has.
7. Causal-bound check (6.5) passes at μ = 0.36 and 0.40.
8. Yoo replica: the e = 0, μ = 0.3 horizon is found earlier than 12 t_H and converges toward 8.8 t_H.

6.7 Diagnostics and classification #

The brief's Section 5 observables and classification rule apply with one change: a horizon counts once θ₊ = 0 and θ₋ < 0 are resolution-stable over a short resolved continuation or validated excision; survival for 0.5 t_AH is recorded, not required, and RESOLVED NON-COLLAPSE is always a finite-deadline statement. Additions: the late-time 2m/R(R) profile, its plateau value and inner logarithmic slope; the velocity-dispersion profile and virial ratio; the number of central passages of the innermost shells; minimum shell radius; outgoing mass and energy flux through R_out; M_AH(t) and its growth exponent.

6.8 Deliverables #

Tables of μ_th(D = 100, 300, 1000; σ_ent; P1–P3) with the Section 9 error budget; t_AH(μ) curves; M_AH(t); the causal-bound test; the LTB–Vlasov gap; the plateau 2m/R against μ and its regulator- and N-convergence; a technical note on the gauge and initial data (publishable on its own if gate 2 of the brief's Section 4 is also cleared by V1b).

6.9 Risks #

• The regulator → 0 limit may not converge for P2 (unbound outer shells, steep interior). Report μ_th against the regulator width and stop when the trend is not monotone.
• Shot noise in J feeds K^θ_θ and the shift; use a triangular-shaped-cloud deposition and treat particle number as a convergence axis (pbhgr's momentum-constraint diagnostic already scales as N^{−1/2}).
• Lapse collapse inside the horizon freezes the interior; excision after 0.5 t_AH handles the long runs.
• Non-commuting limits (σ_v → 0, D → ∞): the plan reports the surface μ_th(σ_v, D) rather than a single number, which is the honest product whatever the limits do.

6.10 V1 status after the first three days of runs (25–27 Sep 2026) #

Code: pbhgr/pbhgr/cosmo_ev.py (areal–CMC solver of Appendix B), cosmo_id.py (Yoo eqs 3.4–3.11 → shells), ltb.py; driver scripts/v1_ladder.py; log pbhgr/docs/v1_devlog.md. Runs on the 2-vCPU VPS at K = 2000 cells, 4–16 shells per cell, x_out = 5 r_m, with the comoving grid re-stretched as a grows so that the physical spacing at the centre stays at 0.3/k (or 0.15/k); 3–45 min per run.

Validation ladder (6.6): FLRW exact to 10⁻¹²; LTB cycloid to 2×10⁻⁴ with second-order convergence; apparent-horizon invariant check exact at μ = 0.5; causal-bound test at μ = 0.36 passed (every shell traps at its LTB t_AH to 10⁻⁴, including r > r₁); static Einstein–Vlasov polytropes (2M/R = 0.11 and 0.13, multi-stream, same gauge with K = 0) static to 10⁻³ over 3 t_dyn; isolated cold-sphere caustic benchmark against a Newtonian shell code agrees to 5–10 % in the weak field (2M/R₀ = 0.002) and, in the strong field (2M/R₀ = 0.02), the GR runaway to a black hole is confirmed by the independent polar-areal code. Three bugs were caught by these tests and fixed (sign of crossing shells on the mirror side; tangential stress double-counted in the lapse source; unresolved pericentres of warm shells, now integrated with per-shell substeps in C).

Results (Yoo profile P1, cold): a horizon forms at every μ tried, from 0.10 down to 0.017, at a universal time after the r_m-scale horizon entry: the first trapped surface at 1.09–1.13 t_C(0) (coordinate time), the earliest per-shell trapping approaching the LTB apparent-horizon curve (t_AH of the labels 0.1–0.2 r_m, 1.01–1.06 t_C(0)), and 1 % of M_H inside the horizon by 1.11–1.15 t_C(0); afterwards the mass grows along the LTB curve. At μ = 0.10 the 0.01 M_H milestone is converged in core spacing (1.149 → 1.152 t_C(0) from 0.3 to 0.15/k); at μ = 0.03 it still moves (1.139 → 1.106) toward the LTB value 1.06, so at low μ the cold horizon time is bracketed between the LTB curve and 1.13 t_C(0). Mechanism: once the central caustic is resolved, the shells that cross the centre are decelerated by the pile-up and do not re-expand past the incoming ones, so the LTB-naked profiles trap essentially as the censored ones do; on grids whose physical core spacing exceeds ≈ 1/k (the fixed comoving grids used at first at low μ) the core is numerically warm and trapping waits for shells with larger turnaround compaction — that produced an apparent rule |2E| ≈ 0.022 and a floor at μ ≈ 0.014, both artefacts now withdrawn.

Horizon times on the re-stretched grid (0.3/k; per-shell proper time of the earliest trapping / CMC coordinate time of the first trapped surface, t_H): μ = 0.10: 25.6 / 28.2; 0.04: 88.1 / 92.8; 0.03: 134.3 / 140.3 (130.4 / 135.0 at 0.15/k); 0.025: 173.5 / 179.6; 0.02: 241.4 / 248.4; 0.017: 306.7 / 313.9. Deadline map from 1.12 t_C(0) = D with t_C(0) = 0.589 e^{3μ} μ^{−3/2} t_H: μ_th(100 t_H) ≈ 0.039 (established: shell-number independent at μ ≥ 0.05, saturated at μ = 0.04 where the converged horizon is at 95.5–96.3 t_H); the 300 t_H entry is ≈ 0.019 (μ = 0.02 traps at 284 t_H coordinate, 267 proper, at converged shell number) and the 1000 t_H entry ≈ 0.0086 (the converged factor is flat at 1.30–1.31 from μ = 0.017 to 0.010; μ = 0.010 traps at 793 t_H coordinate, 770 proper). No cold floor down to μ = 0.010. Yoo et al.'s 3D boundary (0.045–0.050) lies above the spherical 100 t_H value; the difference is for V0/V2.

Consequences for Sections 2.10, 3 and the brief: hypothesis H3 (μ_V,sph ≈ 0.3) is contradicted; the spherical cold threshold is a deadline statement with no floor, the brief's regime F4. The GR ingredient is real (the isolated benchmark: a cold sphere with 2M/R₀ = 0.02 runs away to a black hole in two GR codes while Newtonian shells plateau at 2m/r ≈ 0.06), but the spherical cold trapping itself is the LTB collapse surviving the central crossing.

Velocity dispersion (13 warm runs at 8 velocity samples per cell on the re-stretched grid; the earlier fixed-grid set is superseded): σ_ent ≤ 10⁻² moves the horizon to 0.98–1.06 t_C(0) (earlier than the cold 1.12; the radial part of the kick, v_H σ, binds the inner shells more tightly), while σ_ent = 3×10⁻²–10⁻¹ moves it to 1.15–1.33 t_C(0) with first horizons of 0.02–0.04 M_H behind an angular-momentum barrier. Horizon coordinate times (t_H): μ = 0.04: 93 (cold), 87 (10⁻³), 81 (10⁻²), 101 (3×10⁻²), 100 (10⁻¹); μ = 0.03: 140, 128, 142, 151, 142; μ = 0.05: ≈ 69, 65, —, 75, 76; μ = 0.07: ≈ 44, —, 41, —, 52. So the σ axis is mild: under σ_ent = 0.1 (q ≈ 25) μ_th(100 t_H) ≈ 0.040–0.044 against 0.038 cold and μ_th(300 t_H) ≈ 0.019–0.021 — a ≤ 15 % effect in μ over the whole de Broglie range q ≥ 25. At 32 shells per cell in the core the early-trapping branch survives but shrinks — σ_ent ≤ 10⁻² traps 6–13 % of t_C(0) earlier than cold (μ = 0.04: 87.9 and 89.5 t_H against 95.5; μ = 0.03: 130.4 and 132.8 against 150) and σ_ent = 3×10⁻² is within 1 % of cold (94.3 t_H) — so the axis is settled: the de Broglie range q ≥ 25 moves the thresholds by ≤ 15 % in μ, upward only for q ≲ 30. Sampling converges from 8 velocity samples per cell. The injection (isotropic kick at entry, momenta redshifting as 1/a) is the right model for the scalar-field source: the de Broglie momentum k/a gives σ_ent ≈ √6/q, so σ_ent = 10⁻³–10⁻¹ is q = 2450–25 and the V1b ladder q = 50–1600 is σ_ent = 0.05–0.0015. The substructure maps of 6.5 — Harada et al. 2023 eq. (42), σ_v² ≃ h δ_ent(k) (h = 3/5) released at t* when the collapsing region reaches ρ_vir = 8 ρ_ent δ_ent³, and Ebrahimian et al. eq. (39), σ_v,k² ∼ (GM/R)(k₀/k)² when mode k goes nonlinear during the contraction — are release-during-collapse sources and are now implemented as a per-shell "release at t*" mode; Those sources are run as a per-shell "release at t*" mode and scanned in amplitude and timing (17 runs): a kick of ≈ 0.03 c at turnaround (Ebrahimian k/k₀ = 3, ζ_k = 0.05) prevents trapping to 2.5–3.6 t_C(0) at μ = 0.05, 0.04 and 0.03; ≈ 0.01 c (k/k₀ = 10) delays it to 76, 125 and 260 t_H (1.24, 1.51, 2.09 t_C(0)) and gives no horizon by 450 t_H at μ = 0.025; ≈ 0.003 c (k/k₀ = 30) to 101 and 167 t_H; releasing the same kick later (0.6, 0.3 R_max) shortens the delay (118, 110 t_H at μ = 0.04); a release at ρ_vir (Harada) is irrelevant when it fires at the caustic (k/k̃ = 10: 94 t_H at μ = 0.04 for σ_v = 0.17 and 0.35 c) and matters only when ρ_vir is low enough to fire near turnaround (k/k̃ = 5: 115 t_H). Under substructure released at turnaround the map is μ_th(100 t_H) ≈ 0.040 (k/k₀ = 30), 0.045 (k/k₀ = 10, coinciding with Yoo et al.'s 3D boundary of 0.045–0.050) and no spherical PBH at μ ≤ 0.05 within 300 t_H (k/k₀ = 3); μ_th(300 t_H) ≈ 0.025–0.028. The cold "no floor" result is therefore fragile against realistic small-scale power, in line with Ebrahimian et al.'s Newtonian claim, and the decisive test is 3D with a real spectrum (6.11, item 3b).

Open before these become results: (i) the low-μ cold horizon time is not yet the collisionless answer — at μ = 0.03 it moves later with more shells per cell (1.131 → 1.154 → 1.184 t_C(0) for 4, 8, 16) and earlier with finer cells (1.089 at 0.15/k), the signature of shell-crossing discreteness relaxation in a core that lingers for many dynamical times with only 30–120 shells inside (at μ = 0.10 the result is shell-number independent); the shell-number scan at μ = 0.03 (4 → 8 → 16 → 32 → 64 shells per cell in the core) raises the horizon time 1.131 → 1.154 → 1.185 → 1.209 → 1.209 t_C(0): shell-crossing discreteness relaxation, which saturates once ≈ 1000 shells sit in the core; the converged cold value at μ = 0.03 is 1.21 t_C(0) = 150 t_H (COLLAPSE for D = 300), 8 % later than the 4-shell number; its cell-spacing dependence at converged shell number is 1.5 % per halving (0.15/k: 1.191 t_C(0)), so μ = 0.03 is converged at the 2 % level at 1.18–1.21 t_C(0) = 146–150 t_H; μ = 0.02 at 32 shells per cell gives 1.28 t_C(0) = 284 t_H (267 t_H proper), so the converged cold horizon time rises smoothly toward low μ — 1.115, 1.15, 1.2, 1.28, 1.30 t_C(0) at μ = 0.05, 0.04, 0.03, 0.02, 0.017 (μ = 0.017: 364 t_H; μ = 0.010: 1.31 t_C(0) = 793 t_H, 770 proper) — so item (i) is closed at the 2 % level in cells and shells, and at μ ≥ 0.05 the trapping is fast and LTB-like (shell-number independent at 0.10 and 0.05: 1.115 t_C(0) at 0.05); μ = 0.04 is weakly dependent and saturates at 1.15 t_C(0) = 95.5 t_H (0.01 M_H by 96.3 t_H), so μ_th(100 t_H) ≈ 0.039; (ii) the warm campaign and the release scan on the re-stretched grid; (iii) excision inside the horizon for the M_AH(t) curves; (iv) profiles P2/P3; (v) the strong-field benchmark's onset time (resolution-sensitive, cross-code agreement pending).

6.11 Next steps after the V1 cold map (agreed 26 Sep 2026) #

Order: the velocity-dispersion axis and the profile shapes are finished in 1D on the present VPS (2–3 weeks, no new hardware); the 3D track starts in parallel from day one with semi-analytics and Newtonian 3D, and GR-3D goes to the large machine only once a code with mesh refinement is chosen.

1. Velocity dispersion as an axis of its own (1D, ≈ 1 week; running). Established so far: mild at entry (≤ 15 % in μ up to σ_ent = 0.1), negligible when released at ρ_vir (Harada map), decisive when released right after turnaround (Ebrahimian map: μ_th(100 t_H) ≈ 0.045, floor-like near 0.027). To close: (a) the σ_ent ≤ 10⁻² branch converged at 32 shells per cell in the core on the re-stretched grid (μ = 0.04, 0.03; σ = 10⁻³, 10⁻², one 3×10⁻² control); (b) the release source scanned in amplitude and timing separately — Ebrahimian kick ∝ (k₀/k) with k/k₀ = 3, 10, 30 at turnaround, and a fixed release radius b_nl = 0.95, 0.6, 0.3 R_max at k/k₀ = 10, at μ = 0.05, 0.04, 0.03; Harada at (k/k̃, δ_ent) = (10, 0.2) and (5, 0.05) at μ = 0.04; (c) the report as two surfaces, μ_th(D; q) for the de Broglie source (σ_ent = √6/q, the V1b range q = 25–2450) and μ_th(D; k/k₀, b_nl) for substructure. Decision to take: whether the Ebrahimian source (release when δ_k → 1 in the contracting phase) is the physically preferred one; the recommendation is yes, since it carries a derived release time while the Harada map defers release to ρ_vir by construction.

2. Profile shapes (1D, ≈ 1 week). Generalize cosmo_id.py and ltb.py to a tabulated ζ(r) with spline derivatives (one day plus the LTB validation for each profile), then run the mean peak profiles of the log-normal spectra Δ = 0.1 and 0.5 at ν = 4 (P2, P3) and a broad, near-top-hat control: the cold deadline map at six μ with converged shell number (1–5 h per run) and one σ = 0.03 point each. Running (27 Sep, pbhgr/profiles.py; a shifted-tanh top-hat was rejected — cusp at r = 0, negative energy density in the CMC data — and replaced by the super-Gaussian exp(−(r/2.5)⁴); all three pass the LTB test). Cold map complete (32 shells per cell, first horizon in each profile's own t_H): the narrow-spectrum profile (Δ = 0.1) behaves like the Yoo profile — 1.09, 1.10, 1.16 t_C(0) at μ = 0.05, 0.03, 0.02; the broad-spectrum profile (Δ = 0.5, steep, t_C(r_m)/t_C(0) = 8) traps much later — 1.48 and 1.80 t_C(0) at μ = 0.05 and 0.03, and not by 1.8 t_C(0) at μ = 0.02 (2m/R = 0.08); the flat core traps at 0.97–1.05 t_C(r_m) with 0.77–0.82 M_H at once (top-hat collapse at the edge caustic). With σ_ent = 0.03 the horizon time is unchanged for the narrow-spectrum profile (1.110 t_C(0) at μ = 0.03) and the flat core (0.975 t_C(r_m)), and 12 % later for the broad-spectrum profile (1.663 at μ = 0.05): the steep, slowly fed pile-up is the one sensitive to dispersion. Item 2 is closed: the deadline map is profile-dependent at first order through the synchronization of the collapse (Gaussian-like profiles give the universal 1.1–1.3 t_C(0) seed; steep ones delay trapping by 50–80 % and lose it at low μ within 1.8 t_C(0); flat cores form a near-M_H black hole at once), so any threshold quoted for the early-matter era must state the spectrum shape. The expectation was wrong in direction: a steep profile desynchronizes the collapse, the pile-up is fed slowly and trapping is delayed or absent — the profile shape is a first-order effect on the deadline map, comparable to the substructure axis.

3. Three-dimensional collapse (parallel track, 2–3 months). The first horizon is tiny (R_AH ≈ 0.3/k at a ≈ 10³, 10⁻⁴ of the box in comoving units), so a uniform 3D grid (Yoo's 80³, CosmoGRaPH without refinement) only sees the horizon once it has grown to ≈ 0.05–0.1 M_H, i.e. later than the 1D first-trapping time. Three steps: (a) semi-analytics now (days) — Bond–Myers ellipsoidal collapse for (e, p) drawn from the Zel'dovich distribution at the given ν, collapse times of the three axes relative to t_C(0), and the hoop criterion (Harada et al. 2016) at the last-axis collapse: the expected shift of the deadline map with oblateness before any 3D run, and which (e, p) to simulate. Done (27 Sep, scripts/ellipsoidal_collapse.py): for ν = 4 peaks the Zel'dovich ellipticity has median e = 0.13 and 90th percentile 0.19, for which the last axis collapses at 1.05 and 1.09 t_sph and the pancake forms at 0.91 and 0.85 t_sph (e = 0.3: 1.21 and 0.76; p = ±e/2 changes t₃ by < 1 %). If the horizon needs the third axis, oblateness shifts the deadline map by ≤ 3 % in μ — it does not explain the gap between the spherical 0.039 and Yoo et al.'s 3D 0.045–0.050, which the substructure scan does; whether a pancake of a μ ≈ 0.04 peak already satisfies the hoop criterion is for step (b); (b) Newtonian 3D — done at first order (28–29 Sep, pbhgr/pm3d.py: comoving particle-mesh code, 192³ cells and particles on a 6 r_m box, 2-cell force smoothing to suppress the PM lattice instability found in the shell test; the spherical case collapses 15–20 % before LTB, so ratios to the spherical run are the output). Triaxial peaks with the tidal ellipticity of ν = 4 peaks (checked on the initial deformation tensor) collapse far more anisotropically than the Bond–Myers ellipsoid: at e = 0.13 the pancake forms at 0.77 of the spherical time (Bond–Myers 0.91), the filament at 0.94–0.97, and the third axis collapses at 2.7 × (Bond–Myers 1.05; runs to 2.5 t_C(0)); at e = 0.19: 0.60–0.64 / 0.89 / 1.7–1.76; at e = 0.30: 0.32 / 1.66 / > 1.7 (not run longer). In the Newtonian picture the horizon forms between the filament stage (≈ 0.9–0.97 × the spherical time, if the hoop criterion is met along the two short axes — unresolvable here) and the third-axis collapse (1.7–2.7 ×), so μ_th(100 t_H) lies between the spherical 0.039 and ≈ 0.06–0.08; Yoo et al.'s 3D boundary 0.045–0.050 sits inside this bracket, and closing it needs GR (a filament's compactness), i.e. step (c). A small-scale tail (ζ_rms = 0.05 at k = 2–10 k_p) keeps every Lagrangian set extended (smallest axes 2.5–5 cells) and makes the central mass assembly realization-dependent: one of three seeds fragments the region (13 mass units within 8 cells against 71 spherical), the others still build a dominant clump — the spectrum-driven counterpart of the 1D Ebrahimian-release result. Limits: the code sees collapse times and mass assembly (2M/R ≤ 0.15 at 8 cells), not trapping; one μ (Newtonian ratios are μ-independent); two seeds are not statistics; (c) GR-3D (2–3 months, AX102-1 or CCX43) — the code decision is the main one: CosmoGRaPH (public, particles, uniform grid) against a GRChombo/GRTeclyn build with mesh refinement (particles to be added); without refinement only the late horizon of ≳ 0.05 M_H is measurable, with refinement the early one too. The first mandatory test of any 3D code is to reproduce the spherical map at e = 0 (an independent check of V1 as well); then the grid μ = 0.03–0.06 × e = 0–0.3 × deadlines 100–300 t_H. The cheapest route to 3D remains Yoo et al.'s reply (V0): their code and initial data exist.

6.12 GR-3D code comparison and machine recommendation (29 Sep 2026) #

Requirements from V1/3b: cosmological EdS background with one triaxial peak; long-wavelength CMC initial data; collisionless matter (particles) or, as a first pass, the scalar-field matter of the brief's V1b; mesh refinement (the first horizon is 10⁻⁴ of the box in comoving units, the accreted horizon ≳ 0.05 M_H is 1–5 %); an apparent-horizon finder; runs to 1.5–3 t_C(0).

Recommendation: GRChombo as the main code. Stage 1 (immediate): build it, reproduce de Jong et al.'s matter-domination set-up with scalar matter, and run e = 0 against the spherical 1D map (this is at once the mandatory 3D test and the 3D side of V1b); then e = 0.13, 0.19 at μ = 0.04–0.06 to close the bracket [0.039, 0.06–0.08] from 3b. Stage 2: a particle module (deposit T^μν, geodesic pusher in the 3+1 metric, AMR particle handling — preferably on GRTeclyn/AMReX) for the Vlasov limit and the small-scale tail. In parallel ask for GRAMSES and keep the V0 request to Yoo et al. open.

Cost: with a 128³ base grid on a 6 r_m box the number of coarse steps to 1.5 t_C(0) is only ≈ 700 (Δt ∝ a), so a run is hours on 16 cores even with 4–6 refinement levels; memory tens of GB. Machine: Hetzner AX102-1 (16-core Ryzen 9 7950X3D, 128 GB, €257.30/month + €129 setup) — cheaper per month than CCX43 (16 shared-host vCPUs, 64 GB, €275.99) and about twice as fast per core; it also hosts the 1D scans and a PM-3D seed ensemble at the same time. CCX33 (8 vCPU, 32 GB) is enough only for the 1D and PM work; the present 2-vCPU VPS is not (OOM kills at 1.4 GB per 3D run).

7. V1b — spherical Einstein–Klein–Gordon side-check #

7.1 Objective #

Measure Δ(q) = μ_EKG(q) − μ_V,sph and its exponent (H1, H2), clear gates 2 and 5 of the brief's Section 4, and detect F1/F2. Every EKG run is paired with a V1 run at identical ζ, profile and deadline.

7.2 Codes #

• Code A (in-house): a scalar-field matter module in the same areal–CMC solver as V1, so that gauge, grid, horizon finder and diagnostics are shared and the EKG–Vlasov comparison is run-by-run. Field variables on the grid: φ, Π = (∂_t − β∂_R)φ/α, Φ = ∂_Rφ; sources E = (Π² + Φ²/A²)/2 + V, J = −ΠΦ/A, p = (Π² + Φ²/A²)/2 − V, q = (Π² − Φ²/A²)/2 − V with V = m²φ²/2; fourth-order finite differences, Kreiss–Oliger dissipation, RK4 with mΔt ≤ 0.07 (Section 2.1: 0.1 fails the 10⁻³ energy target at q = 100 to 300 t_H), regridding to keep ≥ 12 points per local de Broglie wavelength (the λ_min monitor of the brief), tested at 16.
• Code B (independent): OllinSphere-BiB (Alcubierre; spherical BSSN, 1+log "cosmocf" slicing with Gamma-driver shift, φ² potential, apparent-horizon finder; the shipped scalarDM_pert.par already evolves a φ² scalar background plus a Gaussian perturbation in an expanding universe). It has no license file; the license must be settled with the author before results are published. GRChombo's cartoon reduction is not public and drops out of the plan; Milligan's code is private but a request for their m/H and a cross-check run costs nothing.

7.3 Initial data and gates #

• Separate-universe release at q_i ≈ 1 (ε_i chosen so that q_i = 1): γ_ij = a²e^{2ζ}δ_ij, K = −3H_i, Ã_ij = 0, homogeneous φ̄ at oscillation phase θ, and Π(x) = √(2[E(x) − V(φ̄)]) with E(x) = 3H_i²/8π + ³R(x)/16π from the Hamiltonian constraint (Section 2.9). Both constraints hold exactly for every phase; existence needs E ≥ V(φ̄) everywhere, which holds at these amplitudes. The kinetic seed φ̄ = 0 (Milligan et al.) is the θ = π/2 member and is used for later starts (q_i ≳ 10).
• Gates (from the brief's Section 4, kept as pass/fail): (1) homogeneous limit ΔH/H ≈ 1/q; (2) decaying-mode fraction below 1% at entry from a fit of the oscillation-averaged δ(t) to a + a^{−3/2}; (3) pre-crossing agreement with the LTB solution with residuals ∝ q^{−2}; (4) constraints converge at three resolutions; (5) μ_th unchanged within the error budget for q_i ∈ {0.5, 1, 2} and four phases θ, with the release slices taken on one and the same homogeneous background trajectory rather than fresh releases from rest. Failure of (2) or (5) at the 1% level stops any 3D EKG plan, as the brief says.
• The background is not matter-like during the first e-fold after release; t_k is defined from the run's own a(t)H(t) (entry when aH r_m = 1) and the deadline counts from there.

7.4 q ladder and cost #

q_k ∈ {25 (soliton diagnostic, excluded from fits), 50, 100, 200, 400, 800} with 1600 only if constraint control holds at 800; three resolutions per point; bisection to 1% brackets; fits with free exponent and an out-of-sample check at the highest q. Time step from the RK4 energy-loss exponent mT x⁵/72 with mT ≈ 199 q at D = 300: for a 10⁻³ target mΔt ≤ 0.08 at q = 100, 0.06 at q = 400 and 0.047 at q = 1600, or a symplectic/exponential integrator for the mass term. Steps to 300 t_H are then ≈ 2500 q (q ≤ 400) to 4200 q (q = 1600): minutes at q ≤ 100, ≈ 2 h at q = 400, ≈ 1.5 days at q = 1600 on one core. The brief's 4×10³–10⁴ core-hours stands as an estimate to be replaced by measured cost per resolved endpoint; the q = 1600 runs are the long pole. Record q γ_AH at first horizon detection at every q.

7.5 Observables #

The brief's Section 5 list, plus the cycle-averaged 2m/R profile (compared with V1's plateau), the soliton-core mass against 0.607 m_Pl²/m and against γq, the gradient-energy fraction, and the minimum points per λ_dB outside the horizon (fewer than 8 makes the run UNRESOLVED).

7.6 Risks #

Constraint control at high q is the known failure mode (Milligan et al. lost it below κ = 0.06 despite N = 10,000). Mitigations: the fully constrained metric solve (m and K^θ_θ from the constraints, α from the CMC condition), the λ_dB-based regridding, and the two-code comparison. The Baumgarte–Clough–Giblin existence problem for CTTK-type data is avoided by the algebraic Π construction.

8. Compute, hardware and schedule #

Ordering: V1 first (weeks 1–8), V1b codes in parallel from week 4, the V0 request now and the CosmoGRaPH spike in weeks 1–2. The current 2-vCPU VPS is enough for development and the validation ladder at low resolution; production needs at least 8 dedicated cores (prices in the project memory: CCX33 €138.49, CCX43 €275.99, AX102-1 €257.30 + setup, ex VAT).

Milestones for V1: week 2, FLRW and LTB tests pass in the new gauge; week 3, initial data and AH finder; week 4, benchmarks 4–5; week 6, causal-bound check and the P1 bisection at σ_v = 0 (this alone tests H3 and the F3 gate); week 8, the deadline, σ_v and profile scans with the error budget.

9. Decisions needed from you #

1. V0 route: send the request to Yoo et al. now (I can draft it), start the CosmoGRaPH spike, or both.
2. Second EKG code: accept OllinSphere-BiB pending its license, or wait for a public cartoon code.
3. Hardware: CCX33 now for V1, then AX102-1 for V1b; or AX102-1 directly.
4. Whether the lognormal profile P3 is in V1's first pass or deferred to after P1 and P2.
5. Whether the σ_v axis should be reported against the Harada 2023 map, the Ebrahimian map, or both (the plan does both).
6. Confirm the memo's P2/P3 prescription (lognormal Δ = 0.1 and 0.5 at ν = 4) and whether COSMOS-S joins OllinSphere-BiB as a cross-check EKG code.

Appendix A. Files in /root/PBH/brief_verification/ #

• check_scalings.py: recomputation of Sections 2, 4, 7, 9, 11 and 12 of the brief.
• ltb_yoo_profile.py and ltb_yoo_profile.out: growing-mode LTB, Hellaby–Lake, nakedness boundary, causal-bound quantities, sinc-profile check.
• ltb_collapse_times.py: exact LTB collapse and trapping times against Yoo's Fig. 7.
• agent_reports/A_yoo2026.md … E_observations_novelty_codes.md: verbatim quotes with section and equation numbers for every literature row above.
• pbh-memo.md: the independent memo of 23 September 2026 extracted from pbh-memo.zip and cross-checked in Section 4a.

Appendix B. Equations of the areal–CMC gauge for V1 #

Metric ds² = −α²dt² + A²(dR + βdt)² + R²dΩ²; normal n = (∂_t − β∂_R)/α; radial unit vector ŝ = ∂_R/A. Matter in the normal frame: E, J = −T(n, ŝ), p = T(ŝ, ŝ), q = T(θ̂, θ̂).

1. Extrinsic curvature: K^θ_θ = β/(αR); K^R_R = −[(∂_tA − β∂_RA)/A − ∂_Rβ]/α; K = K^R_R + 2K^θ_θ. The slicing condition is K(t) = −3H(t), the background value.
2. Momentum constraint with the slicing condition, solved outward from K^θ_θ(0) = K/3: ∂_R K^θ_θ + (3/R) K^θ_θ = K/R − 4πAJ, i.e. K^θ_θ(R) = K/3 − (4π/R³)∫₀^R A J R′³ dR′. Then β = αR K^θ_θ.
3. Misner–Sharp mass from 1 − 2m/R = g^{μν}∂_μR∂_νR = 1/A² − β²/α² and the identity ∂_R m = −4πR² T^t_t: ∂_R m = 4πR² (E − R K^θ_θ A J), m(0) = 0 (or m_in with excision). This is independent of α.
4. Radial metric: 1/A² = 1 − 2m/R + R²(K^θ_θ)². (FLRW: 1 − H²R² + H²R² = 1.)
5. Lapse from ∂_tK − β∂_RK = −D²α + α[K_ijK^ij + 4π(E + p + 2q)] with K = K(t): (1/(AR²)) ∂_R(R² ∂_Rα/A) = α[(K^R_R)² + 2(K^θ_θ)² + 4π(E + p + 2q)] + 3Ḣ, with ∂_Rα(0) = 0 and α(R_out) = 1. FLRW check: 3H² + 4πρ + 3Ḣ = 0 for Ḣ = −4πρ.
6. Particles (unit rest mass, covariant radial momentum p_R, conserved L, Γ = √(1 + p_R²/A² + L²/R²)): dR/dt = α p_R/(A²Γ) − β, dp_R/dt = −αΓ ∂_Rα + p_R ∂_Rβ + (α/Γ)[p_R² ∂_RA/A³ + L²/R³]. Normal-frame radial momentum P = p_R/A. Sources per shell of proper volume 4πR²A dR: E = ΣNΓ, J = ΣNP, p = ΣNP²/Γ, q = ΣN(L²/R²)/(2Γ), each divided by that volume.
7. Horizons: θ₊ ∝ 1/A − β/α, so the apparent horizon is 2m/R = 1 with K^θ_θ > 0; the cosmological anti-trapped surface has K^θ_θ < 0.
8. Unused equation for monitoring: the evolution equation of A (equivalently of K^R_R), whose residual is the Hamiltonian-consistency diagnostic; and mass conservation at R_out.
9. Grid motion: with x = R/a(t) the particle equation gains −H R on the right of dR/dt when written for x; the FLRW exterior is then static in x.

Sign conventions in items 1–7 must be re-derived and unit-tested against FLRW, Schwarzschild in Painlevé–Gullstrand form (α = 1, A = 1, β = −√(2M/R)) and the LTB solution before use; they are written here from the standard 3+1 forms and are a specification, not verified code.