Background, boltzmann, and cls in v3

CMakeLists.txt

@@ -19,7 +19,7 @@ option(FORCE_OPENMP "Forcibly use OpenMP " NO)
 set(CCL_SRC
     src/ccl_background.c src/ccl_core.c src/ccl_error.c
     src/ccl_bbks.c src/ccl_f1d.c src/ccl_f2d.c src/ccl_f3d.c
-    src/ccl_power.c src/ccl_bcm.c
+    src/ccl_power.c
     src/ccl_eh.c src/ccl_musigma.c
     src/ccl_utils.c src/ccl_cls.c src/ccl_massfunc.c
     src/ccl_neutrinos.c

benchmarks/test_bcm.py

@@ -6,22 +6,6 @@
 
 def test_bcm():
     cosmo = ccl.Cosmology(
-        Omega_c=0.25,
-        Omega_b=0.05,
-        h=0.7,
-        A_s=2.2e-9,
-        n_s=0.96,
-        Neff=3.046,
-        m_nu_type='normal',
-        m_nu=0.0,
-        Omega_g=0,
-        Omega_k=0,
-        w0=-1,
-        wa=0,
-        bcm_log10Mc=14.0,
-        baryons_power_spectrum='bcm')
-
-    cosmo_nobar = ccl.Cosmology(
         Omega_c=0.25,
         Omega_b=0.05,
         h=0.7,
@@ -41,12 +25,14 @@ def test_bcm():
     k = data[:, 0] * cosmo['h']
     a = 1
 
-    fbcm = ccl.bcm_model_fka(cosmo, k, a)
+    bar = ccl.BaryonicEffectsBCM(log10Mc=14.)
+    fbcm = bar.boost_factor(cosmo, k, a)
     err = np.abs(data[:, 1]/data_nobar[:, 1]/fbcm - 1)
     assert np.allclose(err, 0, atol=BCM_TOLERANCE, rtol=0)
 
-    ratio = (
-        ccl.nonlin_matter_power(cosmo, k, a) /
-        ccl.nonlin_matter_power(cosmo_nobar, k, a))
+    cosmo.compute_nonlin_power()
+    pk_nobar = cosmo.get_nonlin_power()
+    pk_wbar = bar.apply_baryons(cosmo, pk_nobar)
+    ratio = pk_wbar.eval(k, a)/pk_nobar.eval(k, a)
     err = np.abs(data[:, 1]/data_nobar[:, 1]/ratio - 1)
     assert np.allclose(err, 0, atol=BCM_TOLERANCE, rtol=0)

benchmarks/test_cls.py

@@ -68,21 +68,21 @@ def set_up(request):
 
     # Initialize tracers
     trc = {}
-    trc['g1'] = ccl.NumberCountsTracer(cosmo, False,
-                                       (z1, pz1),
-                                       (z2, bz))
-    trc['g2'] = ccl.NumberCountsTracer(cosmo, False,
-                                       (z2, pz2),
-                                       (z2, bz))
-    trc['l1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1))
-    trc['l2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2))
-    trc['i1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1),
+    trc['g1'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z1, pz1),
+                                       bias=(z2, bz))
+    trc['g2'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z2, pz2),
+                                       bias=(z2, bz))
+    trc['l1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1))
+    trc['l2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2))
+    trc['i1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1),
                                       has_shear=False,
                                       ia_bias=(z1, a1))
-    trc['i2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2),
+    trc['i2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2),
                                       has_shear=False,
                                       ia_bias=(z2, a2))
-    trc['ct'] = ccl.CMBLensingTracer(cosmo, 1100.)
+    trc['ct'] = ccl.CMBLensingTracer(cosmo, z_source=1100.)
 
     # Read benchmarks
     def read_bm(fname):

benchmarks/test_correlation.py

@@ -71,21 +71,21 @@ def set_up(request):
 
     # Initialize tracers
     trc = {}
-    trc['g1'] = ccl.NumberCountsTracer(cosmo, False,
-                                       (z1, pz1),
-                                       (z2, bz))
-    trc['g2'] = ccl.NumberCountsTracer(cosmo, False,
-                                       (z2, pz2),
-                                       (z2, bz))
-    trc['l1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1))
-    trc['l2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2))
-    trc['i1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1),
+    trc['g1'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z1, pz1),
+                                       bias=(z2, bz))
+    trc['g2'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z2, pz2),
+                                       bias=(z2, bz))
+    trc['l1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1))
+    trc['l2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2))
+    trc['i1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1),
                                       has_shear=False,
                                       ia_bias=(z1, a1))
-    trc['i2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2),
+    trc['i2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2),
                                       has_shear=False,
                                       ia_bias=(z2, a2))
-    trc['ct'] = ccl.CMBLensingTracer(cosmo, 1100.)
+    trc['ct'] = ccl.CMBLensingTracer(cosmo, z_source=1100.)
 
     # Read benchmarks
     def read_bm(fname):

benchmarks/test_correlation_MG.py

@@ -48,10 +48,12 @@ def set_up(request):
 
     # Initialize tracers
     trc = {}
-    trc['g1'] = ccl.NumberCountsTracer(cosmo, False, (z1, pz1), (z1, bz1))
-    trc['g2'] = ccl.NumberCountsTracer(cosmo, False, (z2, pz2), (z2, bz2))
-    trc['l1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1))
-    trc['l2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2))
+    trc['g1'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z1, pz1), bias=(z1, bz1))
+    trc['g2'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z2, pz2), bias=(z2, bz2))
+    trc['l1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1))
+    trc['l2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2))
 
     # Read benchmarks
     bms = {}

benchmarks/test_correlation_MG2.py

@@ -48,10 +48,12 @@ def set_up(request):
 
     # Initialize tracers
     trc = {}
-    trc['g1'] = ccl.NumberCountsTracer(cosmo, False, (z1, pz1), (z1, bz1))
-    trc['g2'] = ccl.NumberCountsTracer(cosmo, False, (z2, pz2), (z2, bz2))
-    trc['l1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1))
-    trc['l2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2))
+    trc['g1'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z1, pz1), bias=(z1, bz1))
+    trc['g2'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z2, pz2), bias=(z2, bz2))
+    trc['l1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1))
+    trc['l2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2))
 
     # Read benchmarks
     bms = {}

benchmarks/test_correlation_MG3_SD.py

@@ -50,10 +50,12 @@ def set_up(request):
 
     # Initialize tracers
     trc = {}
-    trc['g1'] = ccl.NumberCountsTracer(cosmo, False, (z1, pz1), (z1, bz1))
-    trc['g2'] = ccl.NumberCountsTracer(cosmo, False, (z2, pz2), (z2, bz2))
-    trc['l1'] = ccl.WeakLensingTracer(cosmo, (z1, pz1))
-    trc['l2'] = ccl.WeakLensingTracer(cosmo, (z2, pz2))
+    trc['g1'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z1, pz1), bias=(z1, bz1))
+    trc['g2'] = ccl.NumberCountsTracer(cosmo, has_rsd=False,
+                                       dndz=(z2, pz2), bias=(z2, bz2))
+    trc['l1'] = ccl.WeakLensingTracer(cosmo, dndz=(z1, pz1))
+    trc['l2'] = ccl.WeakLensingTracer(cosmo, dndz=(z2, pz2))
 
     # Read benchmarks
     bms = {}

benchmarks/test_covariances.py

@@ -42,7 +42,7 @@ def test_ssc_WL():
 
     z, nofz = np.loadtxt(os.path.join(data_dir, "ssc_WL_nofz.txt"),
                          unpack=True)
-    WL_tracer = ccl.WeakLensingTracer(cosmo, (z, nofz))
+    WL_tracer = ccl.WeakLensingTracer(cosmo, dndz=(z, nofz))
 
     ell = np.loadtxt(os.path.join(data_dir, "ssc_WL_ell.txt"))
 

benchmarks/test_hod.py

@@ -69,8 +69,8 @@ def _nz_2mrs(z):
                                     normprof1=True, lk_arr=np.log(k_arr),
                                     a_arr=a_arr)
     # C_ell
-    tr = ccl.NumberCountsTracer(cosmo, False, (z_arr, dndz),
-                                (z_arr, np.ones(len(dndz))))
+    tr = ccl.NumberCountsTracer(cosmo, has_rsd=False, dndz=(z_arr, dndz),
+                                bias=(z_arr, np.ones(len(dndz))))
     cl_hod = ccl.angular_cl(cosmo, tr, tr, l_bm, p_of_k_a=pk_hod)
 
     assert np.all(np.fabs(cl_hod/cl_bm-1) < 0.005)

benchmarks/test_power_nu.py

@@ -4,6 +4,10 @@
 
 import pyccl as ccl
 
+# NOTE: We now check up to kmax=23
+# because CLASS v3 has a slight mismatch
+# (6e-3) at higher k wavenumbers.
+KMAX = 23
 POWER_NU_TOL = 1.0E-3
 
 
@@ -43,12 +47,12 @@ def test_power_nu(model):
         warnings.simplefilter("ignore")
         pk_lin_ccl = ccl.linear_matter_power(cosmo, k_lin, a)
 
-    err = np.abs(pk_lin_ccl/pk_lin - 1)
-    assert np.allclose(err, 0, rtol=0, atol=POWER_NU_TOL)
+    assert np.allclose(pk_lin_ccl[k_lin < KMAX],
+                       pk_lin[k_lin < KMAX],
+                       rtol=POWER_NU_TOL)
 
     data_nl = np.loadtxt("./benchmarks/data/model%d_pk_nl_nu.txt" % (model+1))
     k_nl = data_nl[:, 0] * cosmo['h']
     pk_nl = data_nl[:, 1] / (cosmo['h']**3)
     pk_nl_ccl = ccl.nonlin_matter_power(cosmo, k_nl, a)
-    err = np.abs(pk_nl_ccl/pk_nl - 1)
-    assert np.allclose(err, 0, rtol=0, atol=POWER_NU_TOL)
+    assert np.allclose(pk_nl_ccl, pk_nl, rtol=POWER_NU_TOL)

benchmarks/test_timing.py

@@ -22,11 +22,10 @@ def test_timing():
     cosmo.compute_nonlin_power()
 
     start = time.time()
-    t_g = [ccl.NumberCountsTracer(cosmo, True,
-                                  (z_g, ng),
+    t_g = [ccl.NumberCountsTracer(cosmo, has_rsd=True, dndz=(z_g, ng),
                                   bias=(z_g, np.full(len(z_g), b)))
            for ng, b in zip(dNdz_g, b_g)]
-    t_s = [ccl.WeakLensingTracer(cosmo, (z_s, ns)) for ns in dNdz_s]
+    t_s = [ccl.WeakLensingTracer(cosmo, dndz=(z_s, ns)) for ns in dNdz_s]
     t_all = t_g + t_s
 
     cls = np.zeros([nx, nl])

include/ccl.h

@@ -36,7 +36,6 @@
 #include "ccl_correlation.h"
 #include "ccl_massfunc.h"
 #include "ccl_neutrinos.h"
-#include "ccl_bcm.h"
 #include "ccl_bbks.h"
 #include "ccl_eh.h"
 #include "ccl_halofit.h"

pyccl/__init__.py

@@ -26,7 +26,7 @@
 # Hashing, Caching, CCL base, Mutation locks
 from .base import (
     CCLObject,
-    CCLHalosObject,
+    CCLAutoreprObject,
     Caching,
     cache,
     hash_,
@@ -114,9 +114,10 @@
 )
 
 # Cells & Tracers
-from .cls import angular_cl
+from .cells import angular_cl
 from .tracers import (
     Tracer,
+    NzTracer,
     NumberCountsTracer,
     WeakLensingTracer,
     CMBLensingTracer,
@@ -169,9 +170,14 @@
     nfw_profile_2d,
 )
 
+from .baryons import (
+    BaryonicEffects,
+    BaryonicEffectsBCM
+)
+
 
 __all__ = (
-    'lib', 'Caching', 'cache', 'hash_', 'CCLObject', 'CCLHalosObject',
+    'lib', 'Caching', 'cache', 'hash_', 'CCLObject', 'CCLAutoreprObject',
     'UnlockInstance', 'unlock_instance',
     'CCLParameters', 'physical_constants', 'gsl_params', 'spline_params',
     'CCLError', 'CCLWarning', 'CCLDeprecationWarning',
@@ -200,4 +206,5 @@
     'onehalo_matter_power', 'twohalo_matter_power',
     'massfunc', 'halo_bias', 'massfunc_m2r', 'nfw_profile_3d',
     'einasto_profile_3d', 'hernquist_profile_3d', 'nfw_profile_2d',
+    'BaryonicEffects', 'BaryonicEffectsBCM',
 )

pyccl/background.py

@@ -14,9 +14,10 @@
 """
 import numpy as np
 from . import ccllib as lib
-from .pyutils import _vectorize_fn, _vectorize_fn3
-from .pyutils import _vectorize_fn4, _vectorize_fn5
+from .pyutils import (_vectorize_fn, _vectorize_fn3,
+                      _vectorize_fn4, _vectorize_fn5)
 from .parameters import physical_constants
+from .base import warn_api
 
 species_types = {
     'critical': lib.species_crit_label,
@@ -248,7 +249,8 @@ def omega_x(cosmo, a, species):
                           lib.omega_x_vec, cosmo, a, species_types[species])
 
 
-def rho_x(cosmo, a, species, is_comoving=False):
+@warn_api
+def rho_x(cosmo, a, species, *, is_comoving=False):
     """Physical or comoving density as a function of scale factor.
 
     Args:
@@ -280,7 +282,8 @@ def rho_x(cosmo, a, species, is_comoving=False):
         species_types[species], int(is_comoving))
 
 
-def sigma_critical(cosmo, a_lens, a_source):
+@warn_api
+def sigma_critical(cosmo, *, a_lens, a_source):
     """Returns the critical surface mass density.
 
     .. math::