From 3524b9e58d23a34f3dd5fd4aedbbd9fdc3478e73 Mon Sep 17 00:00:00 2001 From: jlashner Date: Wed, 20 May 2020 13:14:08 -0700 Subject: [PATCH 1/3] First shot at Classy transfer_model. --- hmf/density_field/transfer_models.py | 67 ++++++++++++++++++++++++++++ 1 file changed, 67 insertions(+) diff --git a/hmf/density_field/transfer_models.py b/hmf/density_field/transfer_models.py index b5d68ce..da1c44d 100644 --- a/hmf/density_field/transfer_models.py +++ b/hmf/density_field/transfer_models.py @@ -20,6 +20,11 @@ except ImportError: HAVE_CAMB = False +try: + import classy +except ImportError as e: + pass + _allfits = ["CAMB", "FromFile", "EH_BAO", "EH_NoBAO", "BBKS", "BondEfs"] @@ -713,3 +718,65 @@ class EH(EH_BAO): """Alias of :class:`EH_BAO`.""" pass + + +class CLASS(FromFile): + """ + Transfer function computed by CLASS. + + Parameters + ---------- + + cosmo : :class:`astropy.cosmology.FLRW` instance + The cosmology used in the calculation + + \*\*model_parameters : unpack-dict + Parameters specific to this model. + + **class_params:** + Dict of params to set for CLASS object + **class_object:** + Custom class object. If not set, a new Class object will be + created. + **flip_t:** + If set to True, will return ln(|T|). + """ + _defaults = {"class_params": None, "class_obj": None, + "flip_T": False} + + def __init__(self, *args, **kwargs): + super(CLASS, self).__init__(*args, **kwargs) + self.spline_fn = None + + if self.params['class_obj'] is not None: + self.class_obj = self.params['class_obj'] + else: + h = self.cosmo.h + class_params = { + 'omega_b': self.cosmo.Ob0 * h ** 2, + 'omega_cdm': (self.cosmo.Om0 - self.cosmo.Ob0) * h ** 2, + 'h': h, + } + if self.params['class_params'] is not None: + class_params.update(self.params['class_params']) + class_params['output'] = 'dTk' + + self.class_obj = classy.Class() + self.class_obj.set(class_params) + + def lnt(self, lnk): + if self.spline_fn is None: + trans = self.class_obj.get_transfer(output_format='camb') + if not trans: + self.class_obj.compute() + trans = self.class_obj.get_transfer(output_format='camb') + + _k, _t = trans['k (h/Mpc)'], trans['-T_tot/k2'] + if self.params['flip_T']: + _t = np.abs(_t) + _lnk, _lnt = np.log(_k), np.log(_t) + + if lnk[0] < _lnk[0]: + _lnk, _lnt = self._check_low_k(_lnk, _lnt, lnk[0]) + self.spline_fn = spline(_lnk, _lnt, k=1) + return self.spline_fn(lnk) From 5057eb9c43a5d50e6c4af46c5592de7b628c0c20 Mon Sep 17 00:00:00 2001 From: jlashner Date: Sat, 23 May 2020 23:54:11 -0700 Subject: [PATCH 2/3] Implemented suggestions and added tests --- CONTRIBUTORS.rst | 2 +- hmf/density_field/transfer_models.py | 55 ++++++++++++++++++---------- tests/test_transfer.py | 8 ++++ 3 files changed, 44 insertions(+), 21 deletions(-) diff --git a/CONTRIBUTORS.rst b/CONTRIBUTORS.rst index 6f3002a..d73c45c 100644 --- a/CONTRIBUTORS.rst +++ b/CONTRIBUTORS.rst @@ -5,7 +5,7 @@ Steven Murray: `@steven-murray `_ Contributors ------------ * Jordan Mirocha (UCLA): `@mirochaj `_ -* `@jlashner `_ +* Jack Lashner (USC): `@jlashner `_ * `@liweitianux `_ Comments, corrections and suggestions diff --git a/hmf/density_field/transfer_models.py b/hmf/density_field/transfer_models.py index da1c44d..9635697 100644 --- a/hmf/density_field/transfer_models.py +++ b/hmf/density_field/transfer_models.py @@ -22,10 +22,12 @@ try: import classy -except ImportError as e: - pass -_allfits = ["CAMB", "FromFile", "EH_BAO", "EH_NoBAO", "BBKS", "BondEfs"] + HAVE_CLASS = True +except ImportError: + HAVE_CLASS = False + +_allfits = ["CAMB", "FromFile", "EH_BAO", "EH_NoBAO", "BBKS", "BondEfs", "CLASS"] class TransferComponent(Component): @@ -726,11 +728,10 @@ class CLASS(FromFile): Parameters ---------- - cosmo : :class:`astropy.cosmology.FLRW` instance The cosmology used in the calculation - \*\*model_parameters : unpack-dict + model_parameters : unpack-dict Parameters specific to this model. **class_params:** @@ -741,38 +742,52 @@ class CLASS(FromFile): **flip_t:** If set to True, will return ln(|T|). """ - _defaults = {"class_params": None, "class_obj": None, - "flip_T": False} + + _defaults = {"class_params": None, "class_obj": None} def __init__(self, *args, **kwargs): + if not HAVE_CLASS: + raise ImportError( + "classy failed to import! Cannot use the Class " "transfer model." + ) super(CLASS, self).__init__(*args, **kwargs) self.spline_fn = None - if self.params['class_obj'] is not None: - self.class_obj = self.params['class_obj'] + if self.params["class_obj"] is not None: + self.class_obj = self.params["class_obj"] + if not isinstance(self.class_obj, classy.Class): + raise TypeError("Parameter 'class_obj' must be of type classy.Class.") else: h = self.cosmo.h class_params = { - 'omega_b': self.cosmo.Ob0 * h ** 2, - 'omega_cdm': (self.cosmo.Om0 - self.cosmo.Ob0) * h ** 2, - 'h': h, + "omega_b": self.cosmo.Ob0 * h ** 2, + "omega_cdm": (self.cosmo.Om0 - self.cosmo.Ob0) * h ** 2, + "Omega_k": self.cosmo.Ok0, + "N_ncdm": 1, # Camb/class defaults differ here so set manually + "N_ur": self.cosmo.Neff - 1, + "m_ncdm": sum(self.cosmo.m_nu).value, + "T_cmb": self.cosmo.Tcmb0.value, + "h": h, } - if self.params['class_params'] is not None: - class_params.update(self.params['class_params']) - class_params['output'] = 'dTk' + if self.params["class_params"] is not None: + class_params.update(self.params["class_params"]) + class_params["output"] = "dTk" self.class_obj = classy.Class() self.class_obj.set(class_params) def lnt(self, lnk): if self.spline_fn is None: - trans = self.class_obj.get_transfer(output_format='camb') - if not trans: + if self.class_obj.state == 0: self.class_obj.compute() - trans = self.class_obj.get_transfer(output_format='camb') + trans = self.class_obj.get_transfer(output_format="camb") - _k, _t = trans['k (h/Mpc)'], trans['-T_tot/k2'] - if self.params['flip_T']: + _k, _t = trans["k (h/Mpc)"], trans["-T_tot/k2"] + if np.any(_t < 0): + print( + "Warning! The transfer function dips below zero! " + "Returning ln(|T|)." + ) _t = np.abs(_t) _lnk, _lnt = np.log(_k), np.log(_t) diff --git a/tests/test_transfer.py b/tests/test_transfer.py index 108f734..7f12189 100644 --- a/tests/test_transfer.py +++ b/tests/test_transfer.py @@ -52,6 +52,14 @@ def test_bondefs(): assert np.isclose(np.exp(t._unnormalised_lnT[0]), 1, rtol=1e-5) +def test_class_camb_defaults(): + # Sets a small lnk_max to capture main features but avoid extrap error + tclass = Transfer(transfer_model="CLASS", lnk_max=0)._unnormalised_lnT + tcamb = Transfer(transfer_model="CAMB", lnk_max=0)._unnormalised_lnT + percent_err = 200 * (tclass - tcamb) / (tclass + tcamb) + assert np.all(percent_err < 1) + + @pytest.mark.skip("Too slow and needs to be constantly updated.") def test_data(datadir): import camb From 09626870ac71c57c4dc7d741f64858fdf9282377 Mon Sep 17 00:00:00 2001 From: jlashner Date: Tue, 28 Jul 2020 11:26:32 -0700 Subject: [PATCH 3/3] Implemented docstring comments --- hmf/density_field/transfer_models.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/hmf/density_field/transfer_models.py b/hmf/density_field/transfer_models.py index 9635697..a2ecdcd 100644 --- a/hmf/density_field/transfer_models.py +++ b/hmf/density_field/transfer_models.py @@ -723,7 +723,7 @@ class EH(EH_BAO): class CLASS(FromFile): - """ + r""" Transfer function computed by CLASS. Parameters @@ -739,8 +739,6 @@ class CLASS(FromFile): **class_object:** Custom class object. If not set, a new Class object will be created. - **flip_t:** - If set to True, will return ln(|T|). """ _defaults = {"class_params": None, "class_obj": None}