Bolt.jl

using Bolt, Plots, ForwardDiff, LaTeXStrings, ThreadTools, Base.Threads

This notebook walks you through some of the capabilities of Bolt - first through how to compute some of the main observables (and their derivatives) and then through how to look at some internal quantities.

First, check how many threads we are using. See here for how to start Julia with more threads, which will lead to better performance.

nthreads()

2

FRW Background setup

# Assign cosmological parameters
𝕡 = CosmoParams(Ω_c = 0.3) # set kwargs like so to change the default values

function FRW_setup(𝕡)
    # Compute expansion history quantities
    bg = Background(𝕡)
    # Compute ionization history (via RECFAST)
    𝕣 = Bolt.RECFAST(bg=bg, Yp=𝕡.Y_p, OmegaB=𝕡.Ω_b, OmegaG=𝕡.Ω_r)
    ih = IonizationHistory(𝕣, 𝕡, bg)
    return bg, ih
end

bg, ih = FRW_setup(𝕡);

$P_{L}(k)$

# Matter power spectrum
kmin, kmax, nk = 10bg.H₀, 5000bg.H₀, 32
ks = log10_k(kmin, kmax, nk) # k grid
pL = tmap(k -> plin(k, 𝕡, bg, ih), ks)
plot(
    ks, vcat(pL...),
    xscale=:log10, yscale=:log10, label=false,
    xlabel=L"k \ [h/\mathrm{Mpc}]", ylabel=L"P_L(k) \ [\mathrm{Mpc}/h]^3"
)

# Gradient wrt Ω_c
# Define a function that changes 𝕡 - need to recompute background components, as they depend on Ω_c
function pL_Ω_c(Ω_c::T) where T
    𝕡 = CosmoParams{T}(Ω_c=Ω_c)
    bg, ih = FRW_setup(𝕡)
    return tmap(k -> plin(k, 𝕡, bg, ih)[1], ks)
end
∂pL_∂Ω_c = ForwardDiff.derivative(pL_Ω_c, 0.3);

plot(
    ks, abs.(∂pL_∂Ω_c),
    xscale=:log10, yscale=:log10, label=false,
    xlabel=L"k \ [h/\mathrm{Mpc}]", ylabel=L"\vert \partial_{\Omega_c} P_L(k) \vert"
)

$C^{TT}(\ell)$

# CMB Cᵀᵀ(ℓ)
ℓmin, ℓmax, nℓ = 2, 20, 1200
ℓs = ℓmin:ℓmax:nℓ
kmin, kmax, nk = 0.1bg.H₀, 1000bg.H₀, 100
ks = quadratic_k(kmin, kmax, nk)
sf = source_grid(𝕡, bg, ih, ks, BasicNewtonian()) # set up LOS source function interpolator
Cᵀᵀ = cltt(ℓs, 𝕡, bg, ih, sf)
p2 = plot(ℓs, (@. ℓs^2 * Cᵀᵀ), label=false, xlabel=L"\ell", ylabel=L"\ell^2 C^{TT}(\ell)")

# gradient wrt Ω_b
function Cᵀᵀ_Ω_b(Ω_b::T) where T # type-stable wrapper
    𝕡 = CosmoParams{T}(Ω_b=Ω_b)
    bg, ih = FRW_setup(𝕡)
    sf = source_grid(𝕡, bg, ih, ks, BasicNewtonian())
    return cltt(ℓs, 𝕡, bg, ih, sf)
end
∂Cᵀᵀ_∂Ω_b = ForwardDiff.derivative(Cᵀᵀ_Ω_b,0.045);

plot(
    ℓs, (@. ℓs^2 * ∂Cᵀᵀ_∂Ω_b),
    label=false, xlabel=L"\ell", ylabel=L"\ell^2 \partial_{\Omega_b} C^{TT}(\ell)"
)

Some internal quantities

Expansion history:

# Conformal time
plot(bg.x_grid, bg.η, xlabel=L"\log(a)", ylabel=L"\eta", label=false, yscale=:log10)

# Hubble parameter
plot(bg.x_grid, bg.ℋ, xlabel=L"\log(a)", ylabel=L"\mathcal{H}", label=false, yscale=:log10)

Ionization history:

# Free electron fraction
plot(bg.x_grid, ih.Xₑ, xlabel=L"\log(a)", ylabel=L"X_{e}", label=false)

# Visibility function
plot(bg.x_grid, ih.g̃, xlabel=L"\log(a)", ylabel=L"g", label=false)
xlims!(-10, 0)