Skip to main content
Sandbox Physics
Starlight and intergalactic hydrogen

Recompute a fixed-source ionization history

Sandbox Physics ·

Reproducing an ionization history begins with a definite physical question: which sources shine, how much radiation escapes, how hydrogen recombines, and which time belongs to each array. Reading a file is the first step. Summing its ionized volume and recalculating the history through another mathematical route come next.

Use the accompanying fixed-source data and methods. First read how neutral fraction becomes brightness and power, then extract the archive. Its records.json lists the five preset manifests and physical inputs. Each record contains sixty-five native epochs.

From the extracted source directory, set the default record's two paths:

export AC48_LOCAL_MANIFEST="../records/reionization-default-a512-20261112/native-manifest.json"
export AC48_LOCAL_INPUTS="../records/reionization-default-a512-20261112/declared-inputs.json"

Reading fields and projecting modes needs only the Python standard library. The full model replay later uses exactly Python 3.12.12, NumPy 2.3.5 and SciPy 1.17.1; requirements.txt pins the numerical libraries. Separate notices license the original code and synthetic data. No third-party software or paper images are bundled.

Read the physical inputs first

The record separates physical inputs, a clock, cell arrays and calculated statistics. Begin with source positions and intrinsic photon rates, followed by escape fraction and common recombination parameter. Positions are fractions of a periodic box of side 20 cMpc, with no inverse Hubble parameter in its length. The reaction network follows hydrogen; the hydrogen mass fraction remains 0.76.

Moving sources with their total luminosity and common recombination prescription fixed should preserve the global history in this uniform model. The spatial arrangement can change. Lowering escape fraction also changes the supply entering the box. Write down this prediction before comparing two records. The global ionized-volume model

A cell value fjf_j means the fraction of that cell occupied by the ionized phase. Its equal-volume mean,

f‾=1N3∑jfj,\overline f=\frac{1}{N^3}\sum_j f_j,

should agree with the separately advanced ionized volume fraction QQ at the same time. Here NN counts cells along each edge. Matching the time labels is part of the comparison.

Recover the original eight-byte values

Each raw value occupies eight bytes in little-endian Float64 format. The x index changes fastest, then y, then z. A three-dimensional array therefore has axis order z, y, x. This convention connects a source, probe or slice to its physical cell.

The compressed file and decompressed array each have a SHA-256 fingerprint. A matching fingerprint identifies the specified bytes; mathematical and physical checks follow. Decompress exactly once and retain original precision for the calculation. Rendering textures can have a different precision and purpose.

From the scripts' directory, inspect one native epoch:

python -S inspect_saved_record.py --manifest ../records/reionization-default-a512-20261112/native-manifest.json --epoch 24 --with-modes --output inspected-frame.json

The argument selects a native index from zero through sixty-four. An index has no time unit. The report includes file checks, the field-mean/counter difference and explicit Fourier differences. It retains the original record's qualification information. Omit both --epoch and --with-modes to read all sixty-five frames and check their bytes and photon balance first.

Default frame 24 has Q=0.458014596Q=0.458014596. At that same physical time, low escape gives Q=0.229007298Q=0.229007298 and high recombination gives Q=0.370305927Q=0.370305927. The dispersed-source record has the default's global counters and can have a different spatial pattern. This comparison connects the input prediction to actual saved numbers.

Connect four counters to one spatial field

Normalize counts per hydrogen nucleus. Let qq count retained ionizations, rr cumulative recombinations, ee photons entering the box and uu surplus photons after saturation. They obey

e=q+r+u.e=q+r+u.

The model advances these counters directly and then checks the relation. Defining recombinations as the residual e−q−ue-q-u would force the identity to hold and conceal a counting error in the evolution.

For a separate arithmetic example, take 1.2 entering photons, 0.3 recombinations and 0.1 surplus photons per nucleus. Retained ionizations are then 0.8 per nucleus. Uniform density makes this the ionized volume fraction too. Comparing the original cell mean with the directly integrated qq connects the spatial calculation to this photon balance.

Recalculate through an integrating factor

Before complete ionization, the common recombination prescription is

dQdt=s(t)−b(t)Q,b(t)=αBCnˉH,0a−3.\frac{\mathrm dQ}{\mathrm dt}=s(t)-b(t)Q,\qquad b(t)=\alpha_{\mathrm B}C\bar n_{\mathrm H,0}a^{-3}.

The source ss is escaped photons per nucleus per second; bb also has units of inverse seconds. The recombination coefficient αB\alpha_{\mathrm B}, common parameter CC and proper hydrogen density set the loss rate.

Starting from zero ionized volume, add the contribution supplied at each earlier time after multiplying by its surviving fraction:

Q(t)=∫tonts(t′)exp⁡ ⁣[−∫t′tb(t′′) dt′′]dt′.Q(t)=\int_{t_{\mathrm{on}}}^{t} s(t')\exp\!\left[-\int_{t'}^{t}b(t'')\,\mathrm dt''\right]\mathrm dt'.

This integrating-factor solution accumulates the prescribed recombination losses in the exponent. It offers a different computational route from advancing the candidate differential equation step by step. Reaching the ionized cap requires a separate event treatment for maintenance recombinations and surplus photons; the uncapped expression alone can grow beyond a volume fraction.

Use the same declared physical inputs for the independent reference. Save its calculated counters before comparing them with the record:

import json
import os
import warnings
import numpy as np
from scipy.integrate import IntegrationWarning
from independent.scalar_quadrature import QuadratureReference

with open(os.environ["AC48_LOCAL_INPUTS"]) as source:
    inputs = json.load(source)
with open(os.environ["AC48_LOCAL_MANIFEST"]) as source:
    record = json.load(source)
reference = QuadratureReference(
    inputs,
    record["settings"]["escapeFraction"],
    record["settings"]["commonRecombinationClumping"],
)
epochs = np.array([row["scaleFactor"] for row in record["epochs"]])
assert epochs[0] == reference.initial_scale
with warnings.catch_warnings():
    warnings.simplefilter("error", IntegrationWarning)
    counters, cap_scale_factor = reference.history(epochs)
with open("independent-counters.npz", "xb") as output:
    np.savez(output, scaleFactors=epochs, fourCounters=counters)

The reference evaluates its own background density, time relation and integrals. Its import graph is separate from the candidate evolution. An integration warning requires investigation before a result can pass. Compare all four saved columns afterward. The original normalized-counter tolerance is 5×10−85\times10^{-8}; actual errors must come from the real run report. The original recombination-coefficient fit and conventions

To recompute the whole candidate history from its inputs, use the fixed environment:

python recompute_history.py --manifest ../records/reionization-default-a512-20261112/native-manifest.json --physical-input ../records/reionization-default-a512-20261112/declared-inputs.json --output recomputed-history

This advances all four counters and all 512 spatial allocation steps. It saves each returned field, the unused volume fraction of every cell and each excess contribution, then compares the sixty-five native fields word for word. Every eight-byte value must match this fixed record. The scalar reference follows the separate integrating-factor route above. An existing output directory is never overwritten.

This final archive was extracted into a fresh directory and used to replay all five complete records. All 325 native fields matched exactly; the largest independent normalized-counter difference was 4.44×10−154.44\times10^{-15}, below the original 5×10−85\times10^{-8} tolerance. Bundled evidence consists of labelled summaries, with full original receipts kept in the scientific archive. This checks reproducibility for fixed inputs while retaining the physical assumptions of uniform density, common recombination and noncausal allocation.

Check time independently too. Scale factor equals one today, but elapsed proper seconds begin when the sources turn on at redshift ten. The relation is dt=da/[aH(a)]\mathrm dt=\mathrm da/[aH(a)]. At 31,557,600 seconds per year, an interval of 3.15576×1015 s3.15576\times10^{15}\,\mathrm s is 100 Myr. This example specifies an elapsed interval. Each native frame takes its time from the record's own elapsed seconds.

Ask a separate question about shape

Many spatial fields can conserve the same photon total. An incremental allocation can depend on its grid and on the processing order of overlapping excess contributions. Change epoch spacing, processing order and spatial resolution separately, and retain the corresponding errors. The original study of photon-conserving spatial allocation

A single source gives a tractable geometric reference. An independently obtained ionized volume corresponds to a comoving sphere radius

R=L(3Q4π)1/3.R=L\left(\frac{3Q}{4\pi}\right)^{1/3}.

Compare an infinite-space sphere with cell volumes only below half a box radius, before periodic wrapping and under the stated uniform assumptions. First resolve the reference's intersection volume within each cell, then examine the candidate-field error. Refining that reference integral checks the reference; refining the physical grid checks the spatial calculation.

Finally project the same original field onto sine and cosine waves. The script sums individual cells without an FFT routine. Coefficients carry the directions of the thirty-two declared modes. Four family averages discard additional directional and phase information. Keeping both the field and its complex coefficients lets you investigate problems hidden by four power values.

Problems with explained answers

  1. Every cell value stays unchanged, but x and z axes are exchanged. Will the mean reveal the error? Must four symmetric mode-family powers reveal it? What else can you check?
  2. A second program uses the candidate's derivative function with a smaller numerical step. What does this test, and why also use the integrating-factor reference?
  3. Move two sources apart with total intrinsic luminosity, escape fraction and common recombination fixed. Predict the global counters and spatial powers. What must be reconsidered if recombination is calculated from local density and ionization state?

In the first problem, the mean is unchanged and symmetric family powers can also remain unchanged. A cell marker at known coordinates, directional complex coefficients and source/probe locations can reveal the exchange. Array layout is part of the physical association.

In the second, the comparison tests numerical step accuracy but still shares the derivative implementation. The independent integrating-factor calculation also checks how that derivative connects sources, density and time. Its value comes from checking another potentially mistaken connection.

In the third, the current uniform model with the same prescribed recombination losses preserves all four global counters and permits different spatial arrangements and powers. A local recombination calculation can change accumulated losses when sources move, requiring a new complete source/loss comparison. Keeping the assumptions explicit tells you which results an input change should affect.

Sources

  1. The global ionized-volume model
  2. The original recombination-coefficient fit and conventions
  3. The original study of photon-conserving spatial allocation

Keep reading

What this page links to directly

Starlight and intergalactic hydrogen · Reading path

Photons, ionized volume and a radio signal

  1. ArticleHow starlight changes intergalactic hydrogen
  2. ArticleWhat does half ionized mean?
  3. ArticleSupplying photons to an expanding volume
  4. ArticleWhere the photons go when ionized bubbles meet
  5. ArticleReading ionized bubbles through neutral hydrogen
  6. You are hereRecompute a fixed-source ionization history
  7. ArticleReading HERA's 21 cm upper limits
Something broken, confusing, or missing on this page?Tell us
What kind of feedback?

Sent with this page’s address, its Lab ID, the site language and your device type. Please leave out personal details. Feedback is kept for up to 12 months.Privacy policy