Revisiting the bulge-halo conspiracy I: Dependence on galaxy properties and halo mass

TL;DR

This study investigates how the total mass density profile of massive early-type galaxies varies with galaxy properties and halo mass, using semi-empirical models and observational data, revealing non-universality and specific profile preferences.

Contribution

It provides a systematic analysis of the dependence of galaxy mass profiles on properties and halo mass, challenging some assumptions about dark matter profiles and initial mass functions.

Findings

Density slope steepens with galaxy compactness and mass.

Models with Salpeter IMF and uncontracted halos fit data well.

Massive galaxies favor cored dark matter profiles and flatter stellar distributions.

Abstract

We carry out a systematic investigation of the total mass density profile of massive (Mstar>2e11 Msun) early-type galaxies and its dependence on galactic properties and host halo mass with the aid of a variety of lensing/dynamical data and large mock galaxy catalogs. The latter are produced via semi-empirical models that, by design, are based on just a few basic input assumptions. Galaxies, with measured stellar masses, effective radii and S\'{e}rsic indices, are assigned, via abundance matching relations, host dark matter halos characterized by a typical LCDM profile. Our main results are as follows: (i) In line with observational evidence, our semi-empirical models naturally predict that the total, mass-weighted density slope at the effective radius gamma' is not universal, steepening for more compact and/or massive galaxies, but flattening with increasing host halo mass. (ii) Models characterized by a Salpeter or variable initial mass function and uncontracted dark matter profiles are in good agreement with the data, while a Chabrier initial mass function and/or adiabatic contractions/expansions of the dark matter halos are highly disfavored. (iii) Currently available data on the mass density profiles of very massive galaxies (Mstar>1e12 Msun), with Mhalo>3e14 Msun, favor instead models with a stellar profile flatter than a S\'{e}rsic one in the very inner regions (r<3-5 kpc), and a cored NFW or Einasto dark matter profile with median halo concentration a factor of ~2 or <1.3, respectively, higher than those typically predicted by N-body numerical simulations.

Tables

Table 1. Table 1 Stellar mass-halo mass relations (Equation  9 )

IMF	$M_{\mathrm{star}}^0$	$M_{200c}^0$	$\alpha$	$\beta$
Chabrier-SerExp	10.68	11.80	2.13	1.68
Chabrier-Sersic	10.71	11.82	2.13	1.66
Chabrier-deVac	10.63	11.80	2.17	1.77
Salpeter-SerExp	10.98	11.84	2.14	1.72
Salpeter-Sersic	11.03	11.87	2.09	1.66
Salpeter-deVac	10.93	11.84	2.15	1.77
varIMF-SerExp	10.93	11.83	2.03	1.60
varIMF-Sersic	10.96	11.85	2.05	1.60
varIMF-deVac	10.83	11.81	2.04	1.64

Equations

$$\rho (r)=\frac{{\rho }_0}{\left(r/r_s\right){\left(1+r/r_s\right)}^2},$$

$$\rho (r)=\frac{b{\rho }_0}{\left(1+br/r_s\right){\left(1+r/r_s\right)}^2},$$

$$\rho (r)={\rho }_{-2}\exp (2\mu )\exp \left[-2\mu {\left(\frac{r}{r_{-2}}\right)}^{1/\mu }\right].$$

$$M_i(r_i)r_i=M_f(r_f)r_f,$$

$$M_f(r_f)=M_{\mathrm{star}}(r_f)+(1-f_{\mathrm{gal}})M_i(r_i)$$

$$\begin{array}{rl}{\gamma }^{\prime }(r)& \equiv -\frac{1}{M(<r)}{\int }_0^r\frac{d\log \rho }{d\log x}4\pi x^2\rho (x)dx=\\ & =-\frac{1}{M(<r)}\left(\rho (r)4\pi r^3-3{\int }_0^r\rho (x)4\pi x^{2}dx\right)=\\ & =3-{\frac{\text{d}\log \mathrm{M}}{\text{d}\log \text{x}}|}_{\text{x}=\text{r}}.\end{array}$$

$$\Phi (>M_{\mathrm{star}},z)={\Phi }_c(>M_{\mathrm{halo}},z)+{\Phi }_s(>M_{\mathrm{halo}},z)$$

$$\begin{array}{r}\log {\left(M/L\right)}_{\mathrm{dyn}}-\log {\left(M/L\right)}_{\mathrm{Salp}}=-0.0576+\\ +0.364\times \log \left(\frac{{\sigma }_e}{200\mathrm{km}{\mathrm{s}}^{-1}}\right)\end{array}$$

$$M_{\mathrm{star}}=M_{\mathrm{star}}^0{\left(\frac{M_{200c}}{M_{200c}^0}\right)}^{\alpha }{\left[1+{\left(\frac{M_{200c}}{M_{200c}^0}\right)}^{\beta }\right]}^{-1},$$

$$⟨{\gamma }^{\prime }⟩={\gamma }_0+\alpha (z-0.3)+\beta (\log M_{\mathrm{star}}-11.5)+\xi \log (R_{\mathrm{e}}/5),$$

$${\gamma }^{\prime }\propto {M_{\mathrm{halo}}}^p$$

$$c_{200c}=K_{200c}{\left[\frac{H(z_c)}{H(z)}\right]}^{2/3}$$

$$r_H=R_{2D}(1.356-0.0293m+0.0023m^2),$$

$$f(x)=f_0{x}^{-p}\exp (-x^m),$$

$$p=1.0-0.6097m+0.05463m^2$$

$$f_0=\frac{M_{\mathrm{star}}}{4\pi a^3}\frac{m}{\Gamma [(3-p)/m]}.$$

$$M(<r)=0.5M_{\mathrm{star}}\frac{\gamma [(3-p)/m,x^m]}{\gamma [(3-p)/m,(r_H/a{)}^m)]}$$

Full text

1. Introduction

Early-type galaxies constitute a family of objects of remarkable regularity, captured by a series of tight scaling relations such as the fundamental plane (e.g., Dressler et al., 1987; Djorgovski & Davis, 1987), the color-magnitude relation (e.g. Bower et al., 1992; Mei et al., 2012) and the relation between the mass of the central black hole and the global properties of the host galaxy (e.g., Gebhardt & et al., 2000; Ferrarese & Merritt, 2000; Shankar, 2009; Shankar et al., 2016). This regularity is reflected in their mass structure: the total mass density profile of massive111This is not necessarily the case if one considers a broader mass range (e.g., Shu et al., 2015). galaxies is well approximated, around the scale of the half-light radius, by a power-law $\rho(r)\propto r^{-\gamma}$ with slope close to isothermal ($\gamma\approx 2$) and little scatter across the population (e.g., Koopmans et al., 2006; Gavazzi et al., 2007; Koopmans et al., 2009; Barnabè et al., 2011). This compelling observation takes the name of bulge-halo conspiracy: stars and dark matter have separately density profiles significantly different from an isothermal profile, yet they conspire to produce total density profile close to isothermal (e.g., Treu & Koopmans, 2004).

Recent studies based on strong lensing and stellar dynamics have shown how the density profile of massive early-type galaxies is not exactly universal: the density slope $\gamma^{\prime}$ correlates with projected stellar mass density, being steeper for more compact objects, and anti-correlates with redshift, being shallower for higher redshift object with respect to local systems of the same mass and size (e.g., Auger et al., 2010a; Ruff et al., 2011; Bolton et al., 2012; Sonnenfeld et al., 2013; Dye et al., 2014; Tortora et al., 2014a, b). The slope of the density profile appears to also correlate with halo mass: cluster Brightest Cluster Galaxies (BCGs) have on average shallower density profiles than field galaxies of similar stellar mass (e.g., Newman et al., 2013b, a). The residual scatter left once these scaling relations are taken into account is as small as $6\%$ on the mass density slope (Sonnenfeld et al., 2013), highlighting once more the high degree of self-similarity across the population of early-type galaxies.

The origin of such a regularity is not fully understood. If massive, central early-type galaxies are predominantly grown by mergers (Naab et al., 2009; Shankar et al., 2010a, b; Bai et al., 2014; Shankar et al., 2015; Buchan & Shankar, 2016; Bellstedt et al., 2016), stochastic processes may tend to alter pre-existing scaling relations, or at least increase the scatter in the population (e.g., Nipoti et al., 2009a, b; Shankar et al., 2013), unless additional mechanisms are present such as gas dissipation (e.g., Robertson et al., 2006; Oser et al., 2010; Remus et al., 2013; Sonnenfeld et al., 2014; Remus et al., 2016). More recent numerical and semi-analytic studies suggest that even collisionless mergers could create tight correlations close to those observed (e.g., Shankar et al., 2014a; Taranu et al., 2013, 2015). Remus et al. (2013) specifically investigated the evolution of the total mass density slope in a set of galaxies in a cosmological simulation, and were able to naturally produce slopes close to isothermal consistent with observations. Despite recent progress (e.g., Dubois et al., 2013; Remus et al., 2013, and references therein), however, the still unclear sub-grid physics and large computational power required to run high-resolution cosmological simulations do not allow for a proper quantitative comparison between detailed models and observations.

A powerful, complementary approach to cosmological simulations and semi-analytic models is to use semi-empirical models (e.g., Dutton & Treu, 2014). The latter, by design, are an extremely effective way of making use of the known properties of galaxies together with minimal theoretical inputs to make testable predictions on a set of observables, and set unique, independent constraints on galaxy evolution processes.

In this work we use the observed correlations between photometric properties of massive galaxies (i.e., effective radius, stellar mass), together with basic assumptions on their host dark matter halos, as inputs to create large mock galaxy catalogs. We then test the validity of the (few) input assumptions and parameters, by comparison with data at different environments spanning from the field to the cluster halo mass scale. In the following papers of this series we will specifically include velocity dispersion in our mocks as an additional probe for models, and study the evolution of the mass density profile with redshift.

The structure of this work is the following. In Section 2 we describe the method used to construct mock samples of galaxies and relative measurements of the density slope and dark matter fraction. In Section 3 we present the predictions for the dependence of the total mass density profile on galactic properties, and environment in Section 4. In Section 5 we discuss our results and conclude in Section 6.

In the following we will adopt the reference cosmology with parameters $\Omega_{\rm m}=0.30$, $\Omega_{\rm b}=0.045$, $h=0.70$, $\Omega_{\Lambda}=0.70$, $n=1$, and $\sigma_{8}=0.8$, to match those usually adopted in the observational studies on the stellar mass function and lensing considered in this work. We will also by default adopt a Salpeter (Salpeter, 1955) as our reference Initial Mass Function (IMF), though we will carefully discuss the implications on our model predictions when switching to a Chabrier (Chabrier, 2003) or even variable (Cappellari et al., 2015) IMF.

2. Method

2.1. The Semi-Empirical Model

The aim of this work is to compare semi-empirical models based on $\Lambda$CDM halo relations coupled to observed galaxy scaling relations, to mass modeling analysis principally from strong lensing. Our semi-empirical approach makes use of mostly observational quantities with only a few, basic theoretical inputs. In this respect, it is extremely powerful as it does not rely on any physical assumption for evolving galaxies, e.g., via mergers and/or in-situ star formation.

Our procedure to build mock galaxies relies on the following steps:

• •

We take large catalogues of central, early-type galaxies with measured stellar masses, projected effective radii, and Sérsic indices (Sérsic, 1963).

• •

To each galaxy we assign a host dark matter halo via abundance matching relations.

• •

We associate to each central galaxy a 3D de Vaucouleurs (de Vaucouleurs, 1948) or Sérsic stellar profile (Sérsic, 1963), according to which data set we compare to. Each profile depends on its specific projected effective radius and Sérsic index (equal to four for a de Vaucouleurs profile), as detailed in Appendix B (e.g., Lima Neto et al., 1999).

• •

A Navarro et al. (1997, NFW hereafter) or Einasto (Einasto, 1965) profile characterizes instead each host dark matter halo.

• •

For each galaxy we then build the scale-dependent, total mass density profile.

Although we will adopt the above as our reference model, we will also extensively discuss the impact on our predictions induced by relaxing one or more of the previous assumptions. We will consider changes in the both the dark matter and stellar profiles, in the IMF, in the stellar mass-halo mass mapping, and also in dark matter concentration.

In detail, for our mocks we extract galaxies from the Meert et al. (2015) sample derived from the Sloan Digital Sky Survey (SDSS) DR7 spectroscopic sample (Abazajian et al., 2009) in the redshift range $0.05<z<0.2$, and with a probability $P(E)>0.85$ of being elliptical galaxies based on the Bayesian automated morphological classifier by Huertas-Company et al. (2011). The vast majority of the early-type galaxies in this sample are labeled as “central” galaxies according to the Yang et al. (2007) host halo catalog. Stellar masses are obtained from integrated light profiles multiplied by the color-dependent mass-to-light ratios of Bell et al. (2003) renormalized to a Chabrier IMF. In the following, depending on the data set to compare with, we will consider three different light profiles, a de Vaucouleurs+Exponential (the “cmodel”, see Bernardi et al. 2010), a Sérsic+Exponential, and pure Sérsic profiles (see, e.g., Bernardi et al., 2013). Making use of directly measured physical quantities to build the mock has the advantage that no assumptions have to be made about, e.g., intrinsic correlations and scatters (e.g., Gaussian) among galaxy sizes, stellar masses and Sérsic indices. It is important to stress that SDSS is characterized by an average seeing of $\sim 1.5$ arcsec which would correspond to projected scales of $R\sim 3$ kpc at the average redshift of $z\sim 0.1$. This observational limit is small enough not to bias the measurement of the stellar profile on scales $R\gtrsim 0.2\,R_{\rm e}$ for galaxies $M_{\rm star}\gtrsim 4\times 10^{11}\,M_{\odot}$, though it mostly prevents the possible detection of a stellar “core” or flattening below this scale. In the next Sections we will discuss if extrapolations of a power-law Sérsic profile to very small scales is supported by strong lensing and kinematic data.

On the dark matter side, our reference will be the NFW profile with mass density

[TABLE]

where $r_{s}$ is the scale radius, and $\rho_{0}$ is a reference density. Where relevant, we will also consider other analytic forms for the dark matter profile, namely a “cored” NFW profile (Newman et al., 2013b)

[TABLE]

with a reference value of $r_{\rm core}=r_{s}/b=14$ kpc as the one adopted by Newman et al. (2013b), and an Einasto profile, expressed as (see, e.g., Mamon & Łokas, 2005a, their Appendix A)

[TABLE]

In Equation 3 $\rho_{-2}$ is the local mass density at $r_{-2}$, the radius at which the logarithmic slope of the density profile is equal to -2, which is equal to $r_{s}$ for a NFW profile. We will further discuss below the values of $\mu$ adopted in this work.

Following several previous attempts in the literature (e.g., Dutton & Treu, 2014; Oguri et al., 2014), we also consider NFW profiles modified by the inclusion of adiabatic contraction and/or expansion. As in Dutton et al. (2007), assuming spherical collapse without shell crossing, we express the actual relation between the final and initial radius of the mass distribution as $r_{f}=\Gamma^{\nu}r_{i}$. The parameter $\Gamma=r_{f}/r_{i}$ is the contraction factor (Blumenthal et al., 1986), which can be numerically calculated via the conservation equation

[TABLE]

with $M_{i}$ the initial (baryonic plus dark matter) mass distribution, and

[TABLE]

the final one, with $f_{\rm gal}$ the ratio between the galaxy stellar (possibly plus gas) mass and the host halo mass (see, e.g., Barausse, 2012). The parameter $\nu$ is set to zero for no profile modification, while it assumes positive or negative values for contraction or expansion, respectively222Dutton & Treu (2014, and references therein) discuss that $\nu$ may also be correlated with the shift $\Delta\log M_{\rm star}$ from a Chabrier IMF. This effect is rather small and we ignore it in the present study, though we discuss the effects of changing IMF, host halo masses, etc… . It has already been shown that extreme contraction ($\nu\gtrsim 1$) and/or expansion $\nu\lesssim-0.5$ tend to be disfavored by present data (Dutton & Treu, 2014). We will therefore limit our analysis to non-maximal effects, such as the one by Gnedin et al. (2004) for contraction ($\nu=0.8$), and a slightly milder expansion with $\nu=-0.3$.

Once a total mass density profile for a galaxy has been constructed, one direct quantity usually adopted in the literature to compare with lensing data has been $\gamma$, i.e., the local logarithmic slope of the spherically averaged total density profile $\rho(r)\propto r^{-\gamma(r)}$. In more recent years, the quantity $\gamma^{\prime}$ has been more often adopted. The latter is defined as the mass-weighted slope of the total density profile within a radius $r$, and it is computed as (Dutton & Treu, 2014)

[TABLE]

The mass-weighted slope $\gamma^{\prime}(r)$, specifically computed within the effective radius $R_{\rm e}$, has been shown to well approximate the slope measured in joint strong lensing and stellar kinematics studies (Sonnenfeld et al., 2013). For power-law density profiles, $\gamma^{\prime}=\gamma$. In the following, unless otherwise noted, we will always refer to the mass-weighted slope $\gamma^{\prime}$ computed at $R_{\rm e}$. Sonnenfeld et al. (2013, cfr. their Figure 6) and Dutton & Treu (2014) have shown that the mass-weighted $\gamma^{\prime}$ computed as in Equation 6 is a good proxy of $\gamma^{\prime}$ derived from strong lensing and dynamical measurements (see Sonnenfeld et al., 2013, for a full discussion ).

2.2. The stellar mass-halo mass and effective radius-stellar mass relations

In this work, we focus on central massive galaxies only, and thus accordingly adopt the mean central stellar mass-halo mass relation to map galaxies into halos. The complexities behind establishing a secure mapping for satellites have been under intense study for several years (e.g., Neistein et al., 2011; Rodríguez-Puebla et al., 2012; Watson & Conroy, 2013), and are not vital to the present study focused on massive galaxies with $\log M_{\rm star}/M_{\odot}\gtrsim 11.5$, for which the fraction of satellites is very low, at the percent level (e.g., Guo et al., 2014).

The empirical stellar mass-halo mass correlations adopted as inputs of our empirical models are derived from abundance matching between the stellar mass and (sub)halo mass functions

[TABLE]

with $M_{\rm halo}=M_{\rm 200c}$ the halo masses defined within $r_{200c}$, such that the average density within $r_{200c}$ is 200 times the critical density of the Universe at redshift $z$. The $\Phi_{c}(>M_{\rm halo},z)$ term refers to the host halo mass function of Tinker et al. (2008), with the subhalo term $\Phi_{s}(>M_{\rm halo},z)$ with mass at infall $M_{\rm halo}$ taken from Behroozi et al. (2013).

For the concentration $c_{200c}=r_{200c}/r_{s}$ we use the median $c_{200c}-M_{\rm halo}$ relation as derived by Diemer & Kravtsov (2015) for our cosmology, and include a log-normal scatter of $0.16$ dex. As we discuss in Appendix A, the Diemer & Kravtsov (2015) fit is, at least at $z=0$, in very good agreement with the Bullock et al. (2001) model, as revised by Macciò et al. (2008) and Dutton & Macciò (2014). We find that the exact mass-dependent normalization and scatter of the concentration-mass relation (e.g., Dutton & Treu, 2014), have a relatively small impact to most of our conclusions, except possibly at cluster scales, as discussed below. Note that we assume the scatter in halo concentration to be uncorrelated with stellar mass or galaxy size (see discussions in, e.g., Papastergis & Shankar, 2016; Katz et al., 2016). We have anyhow checked that at fixed galaxy stellar mass, size, and host halo mass a systematic variation of $0.16$ dex in halo concentration induces a relatively small variation of $\sim 4\%$ in the implied $\gamma^{\prime}$ computed at the effective radius.

The stellar mass function on the left side of Equation 7 depends on the assumed light profile. We adopt the results by Bernardi et al. (2013), who provide stellar mass functions for all SDSS galaxies characterized by cmodel, Sérsic+Exponential, and pure Sérsic magnitudes (all based on a Chabrier IMF). As recently discussed by a number of groups, adopting Sérsic profiles naturally results in more integrated light than de Vaucouleurs profiles, and thus, at fixed mass-to-light ratio, larger abundances of massive galaxies (see, e.g., Bernardi et al., 2013; D’Souza et al., 2015; Bernardi et al., 2016b, a; Thanjavur et al., 2016; Bernardi et al., 2017b). This in turn yields steeper stellar-to-halo mappings, i.e., more massive mean stellar masses at fixed halo mass, in better agreement with direct measurements of massive brightest cluster galaxies (see, e.g., Kravtsov et al., 2014; Shankar et al., 2014b). For each light profile we explore the impact of varying the input IMF, from a Chabrier to a Salpeter and variable IMF. For the former, we simply add 0.25 dex to each stellar mass (Bernardi et al., 2010). For the latter, we take the velocity dispersion-dependent dynamical/strong lensing mass-to-light ratios from Cappellari (2016, their Figure 19), expressed as a difference with respect to the stellar (Salpeter) mass-to-light ratio

[TABLE]

with an intrinsic scatter of 0.11 dex.

The resulting mapping between stellar mass and host halo mass333Note that Equation 9 is the mean relation for the entire galaxy population. Dutton et al. (2010), however, showed that the stellar mass-halo mass relation of early-type galaxies is practically indistinguishable from the one characterizing the total population for $M_{\rm star}\gtrsim 2\times 10^{11}\,M_{\odot}$. is well reproduced by the following two-power law relation

[TABLE]

with $M_{\rm 200c}=M_{\rm halo}$, valid in the approximate range $10^{10}\lesssim M_{\rm star}/M_{\odot}\lesssim 10^{12}$ (Chabrier IMF).

We always assume an intrinsic Gaussian scatter of 0.15 dex in $\log M_{\rm star}$ at fixed halo mass, which is a valid approximation especially at high stellar masses and at low redshifts, in agreement with a number of diverse observational and theoretical studies (e.g., Guo et al., 2014; Shankar et al., 2014b; Gu et al., 2016; Tinker et al., 2016, and references therein). To include scatter, we follow the methodology outlined in Shankar et al. (2014b). We first fit Equation 9 to the results of abundance matching without scatter, and then vary the parameter $\beta$ until the implied mock galaxy-halo samples (inclusive of scatter) reproduce well the input stellar mass function within the observational uncertainties. The final values of the parameters of Equation 9 are given in Table 1 for a Chabrier, a Salpeter, and a variable IMF, and for each IMF we consider three different light profiles, a de Vaucouleurs, a Sérsic+Exponential, and Sérsic profile. For each combination of light profile and IMF, the implied mean stellar mass-halo mass relation is calculated by inserting in Equation 7 the relevant light profile-dependant stellar mass function (Bernardi et al., 2013) corrected for an appropriate mass-to-light ratio.

In turn, when assigning a host halo mass to a galaxy in our SDSS sample with measured stellar mass, we make use of the mean halo mass at fixed stellar mass relation, i.e., the “inverse” of Equation 9. The latter mapping is computed by first creating a large halo mock catalogue extracted from the halo mass function, assigning galaxies directly applying Equation 9, with a scatter of 0.15 dex and parameters given in Table 1, and then self-consistently computing the implied mean and scatter in halo mass as a function of stellar mass. As a consistency check, we make sure that our inverse relations and related scatters properly reproduce the halo mass function from random galaxy catalogues extracted from the input stellar mass functions.

In Figure 1, we show the mean stellar mass-halo mass (top), effective radius-stellar mass (middle), and velocity dispersion-stellar mass relation (bottom) relations assuming a Salpeter (left) or variable IMF (right). As labeled, each panel includes three models corresponding to three different light profiles, de Vaucouleurs (cmodel, black, solid lines with their 1$\sigma$ dispersions marked by gray areas), Sérsic-Exponential (blue, dot-dashed lines), and Sérsic (red, long-dashed lines).

It is clear from the top panels that all our models, and in particular the ones based on Sérsic profiles, predict mean stellar mass-halo mass relations with a high-mass end slope of $\delta=\alpha-\beta\sim 0.45$ (see Table 1), steeper than the Moster et al. (2013) (purple, long-dashed line) or Behroozi et al. (2013) relations, characterized by slopes $\delta\sim 0.30-0.35$ which would imply up to a factor of $\sim 3$ lower stellar masses at fixed halo mass (see, e.g., Kravtsov et al., 2014; Shankar et al., 2014b; Buchan & Shankar, 2016). All our models are in broad agreement with the Sonnenfeld et al. (2013) data (blue squares), with the Sérsic-based stellar mass-halo mass relations particularly well aligned with the Newman et al. (2013a) and Newman et al. (2015) data (green triangles and red circles, respectively), which are also characterized by Sérsic-type light profiles. All in all, the variable IMF models are characterized by scaling relations that are very similar to the ones extracted from pure Salpeter models. At best, they might predict stellar masses slightly lower by $\sim 0.1$ dex at low halo masses. Such a high degree of similarity is expected in the range of high stellar masses considered in this work.

Our models result being significantly steeper than the Moster et al. (2013) relation mainly because ours are based on the Bernardi et al. (2013) stellar mass functions, which predict up to two orders of magnitude more massive galaxies than the Li & White (2009) stellar mass function, adopted as the reference one in the Moster et al. (2013) analysis. The Li & White (2009) stellar mass function adopted Petrosian magnitudes (Petrosian, 1976), which have long been known to underestimate the light in the most massive galaxies (e.g., Bernardi et al., 2010). Moreover, the Petrosian magnitudes used by Li & White (2009) were based on the SDSS pipeline, which suffers from bad sky estimates, particularly in crowded fields where massive galaxies are typically found (see, e.g., Bernardi et al., 2010, 2014, and references therein). More recent work by D’Souza et al. (2015) and Bernardi et al. (2016a) confirms that the Li & White (2009) photometry was indeed biased faint. Finally, Li & White (2009) adopted stellar mass-to-light ratios based on templates from Blanton & Roweis (2007), which are inappropriate for massive galaxies (e.g., Bernardi et al., 2013, 2016b, 2016a, 2017b). Bernardi et al. (2017a) have also carefully demonstrated that the difference with the SDSS pipeline photometry does not depend on whether or not one includes intracluster light, but rather it is truly missing light from the body of the galaxy.

The middle panels of Figure 1 report instead the two-dimensional effective radius $R_{\rm e}$ versus stellar mass for the Salpeter (left) and variable IMF (right), and for our three reference stellar light profiles, as labeled. It is clear that Sérsic-based profiles provide significantly more extended light profiles at fixed stellar mass, and thus proportionally larger effective radii up to a factor of $\gtrsim 2-3$ higher than what inferred from pure de Vaucouleurs fits. Within uncertainties, our SDSS-cmodel data compare reasonably well with the effective radii and stellar masses of the galaxies from the Sonnenfeld et al. (2013)444The lenses in Sonnenfeld et al. (2013) come from the SL2S and SLACS surveys (Auger et al., 2010a; Gavazzi et al., 2014)., which are based on de Vaucouleurs profiles. Even closer is the match between our Sérsic-based size-mass relations with the data by Newman et al. (2013a) and Newman et al. (2015), whose half-light radii and total luminosities are extracted from pseudo-isothermal and Sérsic profiles, respectively. Note that relatively tight correlations such as those between stellar mass and halo mass, and size versus stellar mass, should naturally induce additional dependencies between galaxies and their dark matter halos. We find, for example, that as a consequence of the relation between galaxy size and stellar mass, and the assumption of a monotonic stellar mass-halo mass relation, our models produce a roughly linear correlation between effective radius and virial radius, in broad agreement with the one independently measured by Kravtsov (2013).

For completeness, the bottom panels of Figure 1 plot the line-of-sight velocity dispersion within half of the effective radius555Our SDSS velocity dispersions are defined within an aperture of $R_{e}/8$ and have been corrected to $R_{e}/2$ following the mean aperture correction given by Cappellari et al. (2006), based on well-resolved two-dimensional spectroscopic data from SAURON. as a function of stellar mass, in the same format as for the other panels, for Salpeter (left) and variable (right) IMF and for the three reference light profiles. In this case we find the velocity dispersions in the Sonnenfeld et al. (2013) data to be slightly higher by $\sim 0.05$ dex, at fixed stellar mass or, alternatively, stellar masses to be lower by $\sim 0.1$ dex at fixed velocity dispersion. Irrespective of this, it is indeed quite remarkable that, despite possible differences and biases in the selections, and in the surface brightness profile fitting and/or background subtraction algorithms, the data by Sonnenfeld et al. (2013), Newman et al. (2013a) and Newman et al. (2015) show with respect to our SDSS catalogs offsets of at most $\lesssim 0.1$ dex in effective radius, stellar mass or velocity dispersion.

It could be argued from the apparent (systematic) discrepancy in the (mean) velocity dispersion-stellar mass relations (lower, bottom panel of Figure 1), that the reference sample of lensed galaxies is biased high in velocity dispersion at fixed stellar mass with respect to SDSS galaxies. We disfavor a substantial bias in the SLACS galaxies as the galaxies do not appear proportionally more compact at fixed stellar mass (middle, left panel in Figure 1). Indeed, this (relatively small) discrepancy could be mostly ascribed to residual differences between the cmodel and de Vaucouleurs light profiles, and/or between mass-to-light ratios (see, e.g., Bernardi et al., 2013, 2017b, and references therein). Auger et al. (2010a) claimed a closer match between SLACS galaxies and SDSS, however their comparison was based on the velocity dispersion-stellar mass relation calibrated by (Hyde & Bernardi, 2009), who adopted mass-to-light ratios systematically slightly lower than the ones adopted by Bernardi et al. (2013). Treu et al. (2006) also did not find any difference between SLACS and SDSS galaxies at fixed velocity dispersion. More relevantly to what follows, we verified that all our main results are robust against the level of mild systematics seen in the left panels of Figure 1. For instance, we checked that reducing stellar masses in our SDSS mocks by $\sim 0.1$ dex does not affect any of our conclusions. Even when selecting only the SDSS galaxies above the mean cmodel velocity dispersion-stellar mass relation to better match the SLACS sample, also yields similar results. The latter behavior is not unexpected given that, on average, we verified that massive galaxies in all our models share a very weak dependence of $\gamma^{\prime}$ on velocity dispersion (see also Poci et al., 2017; Xu et al., 2016).

In the following, we will adopt the scaling relations discussed in this section as inputs for our mock galaxy-halo catalogs. For self-consistency reasons, we will specifically adopt the cmodel and single-Sérsic scaling relations when comparing with the Sonnenfeld et al. (2013) and Newman et al. (2013a, 2015) data, respectively. Unless otherwise stated, from now on all our stellar masses $M_{\rm star}$ will refer to Salpeter stellar masses, to make contact with the observations considered in this work, though the input IMF may not necessarily be a Salpeter IMF. In practice, when computing $\gamma^{\prime}$ in some models we will assume a different input IMF, Chabrier or variable IMF, but when comparing to the data we will convert our mock stellar masses into the values an observer would measure assuming a Salpeter IMF. This, we remind the reader, corresponds to a positive offset of 0.25 dex (e.g., Bernardi et al., 2010) in stellar mass in case the input IMF is a Chabrier IMF, or a velocity-dependent offset as given in Equation 5 for a variable IMF. All other photometric properties, namely effective radii and Sérsic indices, are instead clearly independent of the choice of input IMF but mainly rely on the choice of assumed light profile.

3. Dependence of $\gamma^{\prime}$ on galactic properties

In this section, we compare the predictions of our models to observational data at $z<0.3$. For simplicity we start by analyzing the average predictions of the semi-empirical models in Section 3.1, and in Section 3.2 we proceed by analyzing further dependencies on galaxy properties. In this work we will not discuss the dependence of $\gamma^{\prime}$ on velocity dispersion, which shares with $\gamma^{\prime}$ a strong correlation in the observational errors (e.g., Sonnenfeld et al., 2013). We will instead mainly focus on dependencies of $\gamma^{\prime}$ on independently-measured quantities such as effective radius and stellar mass. We will get back to velocity dispersions when modeling the velocity dispersion-stellar mass relation in separate work. In Section 4 we will discuss environmental trends as quantified by the host halo mass. We expect the latter to be an important variable as massive galaxies at a given stellar mass could live in quite diverse halos.

3.1. Average trend

The top panels of Figure 2 show the predicted 3D mass density profile for galaxies with stellar mass (Salpeter) in the range $11.3<\log M_{\rm star}/M_{\odot}<11.7$ (left panels) and $11.7<\log M_{\rm star}/M_{\odot}<12.2$ (right panels). The bottom panels of Figure 2 show the ratio of the total mass density profiles with respect to a pure power-law, a profile of the type $\rho(r)\propto r^{-2.2}$. The latter choice is motivated by the findings of Cappellari et al. (2015, see also ) who, confirming and extending previous results, infer from two-dimensional stellar kinematics a total mass density profile of $\rho(r)\propto r^{-\gamma}$ with $\langle\gamma\rangle=2.19\pm 0.03$, valid in the range $0.1R_{\rm e}$ to $4R_{\rm e}$ and $10.2\lesssim\log M_{\rm star}/M_{\odot}\lesssim 11.7$ (their results are reported with pink stripes in the left panels of Figure 2).

It is first of all evident from Figure 2 that the inner regions of massive galaxies are dominated by the stellar component (blue dot-dashed lines), with the dark matter (red long-dashed lines) gradually taking over at larger scales $r\gtrsim R_{\rm e}$, as already pointed out over a decade ago (e.g., Borriello et al., 2003; Mamon & Łokas, 2005a, b). Nevertheless, the total mean mass density profile (black solid lines) is remarkably well approximated by $\gamma\approx 2.2$ (e.g., Koopmans et al., 2009; Cappellari et al., 2015), especially in the region around the effective radius, with a slight tendency towards a somewhat warped profile (Mamon & Łokas, 2005b; Chae et al., 2014). This prediction of the semi-empirical model is in very good agreement with direct strong lensing and dynamical measurements (e.g., Cappellari et al., 2015). It is important to note that the Sérsic effective radii in massive galaxies tends to become progressively larger than the one based on de Vaucouleurs profiles, thus naturally probing regions with flatter density profiles (right panels).

Figure 2 also includes two Hernquist profiles, with their core radii linearly correlated to the 3D half-mass radii (Hernquist, 1990). The latter are in turn derived from the 2D circularized effective radii (vertical black dotted and orange dashed lines) adopting Equation B1. The green dashed and purple dashed lines show, respectively, the resulting Hernquist profiles for the same stellar masses but different core radii derived from a de Vaucouleurs and Sérsic 2D effective radius and index (equal to four for the de Vaucouleurs profile). It is clear that, irrespective of the input effective radius, the Hernquist profiles are disfavored with respect to a de Vaucouleurs or, even better, a pure Sérsic profile (purple dot-dashed lines), tending to drop too quickly in the inner few kiloparsecs with respect to the data by Cappellari et al. (2015, left bottom panel).

3.2. Dependence on stellar mass and effective radius

While probing the total mass density profile of massive galaxies has been a subject of intense study by a number of teams (e.g., Koopmans et al., 2009; Barnabè et al., 2011; Ruff et al., 2011; Bolton et al., 2012), only recently the galaxy lensing samples have become large enough to enable the study of its dependence on galactic properties (e.g., Auger et al., 2010a; Cappellari et al., 2013; Chae et al., 2014; Shu et al., 2015; Serra et al., 2016). Using strong lensing, Sonnenfeld et al. (2013) have empirically parameterized the dependence of the mass-weighted $\gamma^{\prime}$ on galactic properties as (their Table 4, second column)

[TABLE]

with $\gamma_{0}=2.08^{+0.02}_{-0.02}$, $\alpha=-0.31^{+0.09}_{-0.10}$, $\beta=0.40^{+0.16}_{-0.15}$, and $\xi=-0.76^{+0.15}_{-0.15}$, with an average dispersion around the median of $\sigma_{\gamma^{\prime}}=0.12^{+0.02}_{-0.02}$.

The data tend to favor a positive increase of $\gamma^{\prime}$ with increasing stellar mass, at fixed effective radius. This can be broadly understood in terms of the progressively less significant impact of the (steeper) stellar component when moving to lower mass galaxies, inducing a parallel flattening of the overall profile and thus of $\gamma^{\prime}$. When fixing the stellar mass and varying the effective radius, $\gamma^{\prime}$ is instead expected to steepen with decreasing size, as the profile gets progressively dominated by the steeper stellar component. Similarly, $\gamma^{\prime}$ is expected to steepen even faster with mean (stellar) surface density at fixed stellar mass, since $\Sigma=M_{\rm star}/(2\pi R_{\rm e}^{2})$.

In essence, the dependence of $\gamma^{\prime}$ on stellar mass and/or effective radius is connected to probing different portions of the total mass density profile, from the inner steeper profile of the stellar component to the outer flatter profile of the dark matter component, as well as varying the relative contributions of the stellar and dark matter components (cfr. Figure 2). Although the trends of $\gamma^{\prime}$ summarized in Equation 10 are broadly expected, their precise measurements provide invaluable constraints to galaxy evolutionary processes. To this purpose, in this section we will compare detailed predictions of our semi-empirical models to the Sonnenfeld et al. (2013) data and specifically to Equation 10.

Figure 3 reports the probability distribution function (PDF) measured by Sonnenfeld et al. (2013) evaluated at the average redshift of our mock SDSS sample, $z=0.1$. The red, yellow, and pink regions bracket, respectively, the 1$\sigma$, 2$\sigma$, and 3$\sigma$ contour levels for each pair of the four parameters $\xi$, $\beta$, $\gamma_{0}$, and $\sigma_{\gamma^{\prime}}$. In order to make the closest comparison to the data we treat the outputs of each rendition of our semi-empirical model as in the observations. We apply the same fitting procedure used in Sonnenfeld et al. (2013) to infer the distribution of $\gamma^{\prime}$ across the population of galaxies, which we model as a Gaussian distribution with a mean given by Equation 10 and with dispersion $\sigma_{\gamma^{\prime}}$. The fit produces a posterior PDF for the parameters $\gamma_{0}$, $\beta$, $\xi$ and $\sigma_{\gamma^{\prime}}$. For each model realization we plot the peak of the posterior PDF in Figure 3. Due to the large sample size of our mock datasets, statistical uncertainties on these parameters are very small compared to the uncertainty on the same parameters in the Sonnenfeld et al. (2013) observations, and are therefore omitted in the figure. The circles, triangles, and squares in Figure 3 refer, respectively, to models with standard profiles, with adiabatic expansion, and with adiabatic contraction, as labeled, computed as discussed in Section 2.1.

The first apparent feature characterizing the Salpeter IMF model is that a standard, uncontracted NFW dark matter profile (black circles) is consistent with the data at the $1-1.5\sigma$ level in all the parameters’ spaces identified by Sonnenfeld et al. (2013). Adiabatic contraction will instead tend to decrease $\beta$ and increase $\xi$, significantly worsening the match to the data (black squares). The endpoint of adiabatic contraction is in fact to attract dark matter towards more central regions (cfr. Figure 9). By increasing the relative contribution of dark matter at all scales weakens the dependence of $\gamma^{\prime}$ on stellar mass, thus effectively decreasing the value of $\beta$ in Equation 10. At the same time, adiabatic contraction will also weaken the dependence of $\gamma^{\prime}$ on effective radius, rendering $\xi$ less negative, as the profile tends to steepen at larger scales increasing $\gamma^{\prime}$. Conversely, adiabatic expansion (black triangles) tends to progressively lower the dark matter contribution in the inner regions allowing for a more dominant role of the steeper stellar component, thus increasing the zero point $\gamma_{0}$, increasing $\beta$, and decreasing $\xi$. Figure 3 shows that even milder adiabatic contractions or expansions, as the ones adopted here, tend to be ruled out at the $\gtrsim 3\sigma$ level in at least two or more pairs of parameters. Adiabatic expansion, in particular, appears to be excluded at high significance (black triangles).

The second important point to mention is that a Chabrier IMF (green symbols in Figure 3) is also disfavored by the data, irrespective of any adiabatic contraction/expansion included in the models. The mean $\gamma_{0}$, in particular, appears to be $>3\sigma$ away the values constrained by the data, irrespective of the exact dark matter profile assumed (green symbols in the middle panels). The model with a Chabrier IMF predicts in fact lower stellar masses at fixed halo mass, which in turn reduces the contribution of the steeper stellar component and inducing a flatter total mass density profile (see Figure 2), especially around the effective radius, in tension with the strong lensing data by Sonnenfeld et al. (2013, middle panels in Figure 3).

In Appendix C we also show a series of other systematic variations in the inputs of our reference semi-empirical model (Section 2.1). First off, the PDFs for a Salpeter and variable IMF are quite similar and thus they are both consistent with the Sonnenfeld et al. (2013) strong lensing data (Figure 13). This is partly expected given that the input scaling relations appear to be quite similar for both choices of IMF (Figure 1). In future work we will explore additional observational proxies such as the velocity dispersion-stellar mass relation to hopefully break these type of degeneracies. Assuming a Hernquist stellar profile (Figure 14) also has a relatively mild effect on the model predictions, as expected given that it mainly flattens out the mass density profile within the regions inside the effective radius (Figure 12), leaving the total mass density profile almost unaltered around $R_{\rm e}$. Nevertheless, the model characterized by a Hernquist profile in the stellar component tends to be only marginally consistent with the data at the $\sim 2-3\sigma$ level. Moreover, a shallower stellar profile in the inner regions of intermediate-mass early-type galaxies tends to be in tension with strong lensing and 2D dynamical measurements which point to steep total mass density profiles (e.g., Auger et al., 2010a; Cappellari et al., 2015), as we showed in Figure 2. Introducing a core in the dark matter (Figure 15) severely worsens the match with the data in all parameters, mainly inducing steeper mass densities profiles thus larger $\gamma^{\prime}$ values (middle panels). In this case adiabatic contraction may help bringing model predictions closer to the data, though discrepancies remain at the $\sim 3\sigma$ level, especially in terms of the very weak dependence on stellar mass ($\beta\sim 0$) and also significantly weaker in size (green squares). As expected, also models characterized by combinations of Hernquist plus cored dark matter profiles tend to be significantly less aligned with the data, being excluded at the $\sim 3\sigma$ level, irrespective of any adiabatic contraction or expansion (Figure 16). Switching to the Moster et al. (2013) stellar mass-halo mass relation (Figure 17) has a relatively small impact to model predictions, as the Sonnenfeld et al. (2013) data are dominated by intermediate-mass galaxies, and not brightest cluster galaxies, where the difference becomes most prominent (Figure 1). Nevertheless the adiabatic contraction or expansion tend to still be excluded at $\sim 3\sigma$ level.

4. Dependence of $\gamma^{\prime}$ on Environment

So far, we have been discussing trends of $\gamma^{\prime}$ with galaxy properties, ignoring the possible effects of the larger-scale environment, as labeled by the host dark matter halo mass. In order to add this important dimension to our analysis we need to extend our data sets adopted as comparison. As discussed by Treu et al. (2009), the SLACS survey, which is the foundation of the Sonnenfeld et al. (2013) sample, is dominated by massive galaxies not in cluster-scale halos666We point out that when comparing with Sonnenfeld et al. (2013) in Figure 3 or in the Appendices, we do not apply any halo mass cut in our mock galaxies, though we verified that retaining only galaxies with, say, $\log M_{\rm halo}/M_{\odot}<13.5$, does not alter our main results. (only $\sim 17\%$ of the original sample resides at the center of massive groups and clusters). A stacked weak lensing analysis of a subset of those galaxies also indicated a mean host halo mass of $M_{\rm halo}\sim 1.6\times 10^{13}\,M_{\odot}$ (Gavazzi et al., 2007), with an intrinsic scatter of a factor of $\sim 2$ (Newman et al., 2015). The Sonnenfeld et al. (2013) sample also tends to run out of galaxies above $M_{\rm star}\sim 9\times 10^{11}\,M_{\odot}$. We thus complemented our reference observational sample with 10 galaxies from Newman et al. (2015), specifically selected to be very massive galaxies with mean $M_{\rm star}\sim 10^{12}\,M_{\odot}$ (Salpeter IMF) at the center of massive groups with an average mass $M_{\rm halo}\sim 10^{14}\,M_{\odot}$, and seven very massive galaxies by Newman et al. (2013a), all at the center of clusters and with stellar mass $M_{\rm star}\gtrsim 10^{12}\,M_{\odot}$ within clusters of mass $4\times 10^{14}\lesssim M_{\rm halo}/M_{\odot}\lesssim 2\times 10^{15}$, with total mass density profiles obtained via a combination of weak and strong lensing, resolved stellar kinematics within the BCG and X-ray kinematics (Table 9 in Newman et al. 2015). All halo masses are always self-consistently defined as $M_{\rm 200c}$ (see Section 2).

4.1. Dependence on host halo mass

To start off with, we compare in Figure 4 the dependence of $\gamma^{\prime}$ versus (Salpeter) stellar mass for the whole Sonnenfeld et al. (2013) lensing sample (gray band), with the “group” galaxy sample by Newman et al. (2015, green triangles), and with the cluster galaxy sample by Newman et al. (2013a, red circles). The Sonnenfeld et al. (2013) data indicate a very weak dependence of $\gamma^{\prime}$ on stellar mass. In particular, at masses around $M_{\rm star}\sim 10^{12}\,M_{\odot}$, the three samples predict very different total mass density slopes, ranging from values of \gamma^{\prime}$$\sim 2.2 to values as low as \gamma^{\prime}$$\sim 1.1-1.2. To better understand the origin of such huge discrepancies in these different data sets, we compare with the predicted mean $\gamma^{\prime}$ for several different models, a standard de Vaucouleurs profile for a Chabrier, Salpeter and variable IMF (blue dot-dashed, red long-dashed, and black solid lines, respectively), and a Sérsic profile with a Salpeter and variable IMF (purple dotted, and green dashed lines, respectively), as labeled. The first point to note is that all models predict a more or less rapid drop of $\gamma^{\prime}$ with increasing stellar mass, irrespective of the input stellar profile or IMF (left panel). However, if in our mocks we select only “field” galaxies (right panel), i.e. we disregard galaxies at the center of massive groups and clusters above, say, $\log M_{\rm halo}=13.5$ (which is around the upper limit probed by Sonnenfeld et al. (2013)), we tend to find steeper profiles in high stellar mass galaxies induced by the lower contribution of the host dark matter halo characterized by flatter profiles. This simple cut in host halo mass is sufficient to reconcile de Vaucouleurs model predictions with the Sonnenfeld et al. (2013) data (gray band). It can be seen that both the de Vaucouleurs models with Salpeter and variable IMF match observations, with the Chabrier IMF model predicting generally lower $\gamma^{\prime}$, in line with what discussed in reference to Figures 3 and 13.

However, even after excluding massive host halos, the Sérsic-based models tend to fall to lower $\gamma^{\prime}$ values at large stellar masses. This key difference with respect to the predictions of de Vaucouleurs-based models could be mostly ascribed to the progressively larger effective radii in more massive galaxies with Sérsic profiles (middle panels of Figure 1). Having larger radii, the Sérsic-based models probe larger and flatter portions of the mass density profile, in which the dark matter starts dominating, as emphasized in the right panel of Figure 2. We conclude that the Newman et al. galaxies show significantly flatter total mass density profiles with respect to Sonnenfeld et al. mainly because they are central galaxies in groups and clusters, thus they are physically more extended and better described by Sérsic profiles (Bernardi et al., 2011a, b, 2014), and also characterized by significantly more massive host dark matter halos.

To further probe the dependence on structural properties, in Figure 5 we plot $\gamma^{\prime}$ against effective radius for the most massive galaxies in the Sonnenfeld et al. (2013) catalog, and the group and cluster galaxy samples of Newman et al. (2015) and Newman et al. (2013a), shown by blue squares, green triangles, and red circles, respectively. We include three semi-empirical models in this figure, characterized by a standard stellar Sérsic and dark matter NFW profiles (black, solid lines), a stellar Hernquist and dark matter NFW profiles (blue, dot-dashed lines), and a stellar Sérsic and dark matter NFW based on the stellar mass-halo mass relation by Moster et al. (2013, purple, dotted lines). For consistency with the data, in Figure 5 (as well as in Figure 6) we show results for the subsample of mock galaxies with $\log M_{\rm star}>11.9$, which is comparable to the lower limit in stellar mass probed in the Newman et al. samples, and still low enough to retain a statistically relevant number of galaxies in our mocks. For the same reason, we also restrict the Sonnenfeld et al. (2013) data to only galaxies with $\log M_{\rm star}>11.9$. Despite the variety of mass profiles and stellar-dark matter mappings the predicted $\gamma^{\prime}$ distributions show quite similar behaviors, differing only in normalization by just a modest amount of $\lesssim 10\%$. More interestingly, the predicted mean $\gamma^{\prime}$ substantially decreases with increasing size, as the (flatter) dark matter component becomes more dominant at larger scales (see discussion in Section 3.2). Our model predictions, at least the ones without adiabatic contraction, are in fact in broad agreement with the data, and the predicted trend is by itself sufficient to explain a large portion of the difference between the $\gamma^{\prime}$ values in the Sonnenfeld et al. (2013) sample and in the group galaxies of Newman et al. (2015).

However, the models appear to still be highly inconsistent with the $\gamma^{\prime}$ measurements by Newman et al. (2013a) for BCGs (red circles). Their median $\gamma^{\prime}$ lays significantly below our estimates, and even below the data by Newman et al. (2015) for galaxies having similar stellar masses and effective radii but residing in group-scale halos. The mean offset in stellar mass between the Newman et al. (2015) group sample and the Newman et al. (2013a) cluster sample is in fact only $\sim 0.17$ dex, which is too modest to explain the strong observed systematic difference in $\gamma^{\prime}$ purely in terms of stellar mass. This finding is consistent with the hypothesis that halo mass may play a non negligible role in shaping the overall mass density profiles of very massive galaxies at the center of clusters (see also Kochanek & White, 2001).

To explore this issue further, we present in Figure 6 the predicted dependence of $\gamma^{\prime}$ on host halo mass from our reference models, compared with data from Sonnenfeld et al. (2013) for field/group galaxies (blue squares), Newman et al. (2015) for group central galaxies (green triangles), and Newman et al. (2013a) for BCGs (red circles). It is first of all apparent, as also pointed out by Newman et al. (2015), that $\gamma^{\prime}$ depends on environment/halo mass. While some of the data appear more mixed in the $\gamma^{\prime}$-$R_{\rm e}$ plane (Figure 5), they tend to separate well when $\gamma^{\prime}$ is expressed as a function of host halo mass. Predictions for a series of models (with Salpeter IMF) are included in Figure 6: a standard stellar Sérsic and dark matter NFW profiles (black, solid lines), a stellar Sérsic and dark matter cored NFW profiles (red, long-dashed lines), a stellar Hernquist and dark matter NFW profiles (blue, dot-dashed lines), a stellar Hernquist and cored NFW profiles (orange, dashed lines), and a stellar Sérsic and dark matter NFW profiles adopting the stellar mass-halo mass relation by Moster et al. (2013, purple, dotted lines).

All models, especially the ones with no adiabatic contraction, predict a drop of $\gamma^{\prime}$ with increasing halo mass in line with the data. Our reference model with uncontracted NFW dark matter profiles, in particular, predicts (left panel) a dependence on halo mass of the type

[TABLE]

with $p=-0.30\pm 0.01$, in good agreement with Newman et al. (2015) who find $p=-0.33\pm 0.07$ when combining different samples of lensed galaxies in the field, in groups, and in clusters (each single data set, by itself, does not show any dependence in halo mass). The value of $p$ we find is strictly valid for $\log M_{\rm halo}/M_{\odot}\gtrsim 13$, with a possible gradual flattening below this mass scale. To make a fairer comparison with the data, we again consider in Figure 6 only galaxies from our semi-empirical models with $\log M_{\rm star}/M_{\odot}>11.9$. However, we have verified that the dependence of $\gamma^{\prime}$ on halo mass is quite independent of the exact cuts in stellar mass or effective radius. As such it could be considered as an addition to Equation 6. We note, however, that observational calibrations of the slope $p$ require uniform and precise measurements of strong and weak lensing galaxies over a large range of host halo masses. Non-negligible errors in halo mass could in fact wash out some of the trends with halo mass (e.g., Shankar et al., 2014a). The dependence of $\gamma^{\prime}$ on halo mass is a genuine effect and only partly induced by the increase of effective radius with stellar mass. This is clearly seen in the bottom panels of Figure 6 in which we plot the same models as in the upper panels but restricted in the narrow bins of stellar mass and effective radius, $12.05<\log M_{\rm star}/M_{\odot}<12.15$ and $1.2<\log R_{\rm e}/{\rm kpc}<1.3$. In this case the trend of $\gamma^{\prime}$ with halo mass gets somewhat flatter, but still highly significant $p=-0.21\pm 0.02$.

While Equation 11 has a relatively weak dependence on galaxy properties, it has a significant dependence on the choice of input NFW dark matter halo profile. The left, middle, and right panels of Figure 6 show predictions for the same models assuming dark matter halo profiles uncontracted, with Adiabatic Expansion, and Adiabatic Contraction, respectively. It is clear that irrespective of the chosen stellar and/or dark matter mass profiles, and/or stellar mass-halo mass mapping, the predicted values of $\gamma^{\prime}$ tend to decrease at broadly the same pace with increasing halo mass, with only a small offset in normalization of $\lesssim 15\%$. This behavior is mainly induced by structural non-homology, with the relative density distribution of dark matter becoming gradually more dominant in cluster-sized halos. Despite decreasing with increasing halo mass, the predicted values of $\gamma^{\prime}$ still tend to be steeper than those actually measured in haloes above $M_{\rm halo}\gtrsim 3\times 10^{14}\,M_{\odot}$ by Newman et al. (2013a). Including adiabatic expansion in the models yields a stronger decrease of $\gamma^{\prime}$ with effective radius, mirroring the flattening of the mass profile at $r\gtrsim R_{\rm e}$ (see Section 3.2), in apparent better agreement with the data. Still, we do not properly reproduce the total mass density profiles derived by Newman et al. (2013a) for dominant galaxies in rich clusters.

Figure 7 reports a closer comparison between model predictions and the total mass density profile measured by Newman et al. (2013a). In all top panels of Figures 7 and 8, the black long-dashed line and gray area mark, respectively, the Newman et al. (2013a) mean total mass density profile and 1$\sigma$ dispersion, while the light green area roughly marks their $1\sigma$ dispersion in the BCG stellar profile. The bottom panels of Figures 7 and 8 plot instead the ratios between model predictions and data for each model. In all panels the predicted stellar and dark matter contributions to the total mass density profiles are reported by blue triple dot-dashed, and red long-dashed lines, respectively. The left panel in Figure 7 shows predictions for the reference stellar Sérsic and dark matter NFW profiles. There are two main issues. The reference model is much steeper than the data at small scales, severely overpredicting them at $r\lesssim 4$ kpc by up to a factor of $\sim 3-4$. At intermediate scales $10\lesssim r\lesssim 300$ kpc the model instead significantly underpredict the data at by a systematic factor of $\sim 2$.

Switching to a Hernquist profile in the stellar component naturally produces profiles at small scales much better aligned with the data (middle panel of Figure 7). The resulting stellar profile is indeed in good agreement with the stellar profile inferred by Newman et al. (2013a, light green areas), predicting a density of $\rho_{\rm star}\sim 3\times 10^{8}\,M_{\rm star}/{\rm kpc^{3}}$ at $r\sim 2$ kpc, and $\rho_{\rm star}\sim 1.5\times 10^{7}\,M_{\rm star}/{\rm kpc^{3}}$ at $r\sim 10$ kpc. A Hernquist profile is a simple but clearly not a unique solution. Adding, for example, an appropriate “core” radius (e.g., Trujillo et al., 2004; Newman et al., 2013a) to the stellar profile could yield similar improvements. A flatter core in the inner regions of very massive galaxies, especially at the center of clusters, has been indeed reported several times in the literature (e.g., Kormendy & Bender, 2009; Postman et al., 2012; Huang et al., 2013; Oldham & Auger, 2016). On the other hand, simply assigning a lower Sérsic index of, say, $n\lesssim 3$ to the most massive galaxies, although flattening the inner density profile (cfr. Figure B), would also decrease the predicted stellar mass density profiles in more external regions around a few to ten kiloparsecs, worsening the match to the galaxy stellar profiles measured by Newman et al. (2013a). Moreover, detailed light profile fits to massive galaxies have revealed their Sérsic indices to be as large as $n\gtrsim 5$ (e.g., Huertas-Company et al., 2013; Bernardi et al., 2014).

Despite the possible improvements at small scales, our reference models still significantly underpredict the mass density observationally inferred by Newman et al. (2013a) at intermediate scales beyond $\gtrsim 10$ kpc, where the mass profile is expected to become fully dominated by the dark matter component. We find that additionally allowing for a boost of a factor of $2$ in the median halo concentration helps resolving the mismatch with the data, as shown in the right panel of Figure 7. Increasing the concentration in fact decreases the core radius $r_{s}=r_{200c}/c_{200c}$ at fixed halo mass and virial radius $r_{200c}$, thus naturally increasing the mass density in the inner regions, with the end result of flattening the total mass density profile, especially around the effective radius (vertical, dotted lines in Figure 7). However, we note that an increased concentration also tends to predict a larger mass density at $r\lesssim 10$ kpc. To avoid the latter effect, when adopting a higher concentration we must also allow for a core in the dark matter profile, for which we adopt Equation 2. We checked that also adiabatic expansion could reduce the dark matter in the inner regions at comparable levels to those of a core.

Is there any observational justification for increasing concentrations over the predictions of N-body simulations? An increase in dark matter concentration at cluster scales is indeed possible (e.g., Comerford & Natarajan, 2007; Auger et al., 2013, and references therein), though at present good evidence for real clusters to be on average more concentrated is still not convincing (e.g., Ettori et al., 2010; Dutton & Macciò, 2014; Meneghetti et al., 2014; Merten et al., 2015; Shan et al., 2015; Umetsu et al., 2016, and references therein). Newman et al. (2013a), and more recently Amodeo et al. (2016) at $z\gtrsim 0.4$, empirically constrained the median concentration at $M_{\rm halo}\sim 10^{15}\,M_{\odot}$ from their samples to be $c\sim 5-6$, in broad agreement with our inferred values of $c\sim 2\times 3\div 4\sim 6\div 8$, which are consistent with the numerical results by Prada et al. (2012, but see ). Evidence for larger concentrations for red, massive galaxies consistent with what presented in Figure 7 also come from the kinematics of satellites, though mostly inferred from isolated galaxies (Wojtak & Mamon, 2013). Giocoli et al. (2014, see also ) discuss how selection effects may bias the measured concentration$-$mass relation at cluster scales, yielding up to 30% higher concentrations than what predicted by N-body simulations (see also Auger et al., 2013; Sereno et al., 2015; Lieu et al., 2017). Several studies have also highlighted the possibility of biases arising from halo triaxiality (e.g., Foëx et al., 2014). An alignment of the major axis and the line of sight may cause an apparent higher concentration, that may be evident in lensed groups and clusters. Although the Newman et al. (2013a) cluster sample is relatively small, the selection is well defined, and in fact Newman et al. (2013a) find that only in one case there is substantial elongation along the line of sight (A383) which could artificially boost the measured concentration.

As numerical simulations and direct observations do not tend to unanimously support dark matter concentrations in cluster-sized halos much higher than what predicted by numerical simulations, we have explored in Figure 8 the effects of switching from a NFW to an Einasto profile, which is characterized by an additional parameter $\mu$ (Equation 3), and it has been shown to be an accurate fit to massive dark matter halos by a number of groups (e.g., Navarro et al., 2004; Stadel et al., 2009; Navarro et al., 2010; Ludlow et al., 2013; Dutton & Macciò, 2014; Klypin et al., 2016). The left panel of Figure 8 shows the results of a standard Sérsic+unperturbed dark matter Einasto profile with $\mu=4$ (or $\alpha=1/\mu=0.25$), which is the average value extracted from numerical N-body simulations of very massive clusters (e.g., Navarro et al., 2004; Gao et al., 2008; Dutton & Macciò, 2014). This basic model tends to present broadly the same shortcomings as the standard NFW-based model presented in the left panel of Figure 7. However, the match is substantially improved at all scales when simply allowing for a Hernquist stellar profile and a median concentration just 30% higher than the reference numerical value (middle panel). A higher value of $\mu=7$ (right panel) further increases the dark matter density at intermediate scales without the need to invoke any increase in median concentration, but the match to the data in this case is poorer.

In summary, the present available lensing data for galaxies at the center of cluster-scale halos (Newman et al., 2013a), tend to favor models with profiles in their inner stellar component at $r\lesssim 3-5$ kpc, flatter than what predicted by extrapolations of their outer Sérsic profiles. The same data also tend to align better with models characterized by a NFW profile with a core and a median concentration in the dark matter component twice the one predicted by dark matter-only N-body simulations, or by an Einasto profile with $\mu\sim 4$ and a dark matter concentration just $\lesssim 30\%$ higher than the simulation results. The latter results are in nice agreement with Auger et al. (2013) who also found 26 group- and cluster-scale strong gravitational lenses to be only $\sim 26\%$ more overconcentrated than similar-mass halos from dark matter simulations. We stress that the relatively small comparison observational sample prevent us to generalize our findings to all galaxy clusters. The idea that Einasto models may be more accurate representations of the mass density profiles of especially the most massive halos is in line with a number of analytic and numerical studies (e.g., Zhao et al., 2003; Gao et al., 2008; Lapi & Cavaliere, 2011; Cen, 2014; Nipoti, 2015; Angulo et al., 2016; Klypin et al., 2016; Ludlow & Angulo, 2017). The common theme in these works is the connection between the final shape of the system, and thus of its parameter $\mu$, with its mass accretion history, in a way that $\alpha=1/\mu$ and “curvature” tend to increase with fast growth, as it might be the case, on average, for the inner regions of very massive halos.

4.2. Dark Matter Fractions

Up to this point, we have mainly discussed the outputs of semi-empirical models with respect to $\gamma^{\prime}$, but clearly additional observables can and must be used to place tighter constraints on galaxy evolution. In particular, scaling relations with velocity dispersions (e.g., Dutton et al., 2013), evolution of the fundamental plane (e.g., Treu et al., 2001; Shankar & Bernardi, 2009; Shankar et al., 2013; Oldham et al., 2017), evolution of the virial relations (e.g., Peralta de Arriba et al., 2014), building up of age/metallicity gradients (e.g., Montes et al., 2014), are all relevant issues to test models.

While we aim to explore at least some of these additional issues in separate work, we here begin by comparing in Figure 9 our predicted 3D ratio of dark to total (dark plus stellar) mass as a function of scale (radius normalized to the effective radius) for all galaxies with $11.3<\log M_{\rm star}/M_{\odot}<11.5$ and $6<R_{\rm e}/{\rm kpc}<7$ (left panel) and $\log M_{\rm star}/M_{\odot}>12$ (right panel), compared with the lensing-based dark matter fractions derived from the Sonnenfeld et al. (2015, thick blue lines, left panel) sample and Newman et al. (2013a, thick solid and dotted red lines, right panel) galaxy samples. Individual dark matter fractions in the Newman et al. (2013a) cluster central galaxies are calculated by integrating the density profiles of their Figure 3. We then take the mean of the dark matter fractions at each radius to produce the solid and dotted thick red lines reported in the right panel of Figure 9. The black solid, orange long-dashed, and purple dot-dashed lines show predictions for the uncontracted NFW profile, with adiabatic contraction, and adiabatic expansion, respectively. Both panels show predictions for our standard stellar Sérsic and dark matter NFW profiles. It is apparent that at cluster scales (right panel) the dark matter fractions rise relatively more rapidly than for lower mass galaxies in the lower mass halos (left panel), irrespective of the details in the stellar and/or dark matter profiles. This trend is mostly induced by the fact that more massive galaxies, characterized by progressively higher effective radii, naturally will include more dark matter mass in their central regions.

Our semi-empirical models, within the uncertainties, align with the observationally-based dark matter fractions, and we verified this broadly continues to hold true even with flat inner stellar profiles, and/or with cores in the dark matter component. In this respect, dark matter fractions appear to be less constraining for models, though adiabatic contraction or expansion continue to be disfavored by the data irrespective of the exact inputs in the semi-empirical models, in line with Oguri et al. (2014). Our predicted dark matter fractions are in full agreement with the $f_{\rm dm}(R_{\rm e})\sim 20-30\%$ values inferred by Cappellari et al. (2013, see also ) for $M_{\rm star}\sim 2-3\times 10^{11}\,M_{\odot}$, or on larger scales with Alabi et al. (2017), who infer in the SLUGGS survey $f_{\rm dm}(<5R_{\rm e})\gtrsim 60\%$ from globular cluster kinematics data. At larger masses, our models predict $f_{\rm dm}(R_{\rm e})\sim 60-70\%$ for $M_{\rm star}\gtrsim 10^{12}\,M_{\odot}$, consistent with the extrapolation of the Cappellari et al. (2013) analytic fit, $f_{\rm dm}(R_{\rm e})\sim 1.3+0.24(\log M_{\rm star}-10.6)^{2}$, to the most massive galaxies in our mocks.

For completeness, to visualize the full dependence of the mass density profiles on the relative contributions between stellar and dark matter components, in Figure 10 we plot $\gamma^{\prime}$ as a function of dark matter fractions $f_{\rm dm}(<R_{\rm e})$ within one effective radius for galaxies with $M_{\rm star}\gtrsim 3\times 10^{11}\,M_{\odot}$ as predicted by our reference semi-empirical model with uncontracted NFW (solid and dotted lines), with adiabatic contraction (orange, dashed line), and with adiabatic expansion (purple, dot-dashed line). All models predict a steady drop of $\gamma^{\prime}$ with increasing dark matter fraction, as expected given that, as discussed with respect to Figure 9, larger $f_{\rm dm}(<R_{\rm e})$ tend, on average, to be linked to larger effective radii with more prominent dark matter contributions, and thus flatter $\gamma^{\prime}$, also in line with Figure 5. However, models with adiabatic expansion tend to predict a substantially faster drop with $f_{\rm dm}(<R_{\rm e})$ than the other two models. Our uncontracted NFW predictions nicely match those from the dynamical modeling on ATLAS3D and SPIDER galaxies by Tortora et al. (2014a, empty squares with error bars) who also assume an underlying NFW profile. On the other hand, while model predictions agree with the dark matter fractions extracted from Sonnenfeld et al. (2015, blue squares), they are highly inconsistent with the dark matter fractions derived from Newman et al. (2013a, red circle) for galaxies in clusters. As already emphasized for Figure 5, galaxies characterized by similar stellar mass, effective radius and even dark matter fractions tend to significantly differ in their total mass density slopes, an additional sign that host halo mass and its (outer) structure is systematically different.

5. Discussion

5.1. Constraints from the dependence of $\gamma^{\prime}$ on size and stellar mass

We have carried out a detailed comparison of a set of semi-empirical models with a diverse set of observational data. The models are based on assigning galaxies, with measured stellar mass, effective radius and Sérsic index, to dark matter halos via abundance matching techniques. By design, our approach is based on a minimal number of parameters and underlying assumptions, mainly related to the underlying dark matter mass distribution, which is relatively well constrained by high-resolution N-body simulations (Prada et al., 2012; Dutton & Macciò, 2014; Diemer & Kravtsov, 2015), and the IMF, which we usually set equal to a Salpeter or even variable IMF, having shown in Figure 3 that a Chabrier IMF is ruled out by the data in massive galaxies, in line with a number of independent studies based on surface gravity-sensitive absorption lines in stellar population synthesis models (e.g., Conroy & van Dokkum, 2012; La Barbera et al., 2013; Spiniello et al., 2014), and/or stellar kinematics also coupled to strong lensing measurements (e.g., Auger et al., 2010a; Treu et al., 2010; Spiniello et al., 2011; Thomas et al., 2011; Barnabè et al., 2013; Cappellari et al., 2013; Posacki et al., 2015).

Our semi-empirical models thus represent ideal tools to probe the consistency among independent data sets, to reliably test the predictions of dark matter numerical simulations, to pin down the role of adiabatic contraction/expansion, and environment in shaping the total mass density profiles of massive galaxies. The semi-empirical models successfully reproduce the total mass density profile of massive galaxies with $M_{\rm star}\gtrsim 3\times 10^{11}\,M_{\odot}$ (Salpeter IMF), along with its dependence on stellar mass and effective radius, at least in non-cluster environments. We find that even moderate modifications to the halo profile in terms of adiabatic contraction/expansion are disfavored at high significance. In the next steps (Shankar et al. in prep.) we will mainly focus on the evolution of $\gamma^{\prime}$ as a function of redshift, as well as the additional constraints that can be gained by a comprehensive Jeans-based calculation of velocity dispersions (e.g., Chae et al., 2014).

Our results on the dependence of $\gamma^{\prime}$ on galaxy properties are in line with several previous works. We discuss some of these here below.

Xu et al. (2016) have recently compared the properties of simulated early-type galaxies in the Illustris hydrodynamic cosmological simulation with the same SLACS data as in this paper. Their simulations broadly align with observations, suggesting a flattening of the total mass density profile with increasing effective radius, though generally shallower than what the data suggest.

Poci et al. (2017) have also recently confirmed from ATLAS3D that the slopes of the total mass density profiles are approximately isothermal, and depend on various galactic properties, most significantly on mass surface density, while the correlation with velocity dispersion becomes negligible above $\log\sigma/{\rm km\,s^{-1}}\gtrsim 2.1$.

Tortora et al. (2014a) have explored the density slope dependence in early-type galaxies over quite a broad range of stellar mass finding clear signs of non-universality. They suggested that the relation with effective radius is the most robust to changes in the underlying modeling. They find, in fact, that the largest galaxies have a density slope close to isothermal, while the smallest ones are characterized by a steeper profile, in line with our findings.

Dutton & Treu (2014) have also addressed the question of dependencies of $\gamma^{\prime}$ on galactic properties, though adopting only one stellar mass-halo mass mapping. In agreement with our work, they disfavor models with significant adiabatic contraction or expansion, and also predict an overall non-universal $\gamma^{\prime}$ dependence on galactic properties.

Auger et al. (2010b) compared data on 53 massive early-type lensed galaxies with a variety of semi-empirical models with varying baryon efficiency, adiabatic contraction, and stellar initial mass function. They highly disfavored light IMFs, such as a Chabrier of Kroupa, as well as strong adiabatic contraction, in nice agreement with the results presented in this work. Their models also suggest that only mass-dependent IMFs could be consistent with no adiabatic contraction, a conclusion not necessarily supported by our findings (see, e.g., Figure 13).

Jiang & Kochanek (2007) also carried out semi-empirical galaxy-halo modeling treating the stellar mass-to-light ratio as a free parameter. They analyzed surveys of gravitational lenses with velocity dispersion measurements assuming combinations of Hernquist and NFW profiles. Their total mass density slopes are peaked around the isothermal values with non negligible variations around the mean. They also found that models without adiabatic contraction favor a stellar fraction of around 3%, in agreement with weak lensing estimates (e.g., Mandelbaum et al., 2006) and the most recent abundance matching results by Shankar et al. (2014b) considered in this work (Figure 1).

Dutton et al. (2011) modeled the scaling relations of early and late-type galaxies, finding that massive galaxies characterized by a Salpeter IMF are not consistent with halo contraction. Napolitano et al. (2010) also emphasized that that the lower central dark matter fractions inferred in the oldest early-type galaxies are inconsistent with adiabatic contraction.

By analyzing a suite of simulated spheroids formed out of binary mergers, Remus et al. (2013) concluded that the exact value of the density slope depends on the amount of gas involved in the merger, with more dissipative gas-rich mergers producing steeper profiles. This would naturally imply some non-homology, as more dissipative mergers would produce more compact galaxies at fixed stellar mass with steeper slopes, in line with what observed. Dissipative processes may however not be the only cause for the dependence of $\gamma^{\prime}$ on galactic properties. Indeed, as discussed in Section 3.2, in basic $\Lambda$CDM models with uncontracted halos, the trend of $\gamma^{\prime}$ with effective radius is simply controlled by the rate of progressive change from the steeper stellar component to the flatter dark matter one.

Remus et al. (2016, see also ) also recently emphasized that the decreasing importance of dissipative processes towards lower redshifts and for more massive systems is mostly responsible in shaping the strong flattening of $\gamma^{\prime}$ with increasing dark matter fractions $f_{\rm dm}(<R_{\rm e})$ within one effective radius as measured by Tortora et al. (2014a). As detailed in Figure 10, galaxies in clusters do however show significantly smaller $\gamma^{\prime}$ distributions at fixed $f_{\rm dm}(<R_{\rm e})$, supporting a scenario in which also the larger scale environment plays a significant role in defining the total mass density profiles of early-type, massive galaxies.

5.2. Constraints from the dependence of $\gamma^{\prime}$ on halo mass

We showed that the mass-weighted total mass density slope at the effective radius $\gamma^{\prime}$ is a strong decreasing function of the host halo mass, ranging from values close to, or higher than, the canonical isothermal value of $\gamma^{\prime}$ $\sim 2$ in the field/low-mass groups, to values closer to $\gamma^{\prime}$ $\sim 1.5$ at the cluster scale. This behavior is mainly induced by the increasing contribution of the dark matter component around the effective radius, thus inducing $\gamma^{\prime}$ to progressively flatten out (cfr. Figure 2).

A closer comparison to the data interestingly reveals that for very massive galaxies with $\log M_{\rm star}/M_{\odot}\gtrsim 11.9$ in cluster-scale environments (cfr. Figure 7), the stellar profile at small scales $r\lesssim 3-5$ kpc should be significantly flatter than what predicted by an extrapolation of a Sérsic profile. This is in line with several independent lines of evidence (Postman et al., 2012; Kormendy & Bender, 2009), as mentioned in Section 4.1.

In general, several physical processes are thought to play some role in shaping the density profiles of massive galaxies, especially in overdense regions. Laporte & White (2015), for example, have presented a series of collisionless $N$-body resimulations to follow the growth of central cluster galaxies and their dark matter halos since redshift $z\lesssim 2$. In line with what discussed in this work, Laporte & White (2015) report a steepening at small scales in the stellar component not mirrored in the data. They propose that the action of repeated super-massive black hole mergers could create enough large cores of a few kiloparsecs to alleviate the tension between models and observations (see also, e.g., Merritt, 2006).

However, other in-situ mechanisms may be able to induce a permanent flattening of the central matter distributions, such as the action of efficient feedback (e.g., Peirani et al., 2008; Martizzi et al., 2013) from central active galactic nuclei (AGN). On the other other hand, an “in-situ” process such as AGN feedback, a phenomenon more and more supported by direct and indirect observations (e.g., Shankar et al., 2016, 2017; Barausse et al., 2017, and references therein), should produce a stellar core in all massive galaxies, irrespective of their environment.

Laporte & White (2015) also discuss that multiple dissipationless mergers can effectively flatten out the inner distributions, with the dark matter component producing a shallower central cusp, and gravitational heating redistributing material in the outer regions (see also, e.g., El-Zant et al., 2004; Nipoti et al., 2004; Tonini et al., 2006; Lapi & Cavaliere, 2009; El-Zant et al., 2016). This type of evolution is in line with our proposed modifications presented in Figure 7 of a cored NFW profile to faithfully align model predictions to observations in the range $3<R/{\rm kpc}<20$.

Finally, we note that other non-baryonic processes mainly linked to the dark sector, such as self-interacting dark matter, may contribute to the diversity in total mass density profiles inferred in this work (e.g., Kaplinghat et al., 2016; Di Cintio et al., 2017).

Chae et al. (2014) fitted two-component models to local SDSS nearly spherical galaxies. Adopting a velocity dispersion-dependent variable IMF, they found a mean slope for the total mass density profiles of $\langle\gamma_{e}\rangle=2.15\pm 0.04$, computed in the radial range $0.1R_{\rm e}<r<R_{\rm e}$, in full agreement with our findings. Chae et al. also went forward in testing the dependence on galactic properties including effective radius and halo mass. Their least-square fit total mass density slope between $0.1R_{\rm e}$ and $R_{\rm e}$ also decreases with effective radius, though more weakly than the one reported in Figure 5. Chae et al. (2014) also infer a dependence of the total mass density slope on host halo mass weaker than the one measured here, though they also show that a much stronger dependence is naturally obtained when the slope of the mass density profile is fitted on scales larger than their effective radius.

6. Conclusions

In this work we have set up a semi-empirical approach to create large catalogs of mock galaxies and host dark matter halos at $z=0.1$. To this purpose, we have selected a large sample of galaxies from SDSS with a probability $P(E)>0.85$ of being ellipticals, and with measured stellar masses, effective radii, and Sérsic indices. We then assigned host dark matter halos to each galaxy in the sample making use of the most recent stellar mass-halo mass relations (Shankar et al., 2014b), as well as exploring the impact of different mappings. We finally computed for each galaxy the total mass density profile and compared with a variety of data sets from strong and weak lensing as cataloged and analyzed by Sonnenfeld et al. (2013), Newman et al. (2015), and Newman et al. (2013a), and also, where appropriate, with dynamical measurements from Cappellari et al. (2015). Our main results, for galaxies with mass above $M_{\rm star}\gtrsim 2-3\times 10^{11}\,M_{\odot}$ (Salpeter IMF), can be summarized as follows:

• •

In line with observational evidence, the semi-empirical models naturally predict a non-universal mass weighted total density profile at the effective radius $\gamma^{\prime}$, independently of the exact input assumptions. In the specific, the models confirm the observational evidence that $\gamma^{\prime}$ decreases for increasing galaxy size at fixed stellar mass, and for lower mass galaxies at fixed effective radius.

• •

The strong lensing data disfavor at $\gtrsim 3\,\sigma$ a Chabrier IMF and are consistent with a Salpeter or variable IMF.

• •

Data also disfavor at $\gtrsim 2-3\,\sigma$ departures from a Sérsic stellar profile and an uncontracted NFW dark matter profile in the inner regions.

• •

Our modeling predicts the total density profile to be roughly isothermal or steeper (\gamma^{\prime}$$\gtrsim 2) in low-mass halos of $M_{\rm halo}\lesssim 10^{13}\,M_{\odot}$, but shallower (\gamma^{\prime}$$\sim 1.5) in cluster environments, with $M_{\rm halo}\gtrsim 3\times 10^{14}\,M_{\odot}$. This is mainly a manifestation of structural non-homology, with the relative density distribution of stars and dark matter varying systematically from isolated galaxies to central galaxies in clusters. Despite this substantial decrease with increasing halo mass, the predicted $\gamma^{\prime}$ are still significantly steeper than the total mass density slopes currently measured in BCGs, which suggest values closer to $\gamma^{\prime}$ $\sim 1.1$.

• •

Adequately reproducing the full total mass density profiles of massive galaxies at the center of massive halos with $M_{\rm halo}\gtrsim 3\times 10^{14}\,M_{\odot}$ requires a stellar mass profile flatter than Sérsic within $r\lesssim 3-5$ kpc. An improved match at intermediate scale $10\lesssim r\lesssim 300$ kpc can instead be obtained by assuming: 1) either a cored NFW dark matter profile, coupled to a median halo concentration a factor of $\sim 2$ higher than the one predicted by N-body numerical simulations; 2) or an Einasto profile with $\mu\sim 4$ ($\alpha\sim 0.25$), coupled to a median halo concentration only a factor of $\lesssim 1.3$ higher than numerical predictions.

Acknowledgments

F.S. acknowledges valuable discussions with A. Newman, C. Kochanek, A. Dutton, and S. Dye.

R.G. acknowledges support from the Centre National des Etudes Spatiales (CNES). T.T. acknowledges support from the Packard Foundation through a Packard Research Fellowship. We thank two anonymous referees for several useful comments and suggestions that greatly improved the quality and presentation of our results. K.-H. C. was supported by Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (NRF-2016R1D1A1B03935804).

Appendix A The concentration-mass relation

The reference concentration-halo mass relation we always use in this paper has been derived by Diemer & Kravtsov (2015). Their results are shown in Figure 11, for our chosen cosmology (Section 1), at $z=0$ (solid line), with the gray area marking the $1\sigma$ uncertainty region of 0.16 dex, and at $z=1$ (dotted line).

We compare the results by Diemer & Kravtsov (2015) with the analytic model put forward by Bullock et al. (2001) and further revised by Macciò et al. (2008). The latter first look for the redshift of collapse $z_{c}$ of a given halo as the redshift at which the characteristic mass is equal to a fraction $F=0.01$ of the halo mass at the observation redshift $z$. It then computes the concentration for halos defined as 200 times the critical density as

[TABLE]

where $H^{2}(z)=H_{0}^{2}[\Omega_{\Lambda}+\Omega_{m}(1+z)^{3}]$, and we set $K_{200c}=4$ in line with the numerical results by Dutton & Macciò (2014).

The Bullock et al. (2001) and Macciò et al. (2008) analytic model (red and purple long-dashed lines) is in very good agreement with the Diemer & Kravtsov (2015) results, but mostly at low redshifts, while it departs from it at increasing redshifts (dotted line; see also Klypin et al. 2016).

Appendix B The Sérsic 3D profile

Given a stellar mass, a projected effective radius $R_{\rm 2D}$, and a Sérsic index $n$ assigned to a SDSS galaxy from photometric image analysis, we work out the full 3D Sérsic profile following, in the specific, the analytic fits by Lima Neto et al. (1999), although other groups have provided useful approximations (e.g., Prugniel & Simien, 1997; Mamon & Łokas, 2005a). We first convert the input, empirical 2D radius to a 3D one using their approximation

[TABLE]

with $m=1/n$. The 3D mass density profile at any radius $r$ reads then as

[TABLE]

where

[TABLE]

and $x=r/a$, $\ln(a)=\ln(R_{\rm 2D})-k$, and $k=[0.6950-\ln(m)]/m-0.1789$. The normalization $f_{0}$ is then fixed by combining Eqs. 9 and 18 in Lima Neto et al. (1999) as

[TABLE]

The cumulative mass $M(<r)$ within a given radius $r$ can be analytically expressed as

[TABLE]

with $\gamma$ here indicating the (lower) incomplete gamma function.

An example of the Sérsic profile for different values of the Sérsic index $n$ at fixed stellar mass and effective radius is given in Figure 12, where we also compare it to the Hernquist profile.

Appendix C Alternative models

We here present a series of measured probability distribution functions as extracted from different semi-empirical models and compared with the data by Sonnenfeld et al. (2013). Each Figure has the same format as Figure 3, with the reference Salpeter IMF model, with stellar de Vaucouleurs and dark matter uncontracted NFW profiles, compared with a model with a variable IMF (Figure 13), a Hernquist stellar profile (Figure 14), a core in the dark matter profile as in Newman et al. (2013a) of size $r_{\rm core}=14$ kpc (Figure 15), a combination of Hernquist plus cored ($r_{\rm core}=14$ kpc) dark matter profiles (Figure 16), or assuming a reference model characterized by the Moster et al. (2013) stellar mass-halo mass relation (Figure 17). Most (but not all) of these systematic variations to the reference model tend to be disfavored by the data at different significance.