Recompute the expansion and growth grid

The physics is derived in Modified gravity in four functions, and when it is stable. Here the job is mechanical: install the exact versions, run the independent checks, reproduce the grid, read its units, and know which failures to expect.

What the package contains

The archive unpacks into three folders. The reference folder holds three scripts: the model, a script of independent checks and the script that freezes the grid. The templates folder holds three hi_class input files. The reference-data folder holds seven JSON files: the grid, the six Lab presets, the synthetic record, the consistency test, the documented failure cases, the independent checks and a run manifest with versions, settings and SHA-256 hashes.

The synthetic record is labelled as synthetic inside the file. It was drawn from a chosen model with a fixed seed. Nothing in the package is an observation.

Set up the same environment

The frozen files were produced with Python 3.9.6, NumPy 2.0.2 and SciPy 1.13.1. The model uses only these and the standard library. Make a clean virtual environment so that nothing else on your machine is imported.

python3 -m venv ac81-env
. ac81-env/bin/activate
python -m pip install numpy==2.0.2 scipy==1.13.1
unzip ac81-consistency-grid.zip -d ac81
cd ac81

Run every script with the isolated flag, which ignores environment variables and the user site directory. The scripts add their own folder to the import path explicitly.

Run the independent checks first

Before trusting the grid, check that the solver reproduces results it did not produce itself.

python -I reference/check_ac81_limits.py checks

On the reference machine all thirteen checks print PASS and the script exits with status zero. A few lines with their relative errors:

PASS einstein-de-sitter-distance-z1
PASS lcdm-heath-integral-a1.0 1.629434869254709e-05
PASS linder-2005-growth-index-sugra-like 0.0017487692157125423
PASS only-run-stability-cM-0.1
PASS phantom-no-go-w0-0.752-wa-0.86

The matter-only checks compare with closed forms. The Heath check compares the growth with an integral that holds only for a cosmological constant. Its error of $1.6\times10^{-5}$ comes from starting the growing mode at $z=30$ as if matter were alone. The Linder checks compare with a published growth-index fit inside 1.5 times its quoted accuracy. The stability checks reproduce two published statements: running alone needs $c_M\ge0$, and dark energy alone cannot cross minus one. If any line prints FAIL, stop: the installation or the script differs from the reference.

Freeze the grid again

python -I reference/freeze_ac81_grid.py rerun
shasum -a 256 rerun/*.json

The run takes about 20 seconds on the reference machine. Compare the hashes with the files entry of the manifest. A second run on the reference machine reproduced all six files byte for byte; every number is written rounded to ten significant figures. On another platform or library version, compare values instead. The Lab's own browser solver is a useful yardstick: it uses a different integrator, yet it reproduces every stability flag and agrees on $f\sigma_8$ within $3\times10^{-5}$ inside the stable domain.

Read one grid point

Each point lists its background, its two coefficients, the stability verdict and the most negative $Dc_s^2$ found. A stable point also carries rows at six redshifts. Here is the row at $z=0.5$ for the cosmological-constant background with $c_B=1.0$ and $c_M=0.5$:

{"z": 0.5, "comoving_distance": 0.4394663128, "comoving_distance_ratio_to_lcdm": 1.0,
 "hubble_ratio": 1.317668395, "mu": 1.173166224, "sigma": 1.109639004, "m2": 1.088071388,
 "growth_rate": 0.7929344984, "sigma8": 0.633617168, "f_sigma8": 0.5024169113,
 "f_sigma8_gr_same_background": 0.4724286158, "e_g": 0.4338165283,
 "e_g_gr_same_background": 0.4094327473}

Units matter here. The comoving distance is in units of $c/H_0$; multiply by $2997.92458/h$ to get megaparsecs. The Hubble ratio is $H/H_0$. Both $\mu$ and $\Sigma$ are ratios to Newton's constant, and the Planck mass is $M_*^2/M_{\mathrm{Pl}}^2$. The amplitude $\sigma_8$ is relative: every model starts with the same early amplitude, normalised so that Einstein gravity with a cosmological constant reaches $0.81$ today. Compare $f\sigma_8$ between models, not with a survey's absolute value, unless you rescale the amplitude.

Each row carries its own Einstein-gravity growth on the same background. That pairing is the test: $0.5024$ against $0.4724$ is a 6.3% excess at fixed distances.

Across the grid, 140 of 182 points on the cosmological-constant background are stable, counting pure Einstein gravity. On the DESI-like crossing background 118 of 182 are. With no braiding on the crossing background, only running of $c_M\ge1$ survives. That shows the two coefficients can each rescue the crossing within this parametrisation.

Recompute with hi_class

hi_class is the reference Boltzmann code for this basis. Zumalacárregui et al. 2017, hi_class The templates were written for hi_class v3.4.0.1, git commit b425aef17d566429e8cf189ed47d81f5fa84594e, whose keys were read from its parameter file. Build the C executable and run a template:

git clone https://github.com/miguelzuma/hi_class_public.git hi_class
cd hi_class
git checkout b425aef17d566429e8cf189ed47d81f5fa84594e
make class
mkdir -p output
./class ../ac81/templates/ac81-lcdm-background-braided.ini

The repository README points to the same code under the hiclass-code organisation and also offers a Python wrapper. The key lines of a template are:

Omega_Lambda = 0
Omega_fld = 0
Omega_smg = -1
gravity_model = propto_omega
parameters_smg = 0.1, 0.9, 0.6, 0., 1.
expansion_model = lcdm
skip_stability_tests_smg = no
method_qs_smg = fully_dynamic

The five numbers are $c_K$, $c_B$, $c_M$, $c_T$ and the initial Planck mass. Setting the scalar density to minus one tells the code to fix the dark-energy density by requiring a flat universe. The crossing templates switch the expansion model to the $w_0$, $w_a$ form.

The background file then contains columns titled G_eff_smg and slip_eff_smg. From the hi_class source at that commit, the first is the coupling of $\Phi$ and the second is the ratio $\Phi/\Psi$. Recover the Lab's quantities as

with $\eta$ the slip column. For the default model today this mapping turns the Lab's $\mu=1.352$ and $\Sigma=1.192$ into $G_{\mathrm{eff}}=1.031$ and $\eta=0.763$. Those two numbers come from the Python solver alone, so your first hi_class run is their test.

Differences you should expect

hi_class is a different calculation, and three of its choices move numbers.

1. It includes radiation and neutrinos and starts deep in the radiation era. The Python solver neglects radiation and starts the growth at $z=30$. The templates use the same physical matter densities, so the two backgrounds differ only by the small radiation term. The size of the resulting growth difference was not measured for this course.
2. Its fully dynamic perturbations do not assume the quasi-static limit. Where $Dc_s^2$ is close to zero the sound horizon is large and the quasi-static $\mu$ is unreliable. The failure-case file lists running-only examples whose $Dc_s^2$ today falls from 0.087 to 0.0017.
3. It tests stability from the earliest time it evolves, while the Python solver tests from $z=30$ to today. Points whose most negative $Dc_s^2$ sits near zero, at either end, can change verdict.

Errors and expected failures

The template ac81-crossing-no-braiding-expected-failure.ini sets braiding and running to zero on the crossing background. hi_class should stop with a message that begins as below, with the numbers filled in. This outcome is predicted from the source, not observed.

Gradient instability for scalar field perturbations with minimum c_s^2=

The Python solver flags the same model as unstable, with $Dc_s^2=-0.121$ near $z=0.94$. Two more failures need recognising. In the Lab, choosing $w_0+w_a$ above about $-0.52$ leaves more than 1% dark energy at $z=30$, and the matter-era start is wrong; the Lab shows a notice, and the grid stays out of that region. In the grid, an unstable point has no rows at all. That is deliberate: inside an instability $\mu$ can diverge, and a number there would be misleading.

Do not relax the safety margins in hi_class to make a model run. A run that passes only because the instability test was skipped has not computed a valid cosmology. The DESI analysis cuts the unstable region out of its posterior for the same reason. Ishak et al. 2024, §4.2.1

Checklist

1. Same Python, NumPy and SciPy versions; scripts run with the isolated flag.
2. Thirteen independent checks pass before anything else is used.
3. Hashes match the manifest, or values agree within tolerance when the platform differs.
4. Distances converted from $c/H_0$; $\sigma_8$ treated as relative.
5. Each model compared with Einstein gravity on its own background, not with a cosmological-constant model.
6. hi_class run at the stated commit with stability tests on; $\mu$ and $\Sigma$ recovered from the effective coupling and the slip.
7. Disagreements near the stability edge read as a difference in method before they are read as physics.