Exponential $f(R)$ cosmology with massive neutrinos as a dynamical dark energy framework

Technical Summary: Exponential $f(R)$ Cosmology with Massive Neutrinos

Problem Statement
The standard $\Lambda$CDM cosmological model, while successful, faces persistent theoretical and observational challenges. Theoretical issues include the fine-tuning and coincidence problems associated with the cosmological constant $\Lambda$, and the lack of a confirmed dark matter particle. Observationally, significant tensions exist, most notably the "Hubble tension," a $\sim 5\sigma$ discrepancy between local distance-ladder measurements ($H_0 \approx 73.04$ km s$^{-1}$ Mpc$^{-1}$) and early-Universe inferences from the Cosmic Microwave Background (CMB) within the $\Lambda$CDM framework ($H_0 \approx 67.36$ km s$^{-1}$ Mpc$^{-1}$). Additionally, recent Baryon Acoustic Oscillation (BAO) data from the DESI collaboration suggest deviations from $\Lambda$CDM, hinting at dynamical dark energy.

A secondary tension involves the total neutrino mass ($\sum m_\nu$). While neutrino oscillation experiments establish a lower bound of $\sum m_\nu \gtrsim 0.06$ eV, cosmological data within the $\Lambda$CDM model often yield upper bounds below this threshold (e.g., $\sum m_\nu < 0.064$ eV from DESI). Phenomenological models like $w_0w_a$CDM can relax these bounds but lack a solid theoretical foundation. This work investigates whether a theoretically motivated modification to gravity, specifically exponential $f(R)$ gravity, combined with massive neutrinos, can simultaneously address the Hubble tension and the neutrino mass constraints while maintaining consistency with the standard model at high redshifts.

Methodology
The authors employ a generalized exponential $f(R)$ gravity model defined by the Lagrangian:
$$f(R) = R - 2\Lambda \left\{ 1 - \exp\left[ -\beta \left( \frac{R}{2\Lambda} \right)^\alpha \right] \right\}$$
where $\alpha$ and $\beta$ are dimensionless parameters controlling the transition from General Relativity (GR) at high curvature to modified gravity at low curvature. The model includes a matter sector comprising cold dark matter, baryons, radiation, and massive neutrinos. The neutrino energy density is computed via a phase-space integral over the Fermi–Dirac distribution, accounting for the transition from relativistic to non-relativistic regimes.

The cosmological dynamics are derived from the modified Friedmann equations in a spatially flat FLRW metric. To facilitate numerical analysis, the authors map the $f(R)$ parameters to equivalent $\Lambda$CDM parameters ($H_0^*, \Omega_i^*$) using a dimensionless Hubble parameter $E(z) = H(z)/H_0^*$. The system of coupled differential equations for $E(z)$ and the dimensionless Ricci scalar $\mathcal{R}(z)$ is integrated numerically from a high-redshift initial condition (where the model recovers $\Lambda$CDM) down to $z=0$.

A Bayesian Markov Chain Monte Carlo (MCMC) analysis is performed using the emcee package to constrain the parameter space $(H_0^*, \alpha, \beta, \omega_b^*, \omega_{bc}^*, \omega_\nu^*)$. The analysis utilizes four observational datasets:

1. Cosmic Chronometers (CC): 31 data points measuring $H(z)$.
2. DESI DR2 BAO: Measurements of $D_M(z)/r_d$ and $D_H(z)/r_d$ across six redshift bins.
3. CMB: Compressed Planck 2018 likelihoods for the acoustic scale $\theta_*$ and physical densities $\omega_b, \omega_{bc}$.
4. Pantheon+ (Pan+): 1550 Type Ia Supernovae distance moduli.

The results are compared against the standard $\Lambda$CDM model and the phenomenological $w_0w_a$CDM model using the Akaike Information Criterion (AIC) and Bayesian Information Criterion (BIC) to assess statistical preference.

Key Contributions and Results
The study provides the first joint constraint on the generalized exponential $f(R)$ model parameters ($\alpha, \beta$) and the total neutrino mass ($\sum m_\nu$) using the latest DESI DR2 data.

• Model Viability: The exponential $f(R)$ model remains consistent with current observations. The stability conditions $f_R > 0$ and $f_{RR} > 0$ are satisfied across the entire posterior parameter space and redshift range, confirming the absence of ghost and tachyonic instabilities.
• Hubble Tension: The model partially alleviates the Hubble tension. When including Pantheon+ supernovae, the inferred $H_0$ for the $f(R)$ model is $70.89 \pm 0.19$ km s$^{-1}$ Mpc$^{-1}$. This value reduces the tension with local measurements compared to $\Lambda$CDM ($71.43 \pm 0.13$) and $w_0w_a$CDM ($71.13 \pm 1.20$), though a discrepancy of $2.3\sigma$ remains between the $f(R)$ and $\Lambda$CDM values in the full dataset combination.
• Neutrino Mass Constraints: The inclusion of massive neutrinos in the $f(R)$ framework leads to a mild relaxation of the upper bound on $\sum m_\nu$ compared to $\Lambda$CDM. With Pantheon+ data, the $1\sigma$ upper limit is $\sum m_\nu < 0.0290$ eV for $f(R)$, compared to $< 0.0126$ eV for $\Lambda$CDM and $< 0.0613$ eV for $w_0w_a$CDM. While this is a relaxation, the bound remains below the oscillation lower limit ($\sim 0.06$ eV), indicating the tension is not fully resolved.
• Statistical Preference:
  • Without Pantheon+ data, $\Lambda$CDM is statistically preferred.
  • With the inclusion of Pantheon+ data, both $w_0w_a$CDM ($\Delta$AIC = 43.75) and exponential $f(R)$ ($\Delta$AIC = 20.46) show "very strong" and "strong" evidence, respectively, over $\Lambda$CDM according to AIC and BIC metrics.
• Parameter Degeneracies: The analysis reveals strong degeneracies between the modification parameters ($\alpha, \beta$) and cosmological parameters ($H_0, \Omega_m$), indicating that modified gravity can partially compensate for changes in matter density and expansion history. The posterior distribution for $\beta$ is non-Gaussian, reflecting the nonlinear role of the exponential correction.

Significance and Claims
The paper claims that the interplay between modified gravity and neutrino physics offers a viable, theoretically robust framework for investigating cosmological tensions. Unlike phenomenological dark energy parameterizations, the exponential $f(R)$ model provides a solid theoretical basis rooted in gravitational Lagrangian modifications.

The authors emphasize that while the model does not fully resolve the Hubble tension or the neutrino mass discrepancy, it demonstrates that a combined approach—modifying the gravitational sector while including massive neutrinos—can improve the statistical fit to observational data and relax constraints on neutrino mass more effectively than $\Lambda$CDM alone. The work highlights the importance of low-redshift datasets (like SNe Ia) in discriminating between alternative models and the standard cosmology. The authors conclude that future investigations must include perturbative effects (e.g., weak lensing, redshift-space distortions) to fully assess the model's viability, as the current analysis is restricted to background cosmology.