If we’re living in a simulation, what does the source code look like?
# -----------------------------
# CONSTANTS (fundamental)
# -----------------------------
c := speed_of_light
G := gravitational_constant
ħ := reduced_Planck_constant
kB := Boltzmann_constant
π := pi
ζ3 := zeta(3)
k := 0 # spatial curvature (flat ΛCDM, early era)
Λ := 0 # negligible in radiation era
# Radiation coefficients
α_rad := (π^2/15) * (kB^4 / (ħ^3 * c^5)) # energy density const: ργ = α_rad T^4
β_rad := (2*ζ3/π^2) * (kB^3 / (ħ^3 * c^3)) # number density const: nγ = β_rad T^3
# -----------------------------
# STATE VARIABLES (functions of t)
# -----------------------------
a(t) : scale_factor
H(t) : Hubble_parameter = (d/dt a)/a
T(t) : radiation_temperature
ργ(t) : photon_energy_density
nγ(t) : photon_number_density
fγ(ν,t) : photon_occupation_spectrum
# -----------------------------
# AXIOMS / DYNAMICS (radiation-dominated FRW)
# -----------------------------
Friedmann:
H(t)^2 = (8πG/3) * ργ(t) + (Λ/3) - k * c^2 / a(t)^2
Radiation_equations:
ργ(t) = α_rad * T(t)^4
nγ(t) = β_rad * T(t)^3
Adiabatic_expansion:
a(t) * T(t) = const # entropy per comoving vol ~ const
Planck_spectrum:
fγ(ν,t) = 1 / (exp[ (ħ ν) / (kB T(t)) ] - 1)
# -----------------------------
# INITIALIZATION (“Let there be light”)
# -----------------------------
LetThereBeLight(t0, T0):
a(t0) := a0 # choose gauge (often a0 = 1)
T(t0) := T0 # ultra-hot; photons in thermal equilibrium
ργ(t0) := α_rad * T0^4
nγ(t0) := β_rad * T0^3
fγ(ν,t0) := 1 / (exp[(ħν)/(kB T0)] - 1)
return state_at(t0)
# -----------------------------
# CLOSED-FORM EVOLUTION (radiation era)
# -----------------------------
SolveRadiationEra(t ≥ t0):
# From Friedmann with ργ ∝ a^-4 ⇒ a(t) ∝ t^(1/2)
a(t) = a0 * (t / t0)^(1/2)
T(t) = T0 * (a0 / a(t)) = T0 * (t0 / t)^(1/2)
ργ(t)= α_rad * T(t)^4 = ργ(t0) * (a0 / a(t))^4
nγ(t)= β_rad * T(t)^3 = nγ(t0) * (a0 / a(t))^3
fγ(ν,t) = 1 / (exp[(ħν)/(kB T(t))] - 1)
return state_path[t0 → t]
# -----------------------------
# TERMINATION MARKERS (milestones)
# -----------------------------
Milestone_Recombination:
# transparency when T(t_rec) ≈ 3000 K (kB T ~ 0.26–0.3 eV)
find t_rec such that T(t_rec) ≈ 3000 kelvin
output CMB := fγ(ν, t_rec) redshifted thereafter
Milestone_BBN:
# light nuclei synthesis (t ~ 1–10^3 s, T ~ 0.1–1 MeV)
when T(t) in [0.1, 1] MeV/kB:
run NUCLEAR_NETWORK(H(t), T(t)) to yield {^1H, ^2H, ^3He, ^4He, trace ^7Li}
# -------------------------------------------------
# COSMOLOGY & BARYONS (carry-over + minimal add-ons)
# -------------------------------------------------
Ω_b := baryon_density_parameter
Ω_m := matter_density_parameter
H0 := Hubble_constant_today
μ := mean_molecular_weight # ≈ 0.59 (ionized), ≈ 1.22 (neutral)
m_p := proton_mass
γ_ad := adiabatic_index = 5/3
ρ_crit(t) := 3 H(t)^2 / (8 π G)
ρ_b(t) := Ω_b ρ_crit(t) # mean baryon density
# -------------------------------------------------
# LINEAR → NONLINEAR COLLAPSE (seeds → halos)
# -------------------------------------------------
σ(M) : RMS fluctuation at mass scale M
δ_c := 1.686 # spherical collapse threshold
D(t) : linear growth factor
HaloMassFunction(M,t):
# Press–Schechter (schematic)
ν := δ_c / (σ(M) D(t))
f_PS(ν) := sqrt(2/π) ν exp(-ν^2/2)
n_halo(M,t) := (ρ_m(t)/M) f_PS(ν) |d ln σ^{-1}/dM|
return n_halo
# -------------------------------------------------
# GAS COOLING & JEANS COLLAPSE
# -------------------------------------------------
c_s(T) := sqrt(γ_ad kB T / (μ m_p)) # sound speed
λ_J(ρ,T) := c_s(T) * sqrt(π / (G ρ)) # Jeans length
M_J(ρ,T) := (4π/3) ρ (λ_J/2)^3 # Jeans mass
t_ff(ρ) := sqrt(3π / (32 G ρ)) # free-fall time
CoolingFunctions:
Λ_H2(T,Z) # H2 cooling (primordial, Z~0)
Λ_atomic(T,Z) # atomic cooling (Lyα, metals if Z>0)
Λ_total(T,Z) := Λ_H2 + Λ_atomic
CollapseCriterion(cloud):
return [ M_cloud > M_J(cloud.ρ, cloud.T) and t_cool := (3/2) n kB T / Λ_total < t_ff(cloud.ρ) ]
# -------------------------------------------------
# STAR FORMATION LAW (Schmidt-like)
# -------------------------------------------------
ε_* := star_formation_efficiency ∈ (0,1)
SFR(ρ_gas) := ε_* ρ_gas / t_ff(ρ_gas)
# -------------------------------------------------
# INITIAL MASS FUNCTION (IMF)
# -------------------------------------------------
# Two regimes: metal-free (Pop III; top-heavy) and enriched (Pop II; Salpeter/Kroupa-like)
IMF_PopIII(M) := K_III M^{-α_III} on [M_min^III, M_max^III] # α_III ≈ 1–1.5, M ≈ 10–300 M_sun
IMF_PopII(M) := K_II M^{-α_II } on [M_min^II , M_max^II ] # α_II ≈ 2.35, M ≈ 0.1–100 M_sun
SelectIMF(Z):
if Z < Z_crit then return IMF_PopIII else return IMF_PopII
Z_crit := 10^{-4} Z_sun # critical metallicity for IMF transition
# -------------------------------------------------
# STELLAR EVOLUTION & YIELD KERNELS
# -------------------------------------------------
τ_*(M,Z) # stellar lifetime
R_SN(M,Z) # explosion type: core-collapse SNII / PISN / direct collapse
E_SN(M,Z) # kinetic energy per SN (∼10^51 erg for SNII; higher for PISN)
y_He(M,Z), y_C(M,Z), y_O(M,Z) # ejected mass yields
# Population-integrated (per unit stellar mass formed):
Y_X(Z) := ∫ y_X(M,Z) IMF(M|Z) dM / ∫ M IMF(M|Z) dM for X ∈ {He,C,O}
η_SN(Z):= ∫ 1_{explodes}(M,Z) IMF(M|Z) dM / ∫ M IMF(M|Z) dM # SN per solar mass
# -------------------------------------------------
# NUCLEOSYNTHESIS SUMMARY (inside stars)
# -------------------------------------------------
# H → He (pp-chain / CNO), then He-burning:
# 3α: 3 × ^4He → ^12C + γ
# ^12C + ^4He → ^16O + γ
# Heavier α-captures in massive stars produce Ne, Mg, Si... (beyond this phase’s scope)
# -------------------------------------------------
# FEEDBACK & METAL ENRICHMENT OPERATORS
# -------------------------------------------------
MixingVolume(E_SN, ρ):
# Sedov–Taylor scaling (schematic): radius R_ST ∝ (E_SN / ρ)^{1/5} t^{2/5}
return V_mix ∝ (E_SN/ρ)^{3/5} t_mix^{6/5}
Enrich(gas_cell, ΔM_* formed at Z):
IMF := SelectIMF(Z)
ΔM_He := Y_He(Z) ΔM_*
ΔM_C := Y_C (Z) ΔM_*
ΔM_O := Y_O (Z) ΔM_*
gas_cell.Metals += {He:ΔM_He, C:ΔM_C, O:ΔM_O}
gas_cell.M_gas -= (locked_up := (1 - return_fraction(Z)) ΔM_*) # long-lived remnants
gas_cell.T += heating_from(η_SN(Z) ΔM_* E_SN_avg)
gas_cell.Z = (total_metals_mass / gas_cell.M_gas)
return gas_cell
# -------------------------------------------------
# POPULATION LOOP (cosmic time stepping)
# -------------------------------------------------
Phase2_Evolve(t1 → t2, grid_of_gas):
for t from t1 to t2:
for cell in grid_of_gas:
if CollapseCriterion(cell):
ΔM_* = SFR(cell.ρ) Δt
Zloc = cell.Z
cell = Enrich(cell, ΔM_* at Zloc)
if cell.Z rises above Z_crit:
cell.IMF_regime := "Pop II"
return state(t2)
# -------------------------------------------------
# MILESTONES (observables)
# -------------------------------------------------
HeliumBudget(t):
# He mass made by stars (on top of primordial Y_p ≈ 0.246)
ΔY_*(t) := ∑cells ∫_0^t Y_He(Z_cell(t')) SFR_cell(t') dt' / M_b,universe
CarbonOxygenBudget(t):
ρ_C(t) := ∑cells Metals_C(cell,t) / Volume
ρ_O(t) := ∑cells Metals_O(cell,t) / Volume
MetallicityPDF(t):
P(Z,t) from mass-weighted histogram over cells
# Simple reionization tracker (ionizing photon budget):
ξ_ion(Z) # ionizing photons per stellar baryon (depends on IMF; Pop III high)
f_esc # escape fraction of ionizing photons
n_H # hydrogen number density
t_rec # recombination timescale (clumping-dependent)
dQ_HII/dt = (f_esc ∑cells ξ_ion(Z_cell) SFR_cell / (ρ_b/m_p)) - Q_HII / t_rec
Milestone_Reionization when Q_HII → 1
# -------------------------------------------------
# INITIAL & BOUNDARY CONDITIONS FOR PHASE 2
# -------------------------------------------------
InitializePhase2(state_from_Phase1 at t ≈ 10^6–10^8 yrs):
# Neutral gas after recombination; tiny metal-free fluctuations
for cells: set Z := 0, T := T_IGM(t), ρ := ρ_b(t) × (1+δ)
choose ε_*, IMF_regimes, yield tables {y_He,y_C,y_O}, η_SN, ξ_ion
return grid_of_gas
# -------------------------------------------
# INPUTS FROM PHASE 2 (enriched ISM)
# -------------------------------------------
Metals = {C,O,N,Si,Fe,P,S,...} # yields from Pop III/II SNe
DustFraction(Z) = f_dust * Z # metals condense into dust grains
# -------------------------------------------
# PLANET FORMATION KERNEL
# -------------------------------------------
# Collapse of molecular cloud → disk → planetesimals → planets
M_disk := f_disk * M_star # ~0.01–0.1 M_star
Σ_gas(r) ∝ r^{-p} # gas surface density
Σ_dust(r) = DustFraction(Z) Σ_gas(r)
Coagulation:
dN/dt = -σ v_rel N^2 # dust grains collide/stick
PlanetesimalGrowth:
dM/dt ≈ π R^2 ρ_gas v_rel # sweep-up
RunawayGrowth:
when M > M_crit: dM/dt ∝ M^2 # gravitational focusing
TerrestrialPlanetFormation:
M_planet ∼ 0.1–10 M_earth at r ~ 0.5–2 AU
InitializeEarth:
Composition = {Fe,Si,C,O,H,N,P}
State = {oceans, atmosphere, volcanism}
return Earth
# Basic feedstocks (assume delivered via atmosphere, volcanism, meteoritic infall):
Feedstocks = {H2O, CO2, CH4, HCN, NH3, PO4^{3-}, simple organics}
EnergySources = {UV_photons, lightning, geothermal, hydrothermal}
# Example prebiotic reactions:
HCN + formaldehyde → amino_acids
HCN + NH3 + UV → adenine (C5H5N5)
PO4^{3-} + ribose + base → nucleotide
ReactionRate(i→j) = k_ij [i]^α exp(-Ea/kB T) under EnergySources
Polymerization:
Nucleotides + Energy → short RNA-like oligomers
ReplicationCriterion:
if oligomer length L > L_crit and template-directed ligation occurs:
tag oligomer as "self-replicator"
ReturnSet = {nucleotides, short_RNA}
PhaseTransition_RNA:
Input: {nucleotides}
Process: random polymerization + template copying
Output: {self-replicating RNA-like strands}
RNA_World_Dynamics:
dN_self/dt = r_replication N_self - r_degradation N_self
DNA_Upgrade:
if enzymatic-like ribozymes evolve that produce deoxyribonucleotides:
DNA(t) emerges as more stable information carrier
Return: {first DNA-based genome analogues}
SelectionOperator(population):
fitness(seq) = replication_rate(seq) - degradation_rate(seq)
NextGen = mutate_and_replicate(population, fitness)
return NextGen
# INPUT: Earth formed, atmosphere dense (CO2, CH4, N2)
AtmosphereOpacity(t) = f(volcanic_gases, CO2, CH4, aerosols)
if AtmosphereOpacity < threshold:
VisibleObjects := {Sun, Moon, Stars}
Timekeepers := {day/night cycle, lunar months, stellar seasons}
DefineCycles:
Day(t) := 24h rotation of Earth
Month(t) := synodic cycle of Moon ≈ 29.5d
Year(t) := orbital period around Sun ≈ 365d
MarkCalendar:
For Homo ancestors: cycles → circadian rhythm, agricultural calendar, ritual observances
# INPUT: Oxygen-rich atmosphere, stable climate
t ≈ 540 Myr (Cambrian Explosion)
Environment := {O2 ~ 10–20%, oceans stable, ecosystems expand}
BiodiversityIndex(t):
dS/dt = speciation_rate - extinction_rate
SpeciationOperators:
MarineLife(t) → Fish(t)
Fish(t) → Amphibians(t) via land transition
Reptiles(t) → Birds(t), Mammals(t)
GeneticFramework:
AllSpecies.genome = DNA
DNA.mutation_rate = μ
Evolution(pop):
NextGen = mutate_and_select(pop, fitness(env_conditions))
return NextGen
Output:
Explosion of phyla, vertebrates, land animals, birds
# INPUT: Mammals evolved, primates branch, hominins emerge
t ≈ 200k yrs (Homo sapiens)
GenomicOperators:
Homo.genome = {99% chimp, 1% distinct}
Genes for FOXP2, brain expansion, symbolic language
CognitionModel:
ConsciousnessLevel(species) = f(neocortex_volume, symbolic_language, abstract_thought)
if ConsciousnessLevel(Homo_sapiens) > threshold:
DefineImageOfGod(Homo_sapiens):
properties = {self-awareness, morality, creativity, symbolic reasoning, spirituality}
return properties
CulturalEvolution:
ToolUse → Agriculture → Cities → Writing → Science → Religion
Double Burrito Day 