Added neutrino parameters and CMB to interface

More specifically:
- T_cmb
- m_ncdm
- N_ur
- N_ncdm

ccoffe.pxd

@@ -230,6 +230,14 @@ cdef extern from "common.h":
 
         double k_pivot
 
+        double T_cmb
+
+        double N_ur
+
+        double m_ncdm
+
+        int N_ncdm
+
         int only_cross_correlations
 
 

coffe.pyx

@@ -485,6 +485,83 @@ cdef class Coffe:
             self._free_covariance_multipoles()
 
 
+    @property
+    def T_cmb(self):
+        """
+        The average temperature of the CMB.
+        """
+        return self._parameters.T_cmb
+
+    @T_cmb.setter
+    def T_cmb(self, value):
+        _check_parameter('T_cmb', value, (int, float), 0, 10)
+        if not np.allclose(value, self.T_cmb):
+            self._parameters.T_cmb = value
+            self._free_power_spectrum()
+            self._free_integrals()
+            self._free_corrfunc()
+            self._free_multipoles()
+            self._free_covariance_multipoles()
+
+
+    @property
+    def N_ur(self):
+        """
+        The number of ultra-relativistic species.
+        """
+        return self._parameters.N_ur
+
+    @N_ur.setter
+    def N_ur(self, value):
+        _check_parameter('N_ur', value, (int, float), 0, 10)
+        if not np.allclose(value, self.N_ur):
+            self._parameters.N_ur = value
+            self._free_power_spectrum()
+            self._free_integrals()
+            self._free_corrfunc()
+            self._free_multipoles()
+            self._free_covariance_multipoles()
+
+
+    @property
+    def m_ncdm(self):
+        """
+        The sum of masses of non-CDM species (mostly for neutrinos), in units
+        of eV.
+        """
+        return self._parameters.m_ncdm
+
+    @m_ncdm.setter
+    def m_ncdm(self, value):
+        _check_parameter('m_ncdm', value, (int, float), 0, 10)
+        if not np.allclose(value, self.m_ncdm):
+            self._parameters.m_ncdm = value
+            self._free_power_spectrum()
+            self._free_integrals()
+            self._free_corrfunc()
+            self._free_multipoles()
+            self._free_covariance_multipoles()
+
+
+    @property
+    def N_ncdm(self):
+        """
+        The number of (massive!) non-CDM species.
+        """
+        return self._parameters.N_ncdm
+
+    @N_ncdm.setter
+    def N_ncdm(self, value):
+        _check_parameter('N_ncdm', value, int, 0, 10)
+        if not np.allclose(value, self.N_ncdm):
+            self._parameters.N_ncdm = value
+            self._free_power_spectrum()
+            self._free_integrals()
+            self._free_corrfunc()
+            self._free_multipoles()
+            self._free_covariance_multipoles()
+
+
     @property
     def sep(self):
         """

src/common.h

@@ -470,7 +470,7 @@ typedef struct coffe_parameters_t
 
     double sigma8; /* amplitude of primordial power spectrum sigma8 */
 
-    double T_cmb, N_ur, m_ncdm;
+    double T_cmb, N_ur, m_ncdm, Omega0_nu; /* last one = neutrino energy density fraction */
 
     int N_ncdm;
 

src/parser.c

@@ -626,17 +626,19 @@ int coffe_parse_default_parameters(
     /* cosmological parameters */
     par->Omega0_m = 0.3;
     par->Omega0_baryon = 0.05;
-    par->Omega0_cdm = par->Omega0_m - par->Omega0_baryon;
     par->Omega0_gamma = 9e-5;
-    par->Omega0_de = 1 - (par->Omega0_cdm + par->Omega0_baryon + par->Omega0_gamma);
+    par->h = 0.67;
     par->w0 = -1.0;
     par->wa = 0.0;
     par->have_class = 0;
     par->N_ur = 2.0328;
     par->T_cmb = 2.726;
     par->N_ncdm = 1;
-    par->m_ncdm = 0.06;
-    par->h = 0.67;
+    par->m_ncdm = 0.00;
+    /* see eq. (19) of https://arxiv.org/abs/1212.6154 */
+    par->Omega0_nu = par->m_ncdm / 93.14 / par->h / par->h;
+    par->Omega0_cdm = par->Omega0_m - par->Omega0_baryon - par->Omega0_nu;
+    par->Omega0_de = 1 - (par->Omega0_m + par->Omega0_gamma);
     par->k_pivot = 0.05;
     par->sigma8 = 0.8156;
     par->n_s = 0.96;
@@ -911,10 +913,9 @@ int coffe_parser_init(
     /* cosmological parameters */
     parse_double(conf, "omega_m", &par->Omega0_m, COFFE_TRUE);
     parse_double(conf, "omega_baryon", &par->Omega0_baryon, COFFE_TRUE);
-    par->Omega0_cdm = par->Omega0_m - par->Omega0_baryon;
     parse_double(conf, "omega_gamma", &par->Omega0_gamma, COFFE_TRUE);
 
-    par->Omega0_de = 1. - (par->Omega0_cdm + par->Omega0_baryon + par->Omega0_gamma);
+    par->Omega0_de = 1. - (par->Omega0_m + par->Omega0_gamma);
 
     /* mean redshift */
     if (
@@ -1372,10 +1373,14 @@ int coffe_parser_init(
         parse_double(conf, "N_ur", &par->N_ur, COFFE_FALSE);
         parse_double(conf, "m_ncdm", &par->m_ncdm, COFFE_FALSE);
         parse_int(conf, "N_ncdm", &par->N_ncdm, COFFE_FALSE);
+        par->Omega0_nu = par->m_ncdm / 93.14 / par->h / par->h;
+        par->Omega0_cdm = par->Omega0_m - par->Omega0_baryon - par->Omega0_nu;
         parse_external_power_spectrum(par);
     }
     else{
 #endif
+        par->Omega0_cdm = par->Omega0_m - par->Omega0_baryon - par->Omega0_nu;
+
         /* the power spectrum */
         parse_string(
             conf,