GitHub

Cosmologies (full models)

Free radiation, matter and cosmological constant (RMΛ)

SymBoltz.RMΛ — Function
RMΛ(;
    acceleration = false,
    adiabatic = true,
    g = metric(),
    r = radiation(g; adiabatic),
    m = matter(g; adiabatic),
    Λ = cosmological_constant(g; adiabatic),
    K = nothing,
    G = general_relativity(g; acceleration),
    I = harrison_zeldovich(g; name = :I),
    name = :RMΛ, kwargs...
)

Create a simple model with pure non-interacting radiation, matter and cosmological constant.

source
using SymBoltz, Unitful, UnitfulAstro, Plots
M = RMΛ()
pars = Dict(M.r.Ω₀ => 5e-5, M.m.Ω₀ => 0.3, M.g.h => 1.0, M.r.T₀ => NaN) # TODO: don't pass h and T₀ to avoid infinite loop
prob = CosmologyProblem(M, pars)
ks = [1e-3, 1e-2, 1e-1, 1e-0] / u"Mpc"
sol = solve(prob, ks)
p1 = plot(sol, log10(M.g.a), [M.r.ρ, M.m.ρ, M.Λ.ρ, M.G.ρ] ./ M.G.ρ)
p2 = plot(sol, log10(M.g.a), M.g.Φ, ks)
plot(p1, p2, layout = (2, 1), size = (600, 600))
Example block output

Standard ΛCDM

SymBoltz.ΛCDM — Function
ΛCDM(;
    lmax = 6,
    recombination = true,
    reionization = true,
    Hswitch = 1,
    Heswitch = 6,
    acceleration = false,
    g = metric(),
    G = general_relativity(g; acceleration),
    γ = photons(g; lmax),
    ν = massless_neutrinos(g; lmax),
    h = massive_neutrinos(g; lmax),
    c = cold_dark_matter(g; name = :c),
    b = baryons(g; recombination, reionization, Hswitch, Heswitch, name = :b),
    K = nothing,
    Λ = cosmological_constant(g),
    I = harrison_zeldovich(g; name = :I),
    name = :ΛCDM,
    kwargs...
)

Create a ΛCDM model.

source
using SymBoltz, Plots, Unitful, UnitfulAstro
M = ΛCDM()
pars = parameters_Planck18(M)
prob = CosmologyProblem(M, pars)
ks = [1e-3, 1e-2, 1e-1, 1e-0] / u"Mpc"
sol = solve(prob, ks)
p1 = plot(sol, log10(M.g.a), [M.γ.ρ, M.ν.ρ, M.h.ρ, M.b.ρ, M.c.ρ, M.Λ.ρ, M.G.ρ] ./ M.G.ρ)
p2 = plot(sol, log10(M.g.a), M.g.Φ, ks)
plot(p1, p2, layout = (2, 1), size = (600, 600))
Example block output

w₀wₐCDM (CPL parametrization)

SymBoltz.w0waCDM — Function
w0waCDM(; name = :w0waCDM, kwargs...)

Create a ΛCDM model, but with w₀wₐ-parametrized dark energy instead of the cosmological constant.

source
using SymBoltz, Plots, Unitful, UnitfulAstro
M = w0waCDM()
pars = merge(parameters_Planck18(M), Dict(
    M.X.w0 => -0.9,
    M.X.wa => 0.2,
    M.X.cₛ² => 1.0
))
ks = [1e-3, 1e-2, 1e-1, 1e-0] / u"Mpc"
prob = CosmologyProblem(M, pars)
sol = solve(prob, ks)
p1 = plot(sol, log10(M.g.a), M.X.w)
p2 = plot(sol, log10(M.g.a), M.X.δ, ks)
plot(p1, p2, layout = (2, 1), size = (600, 600))
Example block output

Brans-Dicke ΛCDM

SymBoltz.BDΛCDM — Function
BDΛCDM(; name = :BDΛCDM, kwargs...)

Create a ΛCDM model, but with the Brans-Dicke theory of gravity instead of General Relativity.

source

Solve background such that E = G = 1 today, and plot scalar field and Hubble function:

using SymBoltz, Unitful, UnitfulAstro, Plots
M = BDΛCDM()
D = Differential(M.τ)
ks = [1e-3, 1e-2, 1e-1, 1e-0] / u"Mpc"
pars = merge(parameters_Planck18(M), Dict(M.G.ω => 100.0, D(M.G.ϕ) => 0.0)) # unspecified: M.Λ.Ω₀, M.G.ϕ
prob = CosmologyProblem(M, pars, Dict(M.G.ϕ => 0.95, M.Λ.Ω₀ => 0.5), [M.g.ℰ ~ 1, M.G.G ~ 1])
sol = solve(prob, ks; verbose = true)
p1 = plot(sol, log10(M.g.a), [M.g.ℰ, M.G.G], ylims = (0.8, 1.2))
p2 = plot(sol, log10(M.g.a), M.G.δϕ, ks)
plot(p1, p2, layout = (2, 1), size = (600, 600))
Example block output

Brans-Dicke RMΛ

using SymBoltz, Unitful, UnitfulAstro, Plots
M = SymBoltz.BDRMΛ()
D = Differential(M.τ)
pars = Dict(M.r.Ω₀ => 5e-5, M.m.Ω₀ => 0.3, M.g.h => 1.0, M.r.T₀ => 0.0, M.G.ω => 10.0, D(M.G.ϕ) => 0.0) # unspecified: M.Λ.Ω₀, M.G.ϕ
prob = CosmologyProblem(M, pars, Dict(M.G.ϕ => 0.8, M.Λ.Ω₀ => 0.8), [M.g.ℰ ~ 1, M.G.G ~ 1])
k = 1e-0 / u"Mpc"
sol = solve(prob, k; verbose = true)
p1 = plot(sol, log10(M.g.a), M.G.G)
p2 = plot(sol, log10(M.g.a), M.G.δϕ, k)
plot(p1, p2, layout = (2, 1))
Example block output

Quintessence-CDM

SymBoltz.QCDM — Function
QCDM(v; name = :QCDM, kwargs...)

Create a ΛCDM model, but with the quintessence scalar field in the potential v as dark energy instead of the cosmological constant.

source
using SymBoltz, Plots
@parameters V0 N
V = ϕ -> V0 * ϕ^N
M = QCDM(V, I = nothing)
D = Differential(M.τ)
pars = merge(parameters_Planck18(M), Dict(M.Q.ϕ => 1, D(M.Q.ϕ) => 1.0, M.Q.V0 => 1e-2, M.Q.N => 2))
prob = CosmologyProblem(M, pars)
sol = solve(prob) # TODO: shoot so M.g.ℰ ~ 1 today
plot(sol, M.Q.ϕ, M.Q.V, line_z = log10(M.g.a)) # plot V(ϕ(τ))
Example block output