Galaxy clustering, CMB and supernova data constraints on $\phi$CDM model with massive neutrinos

We investigate a scalar field dark energy model (i.e., $\phi$CDM model) with massive neutrinos, where the scalar field possesses an inverse power-law potential, i.e., $V(\phi)\propto {\phi}^{-\alpha}$ ($\alpha>0$). We find that the sum of neutrino masses $\Sigma m_{\nu}$ has significant impacts on the CMB temperature power spectrum and on the matter power spectrum. In addition, the parameter $\alpha$ also has slight impacts on the spectra. A joint sample, including CMB data from Planck 2013 and WMAP9, galaxy clustering data from WiggleZ and BOSS DR11, and JLA compilation of Type Ia supernova observations, is adopted to confine the parameters. Within the context of the $\phi$CDM model under consideration, the joint sample determines the cosmological parameters to high precision. It turns out that $\alpha<4.995$ at 95% CL for the $\phi$CDM model. And yet, the $\Lambda$CDM scenario corresponding to $\alpha = 0$ is not ruled out at 95% CL. Moreover, we get $\Sigma m_{\nu}<0.262$ eV at 95% CL for the $\phi$CDM model, while the corresponding one for the $\Lambda$CDM model is $\Sigma m_{\nu}<0.293$ eV. The allowed scale of $\Sigma m_\nu$ in the $\phi$CDM model is a bit smaller than that in the $\Lambda$CDM model. It is consistent with the qualitative analysis, which reveals that the increases of $\alpha$ and $\Sigma m_\nu$ both can result in the suppression of the matter power spectrum. As a consequence, when $\alpha$ is larger, in order to avoid suppressing the matter power spectrum too much, the value of $\Sigma m_\nu$ should be smaller.


I. INTRODUCTION
Neutrino is one of the important bonds linking nuclear physics, particle physics, astrophysics and cosmology [1]. In the Standard Model (SM) of particle physics, it is anticipated that there are three types, or "flavors", of neutrinos: electron neutrino (ν e ), muon neutrino (ν µ ) and tau neutrino (ν τ ), which are also dubbed as three normal/active neutrinos. Besides that, neutrinos are assumed to be massless in the SM of particle physics [2].
It was first predicted by Bruno Pontecorvo in 1957 that if neutrinos are massive the neutrino flavor should be unstable, that is called neutrino (flavor) oscillations [3]. Briefly put, neutrino oscillation is a phenomenon that a neutrino produced in a definite flavor is observed in a different flavor after traveling some distances. In other words, neutrinos are able to oscillate among the three available flavors while they propagate through space. Nowadays there are compelling evidences for neutrino oscillations from a variety of experimental data on solar, atmospheric, reactor and accelerator neutrinos. The discovery of neutrino oscillations implies that neutrinos have small but non-zero masses, with at least two species being non-relativistic today. However, the present experimental results on neutrino oscillations only measure the difference of two squared masses, such as ∆m 2 21 = m 2 2 − m 2 1 and ∆m 2 32 = m 2 3 − m 2 2 , but give no hint on their absolute mass scales. m 1 , m 2 and m 3 are the neutrino mass eigenstates. For example, the solar neutrino analysis supplemented by KamLAND produces an estimate of ∆m 2 21 ∼ 8×10 −5 eV 2 [4], and the measurement of atmospheric neutrino oscillation by Super-Kamiokande I indicates ∆m 2 32 ∼ 3 × 10 −3 eV 2 [5]. If it is the case of oscillations among three light neutrinos, only two of the three ∆m 2 ij are independent, as ∆m 2 21 + ∆m 2 32 + ∆m 2 13 = 0, where ∆m 2 13 = m 2 1 − m 2 3 . Recent reviews on progress in both theoretical and experimental aspects of neutrino oscillations can be found in [6].
A variety of cosmological tests are sensitive to the absolute scale of neutrino mass, such as the cosmic microwave background (CMB) radiation, galaxy surveys, and the Lyman-alpha forest [7]. In [8], the effect of massive neutrinos on the Sunyaev-Zel'dovich and X-ray observables of galaxy clusters are investigated with a set of six very large cosmological simulations (8h −3 Gpc 3 comoving volume). The analysis of current cosmological observations provides an upper bound on the total neutrino mass m ν (summed over the three neutrino families) of order 1 eV or less. However, the limits on m ν from cosmology are rather model dependent and vary strongly with the data combination adopted. For example, in the framework of one-parameter extensions to the base ΛCDM model the Planck 2015 results [9] give 95% upper limits on the sum of neutrino masses, i.e., m ν < 0.23 eV for a combination of Planck TT+lowP+lensing+ext, and m ν < 0.59 eV for Planck TT, TE, EE+lowP+lensing, where "TT" denotes the combination of the TT likelihood at multipoles l ≥ 30 and a low-l temperature-only likelihood, "TE" denotes the likelihood at l ≥ 30 using TE spectra, and "EE"denotes the likelihood at l ≥ 30 using EE spectra,"lowP" denotes the low-l Planck polarization data, "lensing" is the Planck lensing data, and "ext" represents the external data including the baryon acoustic oscillations (BAO), Type Ia supernovae (SNe Ia), and H 0 . In [10], the power law and exponential types of viable f (R) theories along with massive neutrinos are studied. It shows that the allowed scales of m ν in the viable f (R) models are greater than that in the ΛCDM model. The cases of fixing the effective number of neutrino species as N eff = 3.046 and treating N eff as a free parameters are both considered in [10]. The former corresponds to just consider the active neutrinos without the effect of dark radiation. The latter corresponds to include the contribution of dark radiation (represented by ∆N eff = N eff − 3.046). For more details on dark radiation, we refer the reader to [11]. The model of holographic dark energy with massive neutrinos and/or dark radiation is investigated in [12], but the computed results from this model are not compared with those from the ΛCDM model. Actually, the ΛCDM model with massive neutrinos is discussed broadly with constraints from various cosmological observations [13]. The time evolving of neutrino mass is also explored in the literature [14]. For further details on neutrino cosmology, the reader is referred to recent reviews such as [7,15].
In this paper, we will discuss the constraints on the sum of neutrino masses m ν in the framework of φCDM model by using a combination of the CMB data from Planck 2013 and WMAP9, the galaxy clustering data from WiggleZ and BOSS surveys, and the JLA compilation of SNe Ia observations. The effect of dark radiation is not considered in this work, i.e., N eff = 3.046. We also assume that one of the three active neutrinos is massive, and the other two are massless. The φCDM model -in which dark energy is modeled as a scalar field φ with a gradually decreasing (in φ) potential V (φ) -is a simple dynamical model with a slowly decreasing (in time) dark energy density. This model could resolve some of the puzzles of the ΛCDM model [16], such as the coincidence and fine-tuning problems. Here we focus on the scalar field with an inverse power-law potential V (φ) ∝ φ −α , where α is a nonnegative constant [17,18]. When α = 0 the φCDM model is reduced to the corresponding ΛCDM case. The φCDM model with this kind of V (φ) has been extensively investigated [19][20][21], but without considering the massive neutrinos.
The rest of the paper is organized as follows. In Sec. II we present the background and perturbation evolutions of the φCDM model with massive neutrinos. The impacts of m ν and α on the CMB temperature power spectrum and on the matter power spectrum are also discussed. Constraints from the cosmological data are derived in Sec. III, and the results for φCDM model are compared with those for the ΛCDM model. We summarize our main conclusions in Sec. IV.

II. THE φCDM MODEL WITH MASSIVE NEUTRINOS
A. Background evolution of the φCDM model Quintessence as one of the popular scalar field dark energy models is a hypothetical form of dynamical dark energy to explain the late-time cosmic acceleration. Since quintessence is described by the scalar field φ, the corresponding dark energy model can also be called as φCDM model. In what follows, we will use the terms "quintessence" and "φCDM" essentially interchangeably. We consider the self-interacting scalar field φ minimally coupled to gravity on cosmological scales. The action of this φCDM model is given by where g is the determinant of the metric g µν , R is the Ricci scalar, m p = 1/ √ G is the Planck mass with G being the Newtonian constant of gravitation, L is the Lagrangian density for matter and radiation, and L φ is the Lagrangian density for the field φ, given by where V (φ) is the field's potential. In this work, we take a flat Friedmann-Lemaître-Robertson-Walker (FLRW) metric for the background evolution, which is described by where x i is the comoving coordinate. a(t) is the scale factor usually normalized to unity now a 0 = a(z = 0) = 1 and related to the redshift z as a/a 0 = 1/(1 + z). Throughout, the subscript "0" denotes the value of a quantity today. By the variation of the action in Eq. (1) with respect to φ, one can obtain the Klein-Gordon equation (equation of motion) for the scalar fieldφ For the φCDM model, there are many kinds of V (φ) which can satisfy the requirement of the late-time accelerating expansion of the universe [22]. In 1988, Peebles and Ratra [18] proposed a scalar field that is slowly rolling down with a potential V (φ) = 1 2 κm 2 p φ −α at a large φ, where κ and α are nonnegative parameters. This inverse power-law potential can not only lead to the late-time acceleration of the universe but also partially solve the cosmological constant problems. The larger value of α induces the stronger time dependence of the scalar field energy density ρ φ . When α = 0, this φCDM model is reduced to the ΛCDM case. What is more, the parameter κ depends on α (see [19,23] for its dependence on α).
The Friedmann equation of the φCDM model with massive neutrinos can be written as where ρ b , ρ c , ρ φ , ρ γ and ρ ν denote the energy densities of baryons, cold dark matter (CDM), scalar field dark energy, photons and neutrinos, and H(z) ≡ȧ/a is the Hubble parameter. The energy density and pressure of the scalar field dark energy are given by and Then, one can work out the equation of state (EoS) of the field φ, which is clearly bounded in the range −1 < ω φ < 1 and usually non-constant. One can see that if the scalar field φ rolls slowly enough such that the kinetic energy density is much less than the potential energy density, i.e.φ 2 ≪ V (φ), the pressure P φ of the field will become negative with ω φ → −1.
Based on Eqs. (4) and (5), along with the initial conditions described in Refs. [18,19], one can numerically compute the Hubble parameter H(z). We also introduce the dimensionless density parameter for each component as, Ω X = ρ X /ρ cr , where the index "X" denotes the individual components, such as radiation ("r"), neutrino ("ν") and matter ("m"). The critical energy density is expressed as ρ cr = 3H 2 m 2 p /(8π). Ω m is the energy density of matter including both baryons and CDM. Ω ν is the total neutrino energy density which scales as ∝ a −4 at early times, and thereafter evolves as ∝ a −3 after the non-relativistic transition. One can see that the massive neutrinos behave like the radiation at early times and like the matter later.

B. Cosmological perturbation of the φCDM model
Let us consider perturbations of the flat FLRW metric in the Newtonian Gauge [24]. In this gauge, the linear perturbed metric is given by where the scalar perturbations are dominant over vector or tensor perturbations, and η = a −1 dt is the conformal time. The Newtonian force Ψ gives rise to the dynamics of the perturbed fluids, while the curvature perturbation Φ measures the local energy density fluctuations. The linear perturbation theory is a good tool both for describing the early universe at any scales, and the recent universe on the largest scales. It has been shown in [25] that for self-interacting scalar field dark energy models it is phenomenologically sufficient to regard the dark energy component as a perfect fluid. We treat each component in the universe as perfect fluid, including the baryon, CDM, photon, neutrino and scalar field dark energy. In the perfect fluid approach, the perturbed Einstein equations lead to the following Eqs. (10) - (13) in the Fourier space: and The great advantage of linear theory is to obtain independent equations of evolution for each Fourier mode. All of the perturbed quantities (δ X , θ X , Ψ, Φ, etc.) are functions of space x and time t, where X denotes each perfect fluid composing the universe. In the linear perturbations, the anisotropic stress σ X is negligible for the perfect fluids. Note that a prime represents a derivative with respect to the conformal time η. The spatial variation of density fluctuations is expressed by the density contrast δ X ≡ δρ X /ρ X = (ρ X −ρ X )/ρ X , andρ X is the background energy density of component X. In the approximation of negligible irrotational flow, the divergence of the peculiar velocity v X , θ X = ∇ · v X can be used to describe the fluid motion. In the Fourier space, we have θ X ≡ i k · v X . While ω X ≡P X /ρ X is the equation of state of each component, and c 2 s,X ≡ δP X /δρ X represents the sound velocity. Eq. (10) is called as the (perturbed) continuity equation, that states the conservation of local density. Eq. (11) is called as the Euler equation, that represents the conservation of local energy momentum, and describes dynamics of perturbed fluids originated by the Newtonian force Ψ. The curvature perturbation Φ is constrained to the local inhomogeneity via the Poisson equation Eq. (12). We can get Eq. (13) under the assumption that the perturbed fluid remains a perfect fluid. These equations (10) -(13) completely determine the dynamical evolution of large scale structure (LSS) of the universe, within a given expansion history H.

C. Matter power spectrum and CMB power spectrum in the φCDM model
In the framework of φCDM model, we qualitatively investigate the impacts of parameters α and Σm ν on the matter power spectrum and on the CMB power spectrum. The analyses are performed with the CAMB Boltzmann code [26].
Neutrinos rarely interact with matter after thermal decoupling, so they are treated as free streaming particles. Massive neutrinos are the only particles that present the transition from radiation to matter. Before the nonrelativistic transition the neutrinos behave like radiation. Thus, when the neutrino mass Σm ν increases, the time of radiation/matter equality is postponed gradually, and a eq increases. The value of Σm ν can affect the matter power spectrum and the CMB power spectrum mainly resulting from a change in the time of equality, that provides a potential way to constrain it through CMB and LSS observations [7,12,27]. In Fig. 1, we show the impacts of neutrino mass Σm ν on the matter power spectrum P (k) and the CMB temperature anisotropy spectrum C T T l . The upper panels show the cases for φCDM model with varying values of Σm ν , where α is fixed as α = 1, and other parameters are fixed based on the recent Planck results [9]. For comparison, we also display the cases for ΛCDM model in the lower panels. For both φCDM and ΛCDM models, the matter power spectrum is gradually suppressed with the increase of Σm ν , however, the effect is more significant on small scales than that on large scales. One possible reason is that the neutrino perturbations do not contribute to gravitational clustering on scales smaller than the free-streaming scale, while on the very large scales neutrino perturbations are never affected by free streaming, and they become indistinguishable from CDM perturbations in the non-relativistic regime [7]. The CMB temperature anisotropy spectrum C T T l is insensitive to the variation of Σm ν in both ΛCDM and φCDM models. The parameter α indicates the dynamics of dark energy, and then it can affect the expansion history of the universe and the redshift of matter/dark energy equality. When α increases, the expansion of the universe occurs more rapidly, and the epoch of dark energy domination begins earlier [19]. For these reasons, the variation of α can have signatures in the CMB map and the matter clustering. The impacts of α on P (k) and C T T l are presented in Fig. 2. We choose α = 0, 1, and 10 as examples, where α = 0 corresponds to the ΛCDM scenario. The values of other parameters are kept fixed. We find that P (k) is slightly suppressed with the increase of α, and the effect is a bit significant on large scales than that on small scales. The CMB temperature anisotropy spectrum C T T l is a little sensitive to the variation of α on the low-l tail, which may arise from the late Integrated Sachs-Wolfe (ISW) effect. Anyhow, CMB and LSS observations are efficient to distinguish between ΛCDM and φCDM models.

III. OBSERVATIONAL CONSTRAINTS
The observational data sets used to constrain the cosmological parameters are described as follows, including the galaxy clustering, CMB and SNe Ia measurements.
A. Cosmological data sets

Galaxy clustering measurements
Galaxy clustering distilled from the galaxy redshift survey is powerful as cosmological probe [28], that can allow us to measure the cosmic expansion history through the measurement of BAO, and the growth history of cosmic large scale structure through measurements of redshift-space distortions. The length scale of BAO, the comoving sound horizon at the baryon drag epoch r s (z d ), can be applied as a cosmological ruler and accurately calibrated by observations of the CMB radiation. The position of the BAO peak in the angle-averaged galaxy clustering pattern is usually quantified in term of the volume averaged distance [29] It is common to report the BAO distance measurements as combinations of the angular diameter distance, D A (z), and the Hubble parameter, H(z), such as or The redshifts of galaxies include indistinguishable contributions from both the Hubble recession and the peculiar velocity of the galaxies themselves, so that there are errors in the distances we assign to galaxies. The differences between the redshift-inferred distances and true distances are known as redshift-space distortions (RSD) [30]. In another word, the RSD are introduced in the observed clustering pattern by galaxy peculiar motions. As a consequence, the correlation function and the power spectrum measured in the redshift space are different from those in the real space, which have to be corrected to be expressed in real space. Because the effects of RSD couple the density and velocity fields, the RSD signals within the correlation function are difficult to model. On different scales, peculiar motions produce different types of distortions to the power spectrum. On small scales, i.e., in the cluster cores, the peculiar velocities of galaxies are almost randomly oriented, that cause the structures to appear elongated along the line of sight (LOS) when viewed in redshift space (i.e., the so called "finger of God" effect) [31], leading to a damping of the clustering. On large scales, because of gravitational growth, the galaxies tend to fall towards high-density regions, and flow away from low-density regions, such that the galaxy clustering in redshift space is enhanced in the LOS direction compared to the transverse direction [32]. The RSD effect on large scales can be described by linear theory [32,33], while the "finger of God" effect is a non-linear phenomenon. On large scales where the gravitational growth is linear, measuring the relative clustering in both LOS and transverse directions leads to measurements of the parameter combination f (z eff )σ 8 (z eff ), where z eff is the effective redshift. f is the growth rate of cosmic structure, which is associated with the evolution of matter density perturbations δ m via the relation f ≡ d ln δ m /d ln a. In the linear regime, the linear growth rate can be expressed as f = d ln D(z)/d ln a, where D(z) = δ m (z)/δ m (z = 0) is the linear growth factor normalized such that D(z = 0) = 1. σ 8 (z) is the root-mean-square amplitude of the matter fluctuations in spheres of 8h −1 Mpc, and σ 8 (z = 0)/σ 8 (z) = D(z = 0)/D(z). Thus, one can figure out where σ 0 8 = σ 8 (z = 0). In linear theory, the galaxy bias b and the growth rate f are degenerate with σ 8 , so the RSD measurements are better presented in terms of b(z)σ 8 (z) and f (z)σ 8 (z), rather than f (z). Currently, the biasindependent parameter combination f (z)σ 8 (z) measured by RSD are widely used [20,34,35].
The Alcock-Paczynski (AP) test [36] is proved to be a significant link between BAO and RSD. AP test states that if an astrophysical structure is spherically symmetric or isotropic, then it should possess equal comoving transverse and radial sizes. An AP measurement is carried out by comparing the observed transverse and radial dimensions of objects. While the AP test is equally valid for an isotropic process such as the two-point statistics of galaxy clustering. The apparent anisotropy of the two-dimensional correlation function of galaxies mainly arises from the geometry and expansion of the universe which should be correctly embodied in the fiducial cosmological model, and the RSD effect which is supposed to be marginalized by using appropriate RSD model (see [30] for a recent review of RSD models). According to the requirement of AP test, the signature of BAO should have identical comoving sizes (i.e., r s ) in transverse and radial dimensions. The observed transverse dimension is the angular projection ∆θ = r s /[(1 + z)D A (z)]. The radial one is the redshift projection ∆z = r s H(z)/c. The relative radial/transverse distortion depends on the value of where F (z) is dubbed as the AP distortion parameter. By combining the BAO peak, AP test and RSD effect, one can report the galaxy clustering effectively as joint measurements of (A, F, f σ 8 ) or (d z , F, f σ 8 ). These joint measurements are extremely good at helping to constrain basic cosmological parameters and distinguish between the dark energy models. By using large-scale structure measurements from the WiggleZ Dark Energy Survey [37], Blake et al. (2012) [38] have performed joint constraints of (A, F, f σ 8 ) in three overlapping redshift slices with effective redshifts z eff = (0.44, 0.6, 0.73). Utilizing these data, it is straightforward to put constraints on the model parameters by calculating the corresponding χ 2 WiggleZ , given by The observational data vector is  Table 1 of [38]. The vector of theoretical values is where [z 1 , z 2 , z 3 ] = [0.44, 0.6, 0.73], and the corresponding theoretical values of (A, F, f σ 8 ) can be obtained with Eqs. (15), (18) and (17), respectively. C is a 9 × 9 covariance matrix between parameters and redshift slices, and the value of 10 3 C is listed in Table 2 of [38], that is achieved by generating 400 lognormal realizations for each WiggleZ survey region and redshift slice with the methods described in [34]. In the analysis of [38], the fitting formulae provided by Jennings et al. (2011) [39] have been taken as the fiducial RSD model. The effect of different choices of the RSD model is also considered in Section 3.4 of [38]. It turns out that the systematic error induced from modeling RSD is much lower than the statistical error in the measurement. Joint measurements of (d z , F, f σ 8 ) at an effective redshift of z eff = 0.57 are provided in   [40] by utilizing the observed anisotropic clustering of galaxies in the Baryon Oscillation Spectroscopic Survey (BOSS) Data Release 11 (DR11) CMASS sample [41]. We employ this data set in our analysis with the chi-squared statistic The observational data vector is Y obs = [d z , F, f σ 8 ], i.e., Y obs = [13.85, 0.6725, 0.4412] by using the mean values presented in Eq. (30) of [40].

The vector of theoretical values is
where the corresponding theoretical values of (d z , F, f σ 8 ) can be obtained with Eqs. (16), (18) and (17), respectively. The covariance matrix Cov of measurements is listed in Eq. (31) of [40]. A suite of 600 PTHalo simulations are used to estimate the covariance matrix (see [42] for details of mock generation). In the analysis of [40], the "streaming model"-based approach developed in [43] has been adopted to model the RSD signal, that has been demonstrated to fit the monopole and quadrupole of the galaxy correlation function with better than percent level precision to scales above 25h −1 Mpc, for galaxies with bias of b ≃ 2.
The galaxy clustering (GC) measurements from WiggleZ and BOSS DR11 are both employed in this study. Thus, the corresponding chi-squared statistic is expressed as where χ 2 WiggleZ and χ 2 BOSS are given by Eqs. (19) and (22), respectively.

CMB power spectrum measurements
The CMB radiation deemed as the afterglow of the big bang can supply us with some information of the very early universe. The observations of CMB provide another independent test for the existence of dark energy. The recent precise measurements of the CMB radiation from Planck and Wilkinson Microwave Anisotropy Probe (WMAP) projects can efficiently improve the accuracy of constraining the cosmological parameters. Currently, the Planck 2015 results have come out [9], but the Planck 2015 likelihoods are not yet available. Given this, we use the low multipoles (2 ≤ l ≤ 49) and high multipoles ( 50 ≤ l ≤ 2479) temperature power spectrum likelihoods from Planck 2013 [44], together with the low multipoles (l ≤ 23) polarization power spectrum likelihoods from nine-year WMAP (WMAP9) [45]. To employ the previously mentioned CMB power spectrum data in the analysis, we compute the χ 2 CMB statistic where C obs l is the observational value of the related power spectrum, C th l is the corresponding theoretical value in the framework of the cosmological model under consideration, and M is the covariance matrix for the best-fit data spectrum.

3.
Magnitude-redshift measurements of Type Ia supernovae The first direct evidence for the cosmic acceleration came from SNe Ia observations, which provide the measurement of the cosmic expansion history through the measured luminosity distance as a function of redshift, d L (z) = (1+z)r(z). In the spatially flat universe, the comoving distance r(z) from the observer to redshift z is given by wherein p denotes the parameter space of the considered cosmological model, and E = H/H 0 is the dimensionless Hubble parameter. Here, we use the "joint light-curve analysis" (JLA) compilation of SNe Ia [46], which is a joint analysis of SNe Ia observations including several low-redshift samples (z < 0.1), all three seasons from the SDSS-II (0.05 < z < 0.4), three years from SNLS (0.2 < z < 1), and 14 very high redshift (0.7 < z < 1.4) from the Hubble Space Telescope (HST) observations. It totals 740 spectroscopically confirmed SNe Ia with high quality light curves. In Ref [46], SALT2 light-curve model [47,48] have been used to fit the supernova light curves of the JLA sample. From the observational point of view, the distance modulus of a SN Ia can be yielded from its light curve with an empirical linear relation: The light-curve parameters (m ⋆ B , X 1 , C) result from the fit of the SALT2 light-curve model to the photometric data, where m ⋆ B corresponds to the observed peak magnitude in rest-frame B band, X 1 describes the time stretching of the light-curve, and C describes the supernova color at maximum brightness. α, β and M B are nuisance parameters in the distance estimate, which are estimated simultaneously with the cosmological parameters and then marginalized over when obtaining the parameters of interest, wherein M B is the absolute B-band magnitude. The theoretical (predicted) distance modulus is µ th (z; p, µ 0 ) = 5 log 10 [D L (z; p)] + µ 0 , where µ 0 = 42.38 − 5 log 10 h, which is also treated as a nuisance parameter, and the Hubble-free luminosity distance is given by .
The best-fit cosmological parameters from SNe Ia data are determined by minimizing where Cov is the covariance matrix of data vector µ obs B . The values of the covariance matrix Cov and the SALT2 fit parameters (m ⋆ B , X 1 , C) are available from Ref. [46].

B. Results and analysis
In our analysis, the likelihood is assumed to be Gaussian, thus we have the total likelihood where χ 2 tot is constructed as wherein χ 2 GC , χ 2 CMB and χ 2 SNe are given by Eqs. (23), (24) and (29), respectively, and denote the contributions from galaxy clustering, CMB and SNe Ia data sets described above.
We derive the posterior probability distributions of parameters with Markov Chain Monte Carlo (MCMC) exploration using the February 2015 version of CosmoMC [49]. The parameter space of the ΛCDM model is where Ω b h 2 and Ω c h 2 , respectively, stand for the baryon and CDM densities today, θ MC is an approximation to θ * = r s (z * )/D A (z * ) (i.e., the angular size of the sound horizon at the time of decoupling) that is used in CosmoMC and is based on fitting formulae given in [50], τ refers to the Thomson scattering optical depth due to reionization, n s and A s are the power-law index and amplitude of the power-law scalar primordial power spectrum of curvature perturbations, and Σm ν is the sum of neutrino masses. The parameter space of the φCDM model is which has one more parameter than that of ΛCDM model, where α determines the steepness of the scalar field potential in the framework of φCDM model. The one-dimensional (1D) probability distributions and two-dimensional (2D) contours for the cosmological parameters of interest are shown in Fig. 3 for ΛCDM model and in Fig. 4 for φCDM model. It shows that constraints from the joint sample are quite restrictive, though there are degeneracies between some parameters, such as the degeneracies in the Ω m − H 0 and σ 8 − Σm ν planes. In addition, the differences between the marginalized likelihoods and the mean likelihoods are modest in 1D and 2D plots. It implies that the distributions of the parameters are almost Gaussian. We also present best-fit values and mean values with 95% confidence limits for the parameters of interest in Table I both for ΛCDM and φCDM models. We find α < 4.995 at 95% CL for the φCDM model, while the ΛCDM scenario corresponding to α = 0 is not ruled out at this confidence level. The constraints on Ω b h 2 , Ω c h 2 , 100θ MC , τ , ln(10 10 A s ), n s , Ω m , σ 8 and H 0 are consistent at 95% CL for these two models.
Here, we pay attention to the constraints on Σm ν . Note that Σm ν is in unit of eV. The best-fit vale is Σm ν = 0.038(0.043) with Σm ν < 0.262(0.293) at 95% CL in the framework of φCDM (ΛCDM) model. The allowed neutrino mass scale in the φCDM model is a bit smaller than that in the ΛCDM model. As it is shown in Figs. 1 and 2, the increases of α and Σm ν both can result in the suppression of the matter power spectrum. Therefore, when α is larger, in order to avoid suppressing the matter power spectrum too much, the value of Σm ν should be smaller. Consequently, the case with α > 0, i.e., φCDM model, has smaller Σm ν ; correspondingly, the case with α = 0, i.e., ΛCDM scenario, has larger Σm ν . Additionally, in Ref. [10], they obtained Σm ν < 0.200 eV at 95% CL in the ΛCDM model, that is consistent with our result. With the results presented in Table II of [10], one can see that our constraints on Ω b h 2 , Ω c h 2 , τ and n s are also consistent with theirs at 95% CL for the ΛCDM model.

IV. CONCLUSION
We have concentrated on a quintessence model (or called as φCDM model) of dark energy with massive neutrinos. In the φCDM model under consideration, the scalar field φ is taken as a candidate of dark energy to drive the latetime acceleration of the universe with an inverse power-law potential V (φ) ∝ φ −α (α > 0). The larger value of α corresponds to the stronger time dependence of the scalar field energy density. When α = 0, it is reduced to the corresponding ΛCDM scenario. The linear perturbation theory is employed in the framework of this model. Through qualitative analyses, we find that the increases of the sum of neutrino masses Σm ν and the parameter α both can gradually suppress the matter power spectrum P (k). It implies that when the value of α is bigger, in order to avoid suppressing the matter power spectrum too much, Σm ν should be smaller. It is in accordance with the results from the observational constraints. The variations of these two parameters also can have signatures in the CMB temperature anisotropy spectrum C T T l . In order to make a comparison, the impacts of Σm ν on P (k) and C T T l in the context of ΛCDM model have also been presented.
A combination of the CMB data from Planck 2013 and WMAP9, the galaxy clustering data from WiggleZ and BOSS DR11, and the JLA compilation of the SNe Ia observations is used to constrain the parameters. The results indicate that constraints on the cosmological parameters from this joint sample are quite restrictive. It turns out that Σm ν < 0.262 eV (95% CL) in the framework of φCDM model and Σm ν < 0.293 eV (95% CL) in the ΛCDM model. The allowed neutrino mass scale in the φCDM model is a little shrunk comparing to that in the ΛCDM model. In Ref. [10], it is concluded that the allowed neutrino mass scales in the viable f (R) models are bigger than that in the ΛCDM model. Given this, we can infer that the allowed scale of Σm ν in our φCDM model must be smaller than those in the viable f (R) models. In addition, we get α < 4.995 at 95% CL for the φCDM model, meanwhile, the ΛCDM scenario corresponding to α = 0 is not ruled out. Consequently, the observational data that we have employed here still cannot distinguish whether dark energy is a time-independent cosmological constant or a time-varying dynamical component.     In the 1D plots, the solid lines denote the marginalized likelihoods and the dotted lines correspond to the mean likelihoods. In the 2D plots, the contours refer to the marginalized likelihoods while the colors refer to the mean likelihoods. The contours correspond to 68%, 95% and 99% confidence levels.