Strong Bounds on Sum of Neutrino Masses in a 12 Parameter Extended Scenario with Non-Phantom Dynamical Dark Energy ($w(z)\geq -1$)

We obtained constraints on a 12 parameter extended cosmological scenario including non-phantom dynamical dark energy (NPDDE) with CPL parametrization. We also include the six $\Lambda$CDM parameters, number of relativistic neutrino species ($N_{\textrm{eff}}$) and sum over active neutrino masses ($\sum m_{\nu}$), tensor-to-scalar ratio ($r_{0.05}$), and running of the spectral index ($n_{run}$). We use CMB Data from Planck 2015; BAO Measurements from SDSS BOSS DR12, MGS, and 6dFS; SNe Ia Luminosity Distance measurements from the Pantheon Sample; CMB B-mode polarization data from BICEP2/Keck collaboration (BK14); Planck lensing data; and a prior on Hubble constant ($73.24\pm1.74$ km/sec/Mpc) from local measurements (HST). We have found strong bounds on the sum of the active neutrino masses. For instance, a strong bound of $\sum m_{\nu}<$ 0.123 eV (95\% C.L.) comes from Planck+BK14+BAO. Although we are in such an extended parameter space, this bound is stronger than a bound of $\sum m_{\nu}<$ 0.158 eV (95\% C.L.) obtained in $\Lambda \textrm{CDM}+\sum m_{\nu}$ with Planck+BAO. Varying $A_{\textrm{lens}}$ instead of $r_{0.05}$ however leads to weaker bounds on $\sum m_{\nu}$. Inclusion of the HST leads to the standard value of $N_{\textrm{eff}} = 3.045$ being discarded at more than 68\% C.L., which increases to 95\% C.L. when we vary $A_{\textrm{lens}}$ instead of $r_{0.05}$, implying a small preference for dark radiation, driven by the $H_0$ tension.

Apart from inconsistencies among high and low redshift datasets, there are several internal inconsistencies in the Planck data itself. Parameter estimations in ΛCDM differ when considering small scale (l ≥ 1000) and high or intermediate scale (l < 1000) temperature data separately [30]. This is especially true for the measured value of H 0 which is much lower when obtained from the l ≥ 1000 data than when obtained from the l < 1000 data. Another puzzling inconsistency in ΛCDM with Planck data is that the latest measurement of lensing parameter by Planck 2018, A lens = 1.180 ± 0.065 (68% C.L.) in a ΛCDM + A lens model [1] is 2.8σ level higher than ΛCDM prediction of A lens = 1. See also [31,32] on the A lens problem.
A possible explanation for these tensions is the systematics of the observations. But it is also possible that we need physics beyond ΛCDM and standard model of particle physics. These inconsistencies in ΛCDM model and different datasets have motivated several studies of cosmological scenarios in extended parameter spaces [10,21,. Recent studies have also analyzed models with as large as twelve parameters [33][34][35]. The motivation behind studying such a large parameter space is that ΛCDM currently seems to be an over-simplification. Indeed, there is no reason to fix m ν to 0.06 eV (95% C.L.), since it is only approximately the minimum sum of masses required for normal hierarchy of neutrinos and this mass might not be an accurate one. Massive neutrinos produce distinct effects on CMB and large scale structure data and this has been widely studied [62][63][64][65][66][67]. Again, the discrepancy with Planck and HST might be explained by a dark radiation species contributing to N eff [25]. Similarly, existence of tensor perturbations are theoretically well motivated and there seem no reason to not to include them in a analysis.
Apart from massive neutrinos and tensors, another extension to ΛCDM which has recently received a lot of attention is dynamical dark energy, where the dark energy (DE) equation of state (EoS) is not fixed at w = −1 or some other constant, but is varying with time [51]. Dark energy is one of the biggest puzzles, not only in Cosmology, but in the whole of Physics. Currently available datasets, in this era of precision cosmology, can provide us with much better bounds on DE equation of state than it was previously possible. Thus it seems simplistic and unnecessary to assume dark energy as just a cosmological constant, especially when from the quantum field theoretic point of view, it has been a very difficult thing to explain [68]. Hence, in this work, with massive neutrinos, tensors, and dynamical dark energy included, we consider a largely extended cosmology compared to a standard one.
However, we do not include the full dynamical dark energy range. The w = −1 line divides the dynamics of dark energy in two distinct regions, phantom (w < −1) and non-phantom (w ≥ −1). In this work, we discard the phantom region as first done in [54] in the context of cosmological neutrino mass constraints, and specifically consider a non-phantom scenario, since in a universe with phantom dark energy (w < −1), the dark energy density reaches infinity in a finite time leading to dissociation of all bound states, i.e., the so called Big Rip, and seems unphysical in that sense [69,70]. From field theory perspective, Dark energy models with a single scalar field are not able to go across the w = −1 line (i.e., the phantom barrier) and more general models that allow it demand extra degrees of freedom to supply stability gravitationally [71]. Phantom dark energy accommodating field theories are usually plagued with one or more of the following problems like Lorentz violation, unstable vacuum, superluminal modes, ghosts, non-locality, or instability to quantum corrections. On the other hand, however, it is possible to make theories free of such abnormalities by using effects like photon-axion conversion or modified gravity which leads to an apparent w < −1 (see [72] for a brief review), or vacuum phase transition ( [73]), which produces a phantom behaviour of the DE EoS. Nonetheless, there are single scalar field theories like quintessence [74][75][76] which are relatively well motivated theoretically, and are non-phantom in nature. So, in this work we limit ourselves to such theories. Our main motivation to do this work has been to study how effective the currently available datasets are in constraining the cosmological parameters (especially the sum of neutrino masses) in a non-phantom dynamical dark energy scenario instead of a cosmological constant, with minimal assumptions about other parameters coming from the massive neutrinos and tensor sector.
In this work we have first considered a 12 parameter extended scenario with 6 usual ΛCDM parameters, two dynamical dark energy parameters (w 0 − w a approach, CPL parametrization) with w(z) ≥ −1, two neutrino parameters (N eff and m ν ), and two inflationary parameters (r 0.05 and the running of the spectral index, n run ≡ dn s /dln k). We performed a Bayesian analysis to constrain parameters using different combinations of latest available datasets: (1) Cosmic Microwave Background temperature and polarization data from Planck 2015; (2)the latest data released from the BICEP2/Keck Collaboration for the BB mode of the CMB spectrum (BK14); (3) Baryon Acoustic Oscillation Measurements from SDSS III BOSS DR12, MGS and 6dFGS; (4) Supernovae Type Ia Luminosity Distance Measurements from the newly released Pantheon Sample, (5) Planck 2015 lensing data; and (6) the HST Gaussian prior (H 0 = 73.24 ± 1.74 km/sec/Mpc (68% C.L.)) on Hubble constant. Next we turned off the tensor perturbations (i.e., removed r 0.05 ) and constrained this 11 parameter scenario with the same datasets except BK14. Finally we add a new parameter A lens and again constrain this 12 parameter expended space with the mentioned datasets. We emphasize here that this is the first time someone has evaluated the non-phantom dark energy scenario in a 12 parameter extended space. Our main focus in this paper is on sum of neutrino masses, however we provide the constraints on all the varying parameters. Here we would also like to emphasize that we take the datasets at face value, i.e., any discrepancy or tension between datasets in our model is assumed to have a physical reason and not due to unknown systematics involved in the experiments. Also, it is imperative to point out that the best bounds on sum of neutrino masses that we have presented, are strong and comparable or better to the bounds provided by the recently released Planck 2018 results [1] in the ΛCDM + m ν model. Hence our results remain very much relevant although we have used the Planck 2015 data.
It is imperative that we also mention three recent papers which have helped in building the motivation for this work, and also the difference in our analyses with the said papers. In [35], the authors constrained the dark energy dynamics in an extended 12 parameter model, but they included both the phantom and non-phantom sectors of dark energy,and did not consider any tensor modes. In our analysis, we also use 12 parameters, but we have included tensor perturbations, use newer datasets, and more importantly, we have discarded the phantom DE sector as explained above. We would like to mention that this does affect the bounds on m ν greatly, i.e., they become far stronger compared to the case where phantom DE is included. Bounds on other cosmological parameters also improve. The fact that the neutrino mass bounds from cosmology improve greatly in a nonphantom dark energy scenario, and are stronger even compared to the minimal ΛCDM + m ν case was shown by two recent papers [21,54]. However, analyses in both of these papers were done in smaller parameter spaces, and none of these two papers have N eff and A lens as free parameters as we have. Consequently, they have not touched the issues like the possibility of extra radiation species and the A lens problem. Ref. [54] also uses older datasets. In this paper, we have, for the first time, shown that neutrino mass bounds can indeed be stronger than the minimal ΛCDM + m ν model even in a 12 parameter extended scenario if one considers non-phantom dark energy, even though one expects the bounds to relax in such a large extended space. We have also shown that it is possible to effectively constrain cosmological parameters with some reasonable 1-σ ranges with current cosmological data, in a 12 parameter expended scenario with non-phantom dark energy. This paper is arranged as follows: in section 2 we describe the cosmological models used in this paper and the prior ranges of parameters used, along with a brief description of the CPL parametrization. In section 3 we briefly describe the datasets used in this work. In section 4 we present our analysis results. In section 5, we further discuss how the neutrino mass bounds will change in the three models with new values of τ and A lens obtained by the new Planck 2018 collaboration [1]. We provide a discussion and summary in section 6. The main results are in tables 2, 4, and 5.

Models
In this work we have considered 3 different cosmological scenarios to obtain bounds on the cosmological parameters. Below we list the vector of parameters to vary in each of these cosmological scenarios.
For NPDDE11+r model with 12 parameters: In this analysis, the first model, NPDDE11+r, comprises of six additional parameters on top of ΛCDM model. The six parameters of ΛCDM are: present day cold dark matter energy density ω c ≡ Ω c h 2 , present day baryon energy density ω b ≡ Ω b h 2 , reionization optical depth τ , spectral tilt and amplitude of primordial scalar power spectrum n s and A s (evaluated at pivot scale k * = 0.05hM pc −1 ) and Θ s is the ratio between the sound horizon and the angular diameter distance at decoupling.. For our analysis we are adding the following parameters: two dark energy parameters w 0 and w a , effective number of relativistic species at recombination N eff , total neutrino mass m ν , the tensor-to-scalar ratio r 0.05 (evaluated at pivot scale k * = 0.05hM pc −1 ) and the running of spectral index of primordial power spectrum n run (≡ dn s /dln k). In this model, the gravitational lensing amplitude of the CMB angular spectra A lens is fixed at the ΛCDM predicted value of unity.
We also consider two other scenarios. In the NPDDE11 model, we do not run the tensor perturbations and constrain the parameter space considering scalar only perturbations. In the NPDDE11+A lens model we also allow the A lens parameter to vary. This is since the cause of the A lens -anomaly is unknown and therefore it is important to look into the effect of varying A lens on the constraints of rest of the parameter space. CPL Parametrization: For dark energy dynamics we use the famous Chevallier-Polarski-Linder (CPL) parametrization [77,78] which uses a varying equation of state in terms of the redshift z and two parameters w 0 and w a : This uses the Taylor expansion of the equation of state in powers of the scale factor a = 1/(1 + z) and takes only the first two terms. Here w(z = 0) = w 0 is the dark energy EoS at present day (z = 0), whereas w(z → ∞) = w 0 + w a is the dark energy EoS in the beginning of the universe; and w(z) is a monotonic function between these two times. Therefore, to constrain only the NPDDE region of the parameter space i.e. w(z) ≥ −1 it is enough to apply these hard priors: For the cosmological parameters mentioned in eqs. 2.1-2.3, we have assumed flat priors which are listed in table 1, along with hard priors given in eq. 2.5. We obtain the constraints using the Markov Chain Monte Carlo (MCMC) sampler CosmoMC [79] which uses CAMB [80] as the Boltzmann code and the Gelman and Rubin statistics [81] to estimate the convergence of chains. All our chains reached the convergence criterion of R − 1 < 0.01.

Datasets
Below, we provide a description of the datasets used in our analyses. We have used different combinations of these datasets.
Cosmic Microwave Background: Planck 2015 : We have used measurements of the CMB temperature, polarization, and temperature-polarization cross-correlation spectra from the Planck 2015 data release [82,83]. We use a combination of the high-l (30 ≤ l ≤ 2508) and low-l (2 ≤ l ≤ 29) TT likelihood. Along with that, we also include the high-l (30 ≤ l ≤ 1996) EE and TE likelihood and the low-l (2 ≤ l ≤ 29) polarization likelihood. We refer to this whole dataset as Planck.

Baryon Acoustic Oscillations (BAO) Measurements:
We use measurements of the BAO signal obtained from different galaxy surveys in this work. We include the SDSS-III BOSS DR12 Consensus sample ( [84] which includes LOWZ and CMASS galaxy samples at z eff = 0.38, 0.51 and 0.61). Along with it, we also include the DR7 MGS at z eff = 0.15 [85], and the 6dFGS survey at z eff = 0.106 [86]. We denote this full combination as BAO.
Here z eff is the effective redshift of a survey.

Luminosity Distance Measurements from Type Ia Supernovae (SNe Ia):
We also use Supernovae Type-Ia (SNe Ia) luminosity distance measurements from the Pantheon Sample [87]. It comprises of data from 279 Pan-STARRS1 (PS1) Medium Deep Survey SNe Ia (0.03 < z < 0.68) and distance estimates of SNe Ia from SDSS, SNLS, various low-z and HST samples. This combined sample comprises of data from a total of 1048 SNe Ia with a redshift range of 0.01 < z < 2.3 and is the largest one till date. We refer to this data as PAN from now on. This dataset supersedes the Joint Light-curve Analysis (JLA) sample which comprises of information from 740 spectroscopically confirmed SNe Ia [88].

BB Mode Spectrum of CMB :
We use the latest data available from BICEP2/Keck collaboration for the B mode polarization of CMB, which includes all data (range: 20 < l < 330) taken up to and including 2014 [89]. This dataset is denoted as BK14.

Hubble Parameter Measurements:
We use a Gaussian prior of 73.24 ± 1.74 km/sec/Mpc (68% C.L.) on H 0 . This result is a recent 2.4% determination of the local value of the Hubble parameter by [25] which combines the anchor NGC 4258, Milky Way and LMC Cepheids. We denote this prior as HST.
While we use HST in most cases, we also provide some results with a prior with a lower value of H 0 = 71.6 ± 2.7 km/sec/Mpc, which is based on the determination of the Hubble constant from the H0LiCOW programme [26].We call this prior H071p6. This is to compare what happens when we use a H 0 prior that has less tension with Planck than HST.

Planck Lensing Measurements:
We also use the lensing potential measurements via reconstruction through the four point functions of Planck 2015 measurements of CMB [83]. We simply refer to this data as lensing.

Results
We have split the results in the three smaller sections for the three different models we have studied. The description of models and datasets are given at section 2 and section 3 respectively. We have presented the results in the following order: first the NPDDE11+r model, then the NPDDE11 model and lastly the NPDDE11+A lens model. All the marginalized limits quoted in the text or tables are at 68% C.L. whereas upper limits are quoted at 95% C.L.   Table 3. Bounds on cosmological parameters in the ΛCDM model. Marginalized limits are given at 68% C.L.
whereas upper limits are given at 95% C.L. Note that H0 and σ8 are derived parameters.

NPDDE11+r model
Bounds on the NPDDE11+r model parameters are presented in table 2 while the bounds on the  ΛCDM model parameters are presented in table 3. We do not include the bounds from CMB only data as the bounds are not strong enough in the NPDDE11+r model, a finding that corroborates with a recent study [35] which had varied the dark energy EoS in both phantom and non-phantom regions. However adding either BAO or HST with CMB data seems to provide strong bounds on cosmological parameters. Comparing with the bounds on the parameters in the ΛCDM model however we can see that the 68% C.L. spreads of the relevant parameters have increased to different degrees for different parameters. This is an expected phenomenon given the number of parameters has been doubled. Overall the six ΛCDM parameters have been estimated in the NPDDE11+r model with reasonable spreads, showing that it is possible to constrain cosmology effectively in a large parameter space with current datasets. We also find tight bounds on m ν in this model.  m ν < 0.509 eV (95% C.L.) which is incidentally very close to the bound of m ν < 0.49 eV (95% C.L.) reported by Planck collaboration [83] using the same data in the minimal ΛCDM + m ν model. Recent studies [21,54] in smaller parameter spaces have shown that the models comprising of NPDDE provide stronger bounds on m ν than ΛCDM + m ν , because of a degeneracy present between the dark energy EoS w and m ν [91] which leads to the phantom region of the dark energy parameter space preferring larger masses and the non-phantom region preferring smaller masses. However, cosmological datasets usually prefer the phantom region more when the dark energy EoS is allowed to vary both in the phantom and non-phantom regions, which usually leads to weaker bounds on m ν . This work shows that even as a 12 parameter model, the NPDDE11+r is very efficient in constraining m ν , unlike the 12 parameter model in [35], where the bounds on neutrino mass sum loosens up considerably. Contrary to what happens in lower dimensional parameter spaces, the HST prior does not lead to stronger bounds on m ν , as the magnitude of correlation between H 0 and m ν is very small in this model. This small correlation can be explained with the help of mutual degeneracies present between H 0 , m ν , and the DE EoS w. When w is kept constant in a flat ΛCDM + m ν universe, H 0 and m ν are strongly anti-correlated, to keep the distance to the last scattering surface, χ(z dec ) unchanged. Here z dec is the redshift of photon decoupling. χ(z dec ) is sensitive to any changes in the values of H 0 and m ν , and as shown in [21], any change to χ(z dec ) due to increase in m ν can be compensated by decreasing The use of the H071p6 prior, which has a lower value of H 0 than HST, however, leads to lower values of N eff , due to a smaller tension between Planck and H071p6. In particular, with Planck + BK14 + BAO + H07106, we get a bound of N eff = 3.202 +0.200 −0.202 (68%). Thus, N eff = 3.045 is no longer excluded at 68% in this case.
The SNe Ia luminosity distance measurements provide information about evolution of luminosity distance as a function of redshift (0.01 < z < 2.3 for the Pantheon sample). This can be used to measure the evolution of the scale factor [92] and is helpful in constraining the dark energy EoS. We found that addition of the PAN data did help in constraining the dark energy parameters more tightly. For Planck+BK14+BAO, we have a bound of w 0 < −0.859 (95% C.L.), which shrinks to w 0 < −0.933 (95% C.L.) with the addition of PAN. On the other hand, Planck+BK14+BAO produces a bound of w a = 0.013 +0.065 −0.077 (68% C.L.), whereas Planck+BK14+BAO+PAN leads to w a = 0.033 +0.036 −0.063 (68% C.L.). We see that the 68% spreads of w a have shrunk. This has also been depicted in figure 4. The HST prior also has similar but less strong effect. With Planck+BK14+BAO+HST we have w 0 < −0.908 (95% C.L.) and w a = 0.028 +0.046 −0.065 (68% C.L.).  In all cases we found that the cosmology is compatible with a cosmological constant (i.e., w 0 = −1, w a = 0). As far as values of the tensor-to-scalar ratio is concerned, we find that if we run the chains without the BK14 data, we get a bound of r 0.05 < 0.155 (95% C.L.) with Planck+BAO, which is higher than the bound of r 0.05 < 0.12 (95% C.L.) set by Planck collaboration [83]. However, inclusion of the BK14 data leads to a bound of r 0.05 < 0.075 (95% C.L.), which is close to the r 0.05 < 0.07 (95% C.L.) limit set by the BICEP2/Keck collaboration [89]. The value of r 0.05 remains almost unchanged across all the datasets as long as the BK14 data is included.

NPDDE11 model
In this section we consider the NPDDE11 model where we turn off the tensor perturbations and also do not include the BK14 data. This does not affect the bounds much as can be seen from table 4 and comparing with table 2, which verifies the stability of the results in a smaller parameter space.
The 1-D posteriors for m ν and N eff for selected datasets are given in figure 5. We again find strong bounds on the sum of neutrino masses. We notice that the removal of BK14 data has a small effect on m ν which persists over different datasets. For instance, in NPDDE11+r, for Planck+BAO, we find a m ν < 0.131 eV (95% C.L.), which is reduced to m ν < 0.123 eV (95% C.L.) when we add the BK14 data. In the NPDDE11, this bound is m ν < 0.126 eV (95% C.L.) with Planck+BAO, which is our best bound in this model. This is also stronger than the bound obtained in ΛCDM + m ν with Planck+BAO, as in the previous NPDDE11+r model, and a large improvement compared to the ones presented in [35], which varied dark energy parameters in both in phantom and non-phantom range.
The strengthening of the bound from NPDDE11+r to NPDDE11 with Planck+BAO might simply be due to reduction in the parameter space volume. On the other hand it seems BK14 prefers a lower m ν . However even then the changes are small. BK14 data also seems to prefer slightly larger values of σ 8 , thereby increasing the tension with CFHTLenS. Also, the inclusion of HST prior again seems to discard the standard value of N eff = 3.045 at 68% C.L. but again, not at 95% C.L., and also it doesn't lead to stronger m ν , as before in the NPDDE+r model, due to a large positive correlation between H 0 and N eff but a only small correlation between H 0 and m ν . This can be visualized in figure 6. The PAN dataset provides stricter bounds on w 0 and w a , as   is no longer excluded at 68% in this model also.

NPDDE11+A lens model
We present the limits on the cosmological parameters in table 5. A number of important changes happen with the introduction of the new varying parameter A lens . Considering that our main goal in this paper is to constrain neutrino masses, we see a substantial relaxation in the bounds on m ν . In previous cases we had fixed A lens = 1. However now that A lens is varied we find that the data prefers a large A lens and discards the ΛCDM value of A lens = 1 at more than 95% C.L.
(except in case of inclusion of Planck lensing data, which prefers a much lower A lens , implying a tension between Planck and lensing). The increasing of the lensing amplitude A lens has the same effect as the decreasing of m ν [93]. Increasing A lens leads to smearing of high-l peaks in the CMB temperature and polarization angular power spectra (C T T l , C T E l , C EE l , C BB l ), due to increased grav-itational lensing. On the other hand, massive neutrinos help in reducing this smearing, because it decreases the gravitational lensing of the CMB photons, by suppressing the matter power spectrum in small scales, due to neutrinos having large thermal velocities which prevents them from clustering. Increasing the m ν parameter causes increasing suppression of matter power in the small scales [64], which leads to decreasing gravitational lensing of the CMB photons. This leads to a strong positive correlation between A lens and m ν , such as, to compensate for the increase in A lens , the neutrino masses are also increased. The 1-D plots for m ν and N eff for selected datasets are given in figure 8. In this model, the Planck data is almost insensitive to neutrino masses < 0.6 eV. Our tightest bound of m ν < 0.239 eV (95% C.L.) again comes with Planck+BAO data. This bound, while weaker than the previous models we have discussed, is still close to the m ν < 0.23 eV (95% C.L.) bound provided by Planck collaboration [83], and still a large improvement compared to the ones presented in [35], which varied dark energy parameters in both in phantom and non-phantom range and had found a bound of m ν < 0.557 eV (95% C.L.) with Planck+BAO, demonstrating the large difference between phantom and non-phantom dark energies as far as neutrino masses are concerned. The preferred N eff values are also higher in NPDDE11+A lens compared to the previous cases. The addition of the HST data leads to even higher N eff which leads to the N eff = 3.045 value being disallowed even at 95% C.L. with Planck+HST, for which the 68% and 95% limits are N eff = 3.517 +0.196 −0.216 and N eff = 3.517 +0.424 −0.396 respectively. This signifies the presence of tension between Planck and HST in this model, as it was in previous models.
The use of the H071p6 prior, again leads to lower values of N eff . In particular, with Planck + BAO + H07106, we get a bound of N eff = 3.329 +0.207 −0.227 (68%). Thus, N eff = 3.045 is not excluded at 95% in this model, but excluded only at 68%.
Another important change is the change in bounds on the optical depth to reionization, τ . With Planck+BAO, the NPDDE11 model preferred a value of τ = 0.092 ± 0.018 (68% C.L.), whereas this model prefers τ = 0.059 +0. 21 −0.22 (68% C.L.), which is actually closer to the bound of τ = 0.055 ± 0.009 (68% C.L.) given by Planck 2016 intermediate results [94]. This was previously observed in [35] which did the analysis with varying the dark energy parameters in both the phantom and nonphantom sector. This implies that the main effect is through the degeneracy between τ and A lens and has not much to do with dark energy. Again, while the NPDDE11+r and NPDDE11 models failed to reconcile Planck with weak lensing measurements like CFHTLenS, the NPDDE11+A lens model prefers lower values of σ 8 and the agreement with CFHTLenS is considerable. This can be visualized in figure 9. This was also previously seen in [35] and hence, again we can infer that this happens because of varying A lens . The bounds on the dynamical dark energy parameters are however weaker than in the other two models. The cosmological constant is however compatible with the data even in this model.

τ and A lens : Implications for Planck 2018
Both τ and A lens are correlated with m ν , and with each other. In particular, when A lens is fixed, increase in m ν reduces smearing in the damping tail of the CMB power spectra, and it can be compensated by increasing τ [10,21]. Hence they have a positive correlation. On the other hand, increasing A lens increases the smearing of the damping tail, i.e., negative correlation with τ . The value of τ has been significantly improved from Planck 2015 to Planck 2018. Thus we consider a bound on this optical depth to reionization, τ = 0.055 ± 0.009, taken from [95], in which Planck col- laboration removed previously unexplained systematic effects in the polarization data of the Planck HFI on large angular scales (low-l). We refer to this prior as τ 0p055 hereafter. We use τ 0p055 as a substitute for low-l polarization data, and thus we discard the lowP data whenever we apply the τ 0p055 prior, to avoid any double counting. This prior is very close to the bound, τ = 0.0544 +0.0070 −0.0081 (68%), obtained with Planck 2018 temperature and polarization data [1]. Hence, imposition of τ = 0.055 ± 0.009 would produce bounds on m ν that will be close to the bounds produced with Planck 2018 (instead of Planck 2015) in the models that we have considered.
We find that in the NPDDE11+r model, with Planck + BK14 + BAO + τ 0p055, we get m ν < 0.097 eV (95%) (i.e. improvement over the m ν <0.123 eV limit as in 2, with Planck + BAO). This bound is actually lower than the m ν 0.1 eV , i.e. minimum mass required for inverted mass hierarchy of neutrinos. At the same time, in the NPDDE11 model, with Planck + BAO + τ 0p055 we get m ν < 0.107 eV (95%), which is also an improvement from the result: m ν < 0.126 eV (95%) with Planck + BAO (see table 5. This happens, since in both of these models the mean value of τ hovers around 0.09-0.1. The τ 0p055 prior partially breaks the degeneracy between τ and m ν , and produces lower values of m ν by lowering the preferred τ values. On the other hand, in the NPDDE11+A lens model with Planck + BAO + τ 0p055, we found m ν < 0.237 eV, which is almost similar to the bound m ν < 0.239 eV (95%) with Planck + BAO (see table 3). This happens since all the three parameters, τ , A lens , and m ν are varied together. Now, as the data prefers A lens values higher than the ΛCDM value in this model, the degeneracy between A lens and τ leads to a much lowered value of τ , and thus the correlation between τ and m ν is already much smaller in this model, than the other two. Thus τ 0p055 has little effect on the neutrino mass bounds in this model. Also, we obtained limits of A lens in a ΛCDM + A lens model with Planck 2015 full temperature and polarization data. The value we got is A lens = 1.15 +0.072 −0.082 (68% C.L.). In the Planck 2018 Cosmological Parameters paper [1], for similar data and same model, given value of A lens is: A lens = 1.18 ± 0.065 (68%) (see equation 36b). This shows that there is only a very small change in A lens from Planck 2015 to Planck 2018. Thus, it is likely that there will not be any considerable changes in the limits of other cosmological parameters with the Planck 2018 data, in the context of the value of A lens .

Summary
In this work we have studied three different extended cosmological scenarios with non-phantom dynamical dark energy (NPDDE) with a focus on constraining sum of neutrino masses. We have presented bounds on all the varying parameters in these extended scenarios and described the main effects we observed. In the first model, NPDDE11+r, we consider 12 parameters: the 6 ΛCDM parameters, two dynamical dark energy parameters with CPL parametrization (w 0 and w a ) with hard priors to satisfy the non-phantom requirement, number of effective relativistic neutrino species at recombination (N eff and sum of neutrino masses ( m ν ), and the running of the inflation spectral index (n run ) and the tensor-to-scalar ratio (r 0.05 ). We used different combinations of recent datasets including Planck 2015 temperature and polarization data, CMB B-mode spectrum data from BICEP2/Keck collaboration (BK14), BAO SDSS III BOSS DR12, MGS and 6dFS data, SNe Ia Pantheon sample (PAN), the HST prior (H 0 = 73.24 ± 1.74 km/sec/Mpc (68% C.L.)). We found that CMB only data is not very effective in constraining the cosmological parameters. The 1σ spreads for the parameters were however increased in this model compared to ΛCDM due to the doubling of number of parameters. Our best bound on neutrino masses in this model came from Planck+BK14+BAO: m ν < 0.123 eV (95% C.L.) which is a strong bound close to the minimum mass of 0.1 eV (95% C.L.) required for inverted hierarchy of neutrino masses and is stronger than a bound of m ν < 0.158 eV (95% C.L.) obtained in ΛCDM + m ν with Planck+BAO [21] (see also [54] for a similar conclusion in a smaller parameter space). We also found that inclusion of the HST prior leads to a preference for dark radiation at 68% C.L. but not at 95%, while without the HST prior the data is completely consistent with the standard value of N eff = 3.045. Although this is driven by the more than 3σ tension present between Planck and HST regarding the value of H 0 and should be interpreted cautiously. This model did not improve the σ 8 tension present in the σ 8 − Ω m plane between Planck and CFHTLenS. The Pantheon sample improved the bounds on the dark energy parameters. All combinations of data are also compatible with a cosmological constant (w 0 = −1, w a = 0). However, this is mostly because we are restricting the parameter space to w(z) ≥ −1 and [35] had found that the data mostly prefers the phantom region in such an extended parameter space when both phantom and non-phantom regions are allowed.
We tested the stability of these results in a lower parameter space (model:NPDDE11) where we turned off the tensor perturbations and also did not use the BK14 data. We found that the general conclusions made for NPDDE11+r were also true in this model. The tightest bound of m ν < 0.126 eV (95% C.L.) in this model also came from Planck + BAO.
Finally we studied the NPDDE11+A lens model where we also varied the lensing amplitude. We found that except when Planck lensing data is included, the A lens = 1 value predicted by ΛCDM was rejected at more than 95% C.L. by the datasets. Due to this, the m ν bounds also worsened with our best result in this model: m ν < 0.239 eV (95% C.L.) coming from Planck+BAO again. This result is, however, still close to the m ν < 0.23 eV (95% C.L.) bound by Planck collaboration [83], showing that the cosmological data is effective in constraining neutrino masses in a cosmology with NPDDE. The HST prior also preferred a dark radiation component but this time also at 95% C.L. level, as this model also prefers higher values of N eff . On the other hand, we found that this model helps relieve the σ 8 tension between Planck and CFHTLenS considerably.
The recent Planck 2018 results [1] put the bound of m ν < 0.13 eV (95% C.L.) in ΛCDM + m ν with Planck+BAO. Thus, the aggressive bound of m ν < 0.123 eV (95% C.L.) (Planck + BK14 + BAO) is still stronger than this bound by Planck 2018 and hence, our results are very much relevant albeit the analysis is with Planck 2015 dataset. In fact, when we use the following Gaussian prior on optical depth to reionization: τ = 0.055 ± 0.009 from 2016 Planck intermediate results, and discard the low-l polarization data, the bound on neutrino masses improves to m ν < 0.097 eV (95%), which is less than the 0.1 eV mass sum required for inverted hierarchy of active neutrino masses.
While we have used the CPL parameterization in our paper, it is not the only parameterization that can be used for non-phantom dark energy. Any change in parameterization can lead to change in bounds obtained on the sum of neutrino masses. For instance, if we set the w a parameter to zero, i.e., if we consider only a simple w(z) = w 0 parameterization, we find that bounds on m ν relax slightly. In the NPDDE11 model, with w a = 0 and w(z) = w 0 , and using Planck + BAO data, we found m ν < 0.141 eV (95%), instead of m ν < 0.126 eV (95%) when we vary both w 0 and w a . In the NPDDE11+A lens model also, with w a = 0 and w(z) = w 0 , and using Planck + BAO, we obtained m ν < 0.261 eV (95%), instead of m ν < 0.239 eV (95%). Some other parameterizations that can be considered include Logarithmic parameterization [96] (w(a) = w 0 − w a ln(a)), Jassal-Bagla-Padmanabhan (JBP) parameterization [97] (w(a) = w 0 + w a a(1 − a) etc. Analysis involving these parameterizations is beyond the scope of our work in this paper. However, we would like to point the reader to [41], where the authors found similar limits, with CPL and Logarithmic parameterizations, on m ν for the case of degenerate hierarchy. However, in case of JBP, bound on m ν was found to be significantly stronger. While [41] does not discard the phantom region, it is possible that results from analyses with only non-phantom dark energy will also vary depending on the parameterization used, as far as neutrino masses are concerned.
We would like to add a final remark that we have obtained the bounds while taking the datasets at face value. However unresolved systematics present in the dataset could have affected our results and conclusions. For instance the tension between Planck and HST prior can be due to a dark radiation species, but can also be due to systematics present in both the datasets. Thus there is still a lot to learn about robustness of datasets and also about dynamics of dark energy.
Union's Horizon 2020 research and innovation programme Elusives ITN under the Marie Sklodowska-Curie grant agreement No 674896.