Add option to make einasto alpha a free parameter (#987)

* Add option to make einasto alpha a free parameter

* Add update_parameters method

* Add test

pyccl/halos/profiles.py

@@ -842,27 +842,44 @@ class HaloProfileEinasto(HaloProfile):
     """
     name = 'Einasto'
 
-    def __init__(self, c_M_relation, truncated=True):
+    def __init__(self, c_M_relation, truncated=True, alpha=None):
         if not isinstance(c_M_relation, Concentration):
             raise TypeError("c_M_relation must be of type `Concentration`)")
 
         self.cM = c_M_relation
         self.truncated = truncated
+        self.alpha = alpha
         super(HaloProfileEinasto, self).__init__()
         self.update_precision_fftlog(padding_hi_fftlog=1E2,
                                      padding_lo_fftlog=1E-2,
                                      n_per_decade=1000,
                                      plaw_fourier=-2.)
 
+    def update_parameters(self, alpha=None):
+        """Update any of the parameters associated with this profile.
+        Any parameter set to ``None`` won't be updated.
+
+        Arguments
+        ---------
+        alpha : float
+            Profile shape parameter.
+        """
+        if alpha is not None and alpha != self.alpha:
+            self.alpha = alpha
+
     def _get_cM(self, cosmo, M, a, mdef=None):
         return self.cM.get_concentration(cosmo, M, a, mdef_other=mdef)
 
     def _get_alpha(self, cosmo, M, a, mdef):
-        mdef_vir = MassDef('vir', 'matter')
-        Mvir = mdef.translate_mass(cosmo, M, a, mdef_vir)
-        sM = sigmaM(cosmo, Mvir, a)
-        nu = 1.686 / sM
-        alpha = 0.155 + 0.0095 * nu * nu
+        if self.alpha is None:
+            mdef_vir = MassDef('vir', 'matter')
+            Mvir = mdef.translate_mass(cosmo, M, a, mdef_vir)
+            print(cosmo, M, a, mdef_vir)
+            sM = sigmaM(cosmo, Mvir, a)
+            nu = 1.686 / sM
+            alpha = 0.155 + 0.0095 * nu * nu
+        else:
+            alpha = self.alpha
         return alpha
 
     def _norm(self, M, Rs, c, alpha):

pyccl/tests/test_profiles.py

@@ -112,6 +112,13 @@ def test_cib_2pt_raises():
                          prof2=p_tSZ, mass_def=M200)
 
 
+def test_einasto_smoke():
+    c = ccl.halos.ConcentrationDuffy08(M200)
+    p = ccl.halos.HaloProfileEinasto(c)
+    p.update_parameters(alpha=1.)
+    assert p._get_alpha(COSMO, 1E14, 1., M200) == 1.
+
+
 def test_gnfw_smoke():
     p = ccl.halos.HaloProfilePressureGNFW()
     beta_old = p.beta