Page 1 of 1

Class/Camb $C_\ell$ comparisons

Posted: May 18 2020
by Jack Lashner
For a simple sanity check I'm am trying to make sure that the $C_\ell^{TT}$ values I am obtaining with Class and Camb are the same. However, for some reason it looks like I am getting a small but non-negligible difference between the two boltzman-codes for what I think is the same input parameters. I've attached the plots that show the difference between the two packages, and I'm pasting in the code snippet that generated it.

I'm just wondering if there is there a small difference in the cosmological assumptions between the two models that I'm not accounting for? Or is there something wrong with the way that I'm computing $C_\ell$ with either of the packages? Thanks for your help!

Code: Select all

import classy
import camb
from astropy.cosmology import Planck15 as Planck
import numpy as np
import matplotlib.pyplot as plt

h = Planck.h
cosmo = Planck
class LCDM:
    omega_b=Planck.Ob0 * h**2
    omega_cdm=(Planck.Om0 - Planck.Ob0) * h**2
    A_s = 2.097e-9
    tau_reio = 0.0540
    n_s = 0.9652

# Camb initialization
camb_params = camb.CAMBparams()
    H0=LCDM.H0, ombh2=LCDM.omega_b, omch2=LCDM.omega_cdm,
    omk=LCDM.omega_k, nnu=LCDM.Neff, TCMB=LCDM.Tcmb, tau=LCDM.tau_reio
camb_params.WantTransfer = True
camb_params.WantCls = True
camb_params.InitPower.set_params(ns=LCDM.n_s, As=LCDM.A_s);

camb_res = camb.get_results(camb_params)
camb_trans = camb.get_transfer_functions(camb_params)
camb_trans = camb_trans.get_matter_transfer_data().transfer_data[[0,6],:,0] # Returns [k, T_tot(k)]
print("Initialized CAMB")

# Class Initialization
class_params = {
    'h': LCDM.h,
    'omega_b': LCDM.omega_b,
    'omega_cdm': LCDM.omega_cdm,
    'Omega_k': LCDM.omega_k,
    'N_ur': LCDM.Neff,
    'A_s': LCDM.A_s,
    'n_s': LCDM.n_s,
    'T_cmb': LCDM.Tcmb,
    'tau_reio': LCDM.tau_reio,
    'P_k_max_1/Mpc': 100,
    'output': 'dTk,mPk,tCl'
cl = classy.Class()
print("Initialized CLASS")

l, cl_class = cl.raw_cl()['ell'], cl.raw_cl()['tt']
l_fac = l*(l+1)/(2 * np.pi)
cl_class *= l_fac
cl_camb = camb_res.get_unlensed_total_cls(lmax=l[-1]).T[0]

m = l > 3
x = l[m]
y1 = cl_class[m]
y2 = cl_camb[m]

fig, (ax1,ax2) = plt.subplots(
    2, 1, figsize=(8,7),
    gridspec_kw={'height_ratios': [3,1], 'hspace': 0.05}

ax1.plot(x, y1, label='class')
ax1.plot(x, y2, label='camb')
ax1.tick_params(axis='x', which='both', direction='in', labelbottom=False)

perc_err = 100 * (y1-y2) / ((y1 + y2)/2)
ax2.plot(x, perc_err)
ax2.set(ylabel='% Error', ylim=(-1, 2))


Re: Class/Camb $C_\ell$ comparisons

Posted: May 18 2020
by Antony Lewis
You didn't account for $\Omega_\nu$ when calculating omega_cdm? I also don't know if CAMB and Class's default BBN, Halofit and neutrino mass models are quite the same (you could set Y_He explicitly to avoid the former).

Re: Class/Camb $C_\ell$ comparisons

Posted: May 19 2020
by Jack Lashner
Hi Antony, I'm not exactly sure what you mean by accounting for $\Omega_\nu$ when calculating omega_cdm. Do you think you can elaborate?

I also tried setting the Class/Camb Y_He parameters to be equal but no luck...

Re: Class/Camb $C_\ell$ comparisons

Posted: May 19 2020
by Antony Lewis
If $\Omega_{m0}$ includes massive neutrinos then

Code: Select all

omega_cdm=(Planck.Om0 - Planck.Ob0) * h**2
would not be correct.

Re: Class/Camb $C_\ell$ comparisons

Posted: May 20 2020
by Jack Lashner
I just checked and it looks like

Code: Select all

Planck.Om0 + Planck.Ode0 + Planck.Onu0 + Planck.Ogamma0 + Planck.Ok0
is equal to 1, so I'm guessing the neutrino mass is not included in $\Omega_{m,0}$.

Thanks to Julien Lesgourgues I learned that Class and Camb used different defaults for the number of massive neutrinos and the neutrino mass (class's default is to use all massless neutrinos). Fixing this in both cases solved the problem, which may be what you were suggesting.

Thanks for the help!