Skip to content

Commit

Permalink
Enhancement to the perturbative modelling (#1002)
Browse files Browse the repository at this point in the history
* Option to skip updating ptc

* test for ptc update

* Pass extrapolation params along
  • Loading branch information
tilmantroester committed Dec 19, 2022
1 parent bc129ec commit 2e957aa
Show file tree
Hide file tree
Showing 2 changed files with 34 additions and 6 deletions.
22 changes: 16 additions & 6 deletions pyccl/nl_pt/power.py
Original file line number Diff line number Diff line change
Expand Up @@ -445,7 +445,7 @@ def get_pt_pk2d(cosmo, tracer1, tracer2=None, ptc=None,
nonloc_pk_type='nonlinear',
a_arr=None, extrap_order_lok=1, extrap_order_hik=2,
return_ia_bb=False, return_ia_ee_and_bb=False,
return_ptc=False):
return_ptc=False, update_ptc=True):
"""Returns a :class:`~pyccl.pk2d.Pk2D` object containing
the PT power spectrum for two quantities defined by
two :class:`~pyccl.nl_pt.tracers.PTTracer` objects.
Expand Down Expand Up @@ -503,6 +503,8 @@ def get_pt_pk2d(cosmo, tracer1, tracer2=None, ptc=None,
initialized when this function is called. If `False` (default)
the ptc is not output, whether or not it is initialized as
part of the function call.
update_ptc (bool): if `False`, do not recompute the tracer-independent
perturbative quantitities in the fastpt object.
Returns:
:class:`~pyccl.pk2d.Pk2D`: PT power spectrum.
Expand Down Expand Up @@ -530,6 +532,7 @@ def get_pt_pk2d(cosmo, tracer1, tracer2=None, ptc=None,
ptc = PTCalculator(with_dd=with_dd,
with_NC=with_NC,
with_IA=with_IA)
update_ptc = True
if not isinstance(ptc, PTCalculator):
raise TypeError("ptc should be of type `PTCalculator`")

Expand Down Expand Up @@ -558,8 +561,9 @@ def get_pt_pk2d(cosmo, tracer1, tracer2=None, ptc=None,
ga = growth_factor(cosmo, a_arr)
ga4 = ga**4

# update the PTC to have the require Pk components
ptc.update_pk(pk_lin_z0)
if update_ptc:
# update the PTC to have the require Pk components
ptc.update_pk(pk_lin_z0)

if nonlin_pk_type == 'nonlinear':
Pd1d1 = np.array([nonlin_matter_power(cosmo, ptc.ks, a)
Expand Down Expand Up @@ -680,11 +684,15 @@ def get_pt_pk2d(cosmo, tracer1, tracer2=None, ptc=None,
pt_pk_ee = Pk2D(a_arr=a_arr,
lk_arr=np.log(ptc.ks),
pk_arr=p_pt[0].T,
is_logp=False)
is_logp=False,
extrap_order_lok=extrap_order_lok,
extrap_order_hik=extrap_order_hik)
pt_pk_bb = Pk2D(a_arr=a_arr,
lk_arr=np.log(ptc.ks),
pk_arr=p_pt[1].T,
is_logp=False)
is_logp=False,
extrap_order_lok=extrap_order_lok,
extrap_order_hik=extrap_order_hik)
if return_ptc:
return pt_pk_ee, pt_pk_bb, ptc
else:
Expand All @@ -693,7 +701,9 @@ def get_pt_pk2d(cosmo, tracer1, tracer2=None, ptc=None,
pt_pk = Pk2D(a_arr=a_arr,
lk_arr=np.log(ptc.ks),
pk_arr=p_pt.T,
is_logp=False)
is_logp=False,
extrap_order_lok=extrap_order_lok,
extrap_order_hik=extrap_order_hik)
if return_ptc:
return pt_pk, ptc
else:
Expand Down
18 changes: 18 additions & 0 deletions pyccl/tests/test_pt_power.py
Original file line number Diff line number Diff line change
Expand Up @@ -333,3 +333,21 @@ def test_return_ptc():
assert np.allclose(pk_2.eval(ks, 1., COSMO), pk.eval(ks, 1., COSMO))
assert np.allclose(pee2_2.eval(ks, 1., COSMO), pee2.eval(ks, 1., COSMO))
assert np.allclose(pbb2_2.eval(ks, 1., COSMO), pbb2.eval(ks, 1., COSMO))


def test_pt_no_ptc_update():
pee1 = ccl.nl_pt.get_pt_pk2d(COSMO, TRS['TI'], ptc=PTC)
pee1_no_update = ccl.nl_pt.get_pt_pk2d(COSMO, TRS['TI'], ptc=PTC,
update_ptc=False)

COSMO2 = ccl.Cosmology(
Omega_c=COSMO["Omega_c"], Omega_b=COSMO["Omega_b"], h=COSMO["h"],
sigma8=COSMO["sigma8"]+0.05, n_s=COSMO["n_s"],
transfer_function='bbks', matter_power_spectrum='linear')

pee2_no_update = ccl.nl_pt.get_pt_pk2d(COSMO2, TRS['TI'], ptc=PTC,
update_ptc=False)
pee2 = ccl.nl_pt.get_pt_pk2d(COSMO2, TRS['TI'], ptc=PTC)

assert pee1.eval(0.1, 0.9, COSMO) == pee1_no_update.eval(0.1, 0.9, COSMO)
assert pee2.eval(0.1, 0.9, COSMO) != pee2_no_update.eval(0.1, 0.9, COSMO2)

0 comments on commit 2e957aa

Please sign in to comment.