On the power spectrum of dark matter substructure in strong gravitational lenses
Abstract
Studying the smallest selfbound dark matter structure in our Universe can yield important clues about the fundamental particle nature of dark matter. Galaxyscale strong gravitational lensing provides a unique way to detect and characterize dark matter substructures at cosmological distances from the Milky Way. Within the cold dark matter (CDM) paradigm, the number of lowmass subhalos within lens galaxies is expected to be large, implying that their contribution to the lensing convergence field is approximately Gaussian and could thus be described by their power spectrum. We develop here a general formalism to compute from first principles the substructure convergence power spectrum for different populations of dark matter subhalos. As an example, we apply our framework to two distinct subhalo populations: a truncated NavarroFrenkWhite subhalo population motivated by standard CDM, and a truncated cored subhalo population motivated by selfinteracting dark matter (SIDM). We study in detail how the subhalo abundance, mass function, internal density profile, and concentration affect the amplitude and shape of the substructure power spectrum. We determine that the power spectrum is mostly sensitive to a specific combination of the subhalo abundance and moments of the mass function, as well as to the average tidal truncation scale of the largest subhalos included in the analysis. Interestingly, we show that the asymptotic slope of the substructure power spectrum at large wave number reflects the internal density profile of the subhalos. In particular, the SIDM power spectrum exhibits a characteristic steepening at large wave number absent in the CDM power spectrum, opening the possibility of using this observable, if at all measurable, to discern between these two scenarios.
I Introduction
In our Universe, structure formation based on the Cold Dark Matter (CDM) paradigm Davis et al. (1981); Blumenthal et al. (1982, 1984); Davis et al. (1985) has been extremely successful at explaining the largescale distribution of matter across cosmic times. On subgalactic scales however, assessing whether CDM provides a good fit to observations is significantly more difficult. On the one hand, baryonic processes can play an important role on these scales Brooks and Zolotov (2014); Brooks et al. (2013); Arraki et al. (2014); Oñorbe et al. (2015); Wetzel et al. (2016); Sawala et al. (2016); Creasey et al. (2017); Sawala et al. (2017); GarrisonKimmel et al. (2017), thus significantly affecting the dark matter distribution inside galaxies and their satellites, and making it difficult to compute robust theoretical predictions that can be compared to observations. On the other hand, star formation becomes increasingly inefficient in lowmass CDM halos Fitts et al. (2017); Read et al. (2017), rendering their detection and characterization within the Local Group quite challenging.
Further complicating the picture is the fact that key aspects of the particle nature of dark matter might have important consequences on these subgalactic scales. For instance, significant dark matter free streaming Bond and Szalay (1983); Dalcanton and Hogan (2001); Bode et al. (2001); Boyanovsky et al. (2008); Boyanovsky and Wu (2011) or possible interactions with relativistic species Carlson et al. (1992); Boehm et al. (2002); Ackerman et al. (2009); Feng et al. (2009); Kaplan et al. (2010); van den Aarssen et al. (2012); CyrRacine and Sigurdson (2013); Chu and Dasgupta (2014); BuenAbad et al. (2015); Chacko et al. (2016); Ko and Tang (2016) at early times can substantially reduce the number of lowmass subhalos orbiting a typical galaxy Buckley et al. (2014); Schewtschenko et al. (2015); Vogelsberger et al. (2016). In addition, dark matter selfinteraction Spergel and Steinhardt (2000); Yoshida et al. (2000); Davé et al. (2001); Colín et al. (2002) could modify the density profile of main and satellite halos Vogelsberger et al. (2012); Rocha et al. (2013); Peter et al. (2013); Zavala et al. (2013); Kaplinghat et al. (2014); Vogelsberger et al. (2016); Kaplinghat et al. (2016) away from the standard CDM prediction Navarro et al. (1997). Other dark matter particle candidates such as ultralight axions (Hu et al., 2000; Hui et al., 2017) might also lead to interesting phenomenology on small scales (see e.g. Refs. Du et al. (2017); Mocz et al. (2017)).
Disentangling the impact of dark matter physics on structure formation from that of baryons is key to probing the fundamental nature of dark matter. While it is never entirely possible to neglect the influence of baryonic structures on the evolution of the smallscale dark matter distribution (see, e.g. Ref. GarrisonKimmel et al. (2017)), it can be minimized by focusing our attention on the lowest mass subhalos present in galaxies. As mentioned above, these small subhalos are largely devoid of stars, which makes them less susceptible to baryonic feedback effects, while their abundance and internal structure are quite sensitive to the particle nature of dark matter, making them an important laboratory to test the consistency of the CDM paradigm on small scales.
These dark subhalos could potentially be probed within the Local Group using detailed observations of tidal streams Ngan and Carlberg (2014); Carlberg (2016); Bovy (2016); Erkal et al. (2016); Bovy et al. (2017) or the motion of stars within the Milky Way disk Widrow et al. (2012); Feldmann and Spolyar (2015). Beyond our local neighborhood however, gravitational lensing is the only technique capable of probing lowmass subhalos at cosmological distances from the Milky Way. In particular, galaxyscale strong lensing systems in which a massive foreground galaxy is multiplyimaging a background source (such as another galaxy or a quasar) constitute ideal environments to study the cosmological population of lowmass dark subhalos. For instance, the study of fluxratio anomalies in strongly lensed quasars Mao and Schneider (1998); Metcalf and Madau (2001); Chiba (2002); Metcalf and Zhao (2002); Dalal and Kochanek (2002a); Keeton et al. (2003); Metcalf and Amara (2012) has lead to a measurement of the typical abundance of mass substructures within lens galaxies Dalal and Kochanek (2002b), and has also been used to put constraints on the position and mass of potential individual subhalos within the lens galaxies Fadely and Keeton (2012); Nierenberg et al. (2014, 2017).
Individual mass substructure can also be detected by carefully examining the surface brightness variation of extended lensed arcs and rings. This direct “gravitational imaging” Koopmans (2005); Vegetti and Koopmans (2009a) has lead to the statistically significant detection of a few mass substructures with masses above Vegetti et al. (2010a, b, 2012). A somewhat similar technique using spatially resolved spectroscopic observations of gravitational lenses Moustakas and Metcalf (2003); Hezaveh et al. (2013) has also lead to the direct detection of a subhalo Hezaveh et al. (2016a). Taken together, these measurements can be used to put constraints on the subhalo mass function (see, e.g. Refs. Vegetti and Koopmans (2009b); Vegetti et al. (2014); Li et al. (2016)).
The main limitation of direct subhalo detection efforts is that only the most massive substructures lying within or very close to lensed arcs can be detected with large statistical significance. While not directly detectable, smaller mass substructures or those lying further away from lensed images could still potentially lead to observable effects on the lensing signal, especially on the relative arrival time delay between lensed images Keeton and Moustakas (2009); CyrRacine et al. (2016), but also on extended arcs. For instance, Refs. Birrer et al. (2017); Brewer et al. (2015); Daylan et al. (2017) have recently proposed statistical techniques to harness the constraining power from these marginal detections on the properties and abundance of dark matter subhalos within lens galaxies.
Within the CDM paradigm, the subhalo mass function is expected to rise rapidly toward smaller masses Springel et al. (2008), implying that typical lensed images could be perturbed by a fairly large number of unresolved lowmass substructures. In this limit, it becomes somewhat impractical to phrase the perturbations to lensed images in terms of individual subhalos. A more fruitful approach in this case is to describe the substructure convergence field in terms of its point correlation functions. For CDM, the large number of smallmass subhalos contributing to the total substructure convergence field implies that the statistics of the latter should be nearly Gaussian. In this case, we expect the twopoint correlation function (or its Fourier transform, the power spectrum) to dominate the statistical description of the substructure field. This last point was put forth in Ref. Hezaveh et al. (2016b) to motivate an exploratory study of the detectability of the substructure convergence power spectrum within lens galaxies using the Atacama Large Millimeter/submillimeter Array (ALMA). In practice, given that strong lensing is probing the matter density field deep in the nonlinear regime, we do not expect the substructure density field to be entirely Gaussian. Nevertheless, measuring the substructure power spectrum might still lead to important insights about the abundance and internal structure of subhalos within lens galaxies.
Interestingly, Ref. Hezaveh et al. (2016b) showed that it is possible, in principle, to measure the substructure convergence power spectrum by looking at the correlations of lensed image residuals, once a model image obtained from a purely smooth lens potential is subtracted from the data. They further showed that deep observations of strong gravitational lenses with ALMA could lead to 3 detection of the nonvanishing amplitude of the substructure power spectrum (at least if there is abundant substructure, which is the case in CDM). Given that such measurements might be possible in the near future, the immediate question that comes to mind is: What will we learn about lowmass subhalos from measuring the substructure convergence power spectrum?
In this paper, we present some muchneeded answers to this question. Using the standard halo model Cooray and Sheth (2002) as our framework, we first develop a general formalism to compute the power spectrum of the convergence field on the lens plane due to substructure. We extend the initial approach presented in Ref. Hezaveh et al. (2016b) to include subhalo populations that are not necessarily isotropic and homogeneous, and also take into account the 2subhalo term. This formalism is developed in a way that makes it easy to change the statistical properties of the population as well as the intrinsic properties of subhalos, in order to facilitate its application to different dark matter scenarios. As an example we apply it to two different subhalo populations: one in which subhalos are modeled as truncated NavarroFrenkWhite (NFW) halos as would occur in standard CDM, and another one in which they are modeled as truncated cored halos as would happen in the presence of selfinteracting dark matter (SIDM). We choose the latter because of there is evidence of cored density profiles in at least some of the Local Group satellites (e.g. Refs. Gilmore et al. (2007); Walker and Peñarrubia (2011); Amorisco et al. (2014)). We then use these two examples as a springboard to discuss how the internal structure, statistical properties, and abundance of lowmass subhalos affect the shape and amplitude of the substructure convergence power spectrum.
This paper is organized as follows. In Sec. II we present our halo modelbased formalism to compute the substructure convergence power spectrum from first principles. In Sec. III we apply this formalism to study the 1 and 2subhalo contributions to the substructure power spectrum from a population of truncated NFW subhalos. In Sec. IV we turn our attention to the substructure power spectrum in the presence of a population of truncated cored subhalos, highlighting along the way the differences from the NFW case. We finally discuss our findings and conclude in Sec. V.
Ii Substructure statistics within the halo model
We work within the framework of the halo model Cooray and Sheth (2002), where all the dark matter is bound in roughly spherical halos. Within this model, the dark matter content of a typical lens galaxy is comprised of a smooth dark matter halo containing most the galaxy’s mass, as well as a certain number of subhalos orbiting within the smooth halo. In the following, we will be concerned with these subhalos.
ii.1 Preliminaries: Subhalo statistics
We work in projected twodimensional (2D) space, with denoting the projected 2D vector in the plane of the sky. The total convergence at a given point on the lens plane is
(1) 
where denotes the contribution from the smooth lens model (dark matter + baryons) and denotes that from the subhalos. Note that the convergence is nothing more than the projected mass along the line of sight in units of the critical density for lensing, , where depends on the angular diameter distance between the observer and the source , the observer and the lens and the lens and the source :
(2) 
Here, is the gravitational constant and the speed of light.
The convergence is also related to the projected Newtonian gravitational potential via the Poisson equation: . According to the standard CDM model, a typical lens galaxy will contain a large population of subhalos, all of which contribute to as:
(3) 
where and are the convergence and the position of the th subhalo, respectively, is the total mass of the th subhalo, and the ’s are sets of parameters that determine the internal properties of the th subhalo. is the total number of subhalos contributing to the lensing convergence at position . Note that in Eq. (3) we have taken advantage of the fact that the overall contribution of the subhalo population is equivalent to the sum of the effect of each subhalo, which follows from the linearity of Poisson’s equation. Since the convergence profile of a subhalo is always directly proportional to the subhalo mass , it is useful to define . The advantage of this notation is that obeys a very simple normalization condition
(4) 
independent of the value of . Here, the integral runs over the whole lens plane.
In general, it is impossible to know the mass, position, and internal properties of every subhalo within a lens galaxy. Instead, we would like to determine the “ensembleaveraged” properties of gravitational lensing observables given the statistical properties of subhalos, such as their mass function and spatial distribution. We shall denote by the ensemble average of quantity over all possible realizations of the subhalo density field within a lens galaxy. On the other hand, the notation will be used to denote the “spatial” average of over a given area of the lens plane.
Let us assume that all the statistical properties of subhalos within a lens galaxy are captured by a probability distribution function . It is, in general, a very good approximation (see Refs. Springel et al. (2008); Fiacconi et al. (2016)) to assume that the mass and projected position of a subhalo are uncorrelated. This allows us to write the overall distribution as a product of a mass and position probability distributions as follows:
(5) 
where we have taken into account that the intrinsic properties of a given subhalo likely depend on its mass and position within the lens galaxy. The distribution contains all the information about the projected spatial distribution of subhalos within the host galaxy. Given a projected number density of subhalos, the probability of finding a subhalo within an area centered at position is
(6) 
where is the area of the lens plane where we have sensitivity to substructures (see below). The denominator in Eq. (6) is just the total number of subhalos within the area
(7) 
where is the average number density of subhalos averaged over the whole area . It is useful to write the subhalo number density as
(8) 
where is a stochastic random variable with . Here, the field describes the fractional excess probability (compared to ) of finding a subhalo at position . While any choice of fully specifies the probability density function as per Eq. (6) statistically independent, we will, in general, be interested in ensembleaveraging over realizations of the field.
Numerical studies Springel et al. (2008); Fiacconi et al. (2016) indicate that the 3D spatial distribution of subhalos near the central part of the host has a rather weak radial dependence. Taking into account projection effects and the fact that galaxyscale strong lensing is mostly probing a small region near the projected center of the host, it is usually an excellent approximation to take constant.
The subhalo mass probability distribution can be written as
(9) 
where is the standard subhalo mass function. While our results are easily generalizable to any choice of mass function, we restrict ourselves to a power law mass function, , for . In the following, we assume that is normalized such that
(10) 
As in most lensing calculations in the literature, the calculations presented in the remainder of this paper assume that each subhalo represents an independent draw from the probability distribution. We emphasize though that this does not mean that we neglect spatial correlations between subhalos; these are fully encoded in our choice of . In this case, the probability distribution describing the properties of the whole subhalo population can be factored out as a product of the probability distribution for single subhalos
(11) 
We now have all the ingredients to perform ensemble averages over all possible realizations of a subhalo population.
ii.2 Ensembleaveraged substructure convergence
It is instructive to first compute the mean ensembleaveraged substructure convergence on the lens plane . It is given by
(12)  
where we used the fact that every term in the sum in Eq. (3) contributes equally to . The result is not surprising since it just states that the average convergence for the whole population of (statistically independent) subhalos is just times the average convergence of a single subhalo. Next, we note that the integral above is nothing more than the convolution of the subhalo density profile with the spatial distribution . Using the general result for the integral of a convolution,
(13) 
we obtain,
(14) 
where we used Eq. (4). In the above, we have introduced the notation
(15) 
to denote the average subhalo mass. We note that Eq. (II.2) is useful to relate and to the physically relevant quantities and .
ii.3 The power spectrum of the convergence field
We now turn our attention to the computation of the twopoint correlation function of the substructure density field, or its Fourier transform, the substructure power spectrum. We emphasize that we do not assume here that the substructure convergence field is necessarily Gaussian. As such, we do not expect the power spectrum to characterize the substructure density field completely, and expect higherpoint correlation functions to also contain nontrivial information. Nevertheless, the rapidly rising subhalo mass function toward the lowmass end in CDM models ensures that Gaussianity is a good first approximation CyrRacine et al. (2016). Importantly, the main contributors of nonGaussianities to the substructure field are the most massive subhalos within the lens galaxy Hezaveh et al. (2016b). Since we expect them to be directly detectable Vegetti et al. (2012, 2010b); Hezaveh et al. (2013, 2016a), we can limit their influence on the statistics of the field by absorbing the most massive subhalos within the macrolens mass model .
To obtain a general expression for the substructure power spectrum , we first compute the lens planeaveraged connected twopoint correlation function of the substructure convergence field . To simplify the derivation and avoid clutter, we first focus exclusively on performing the spatial averages encoded in the probability distribution . The averages over the subhalo mass and internal properties will be restored at the end of the calculation. The substructure convergence twopoint function takes the form
(16)  
Substituting Eq. (3) in the above and using the normalization condition given in Eq. (10), we obtain
(17) 
The first term arises from ensembleaveraging over the spatial distribution of a single subhalo (the “1subhalo” term), the second term arises from averaging over pairs of distinct subhalos (the “2subhalo” term), while the last three terms ensure that we are computing only the connected part of the twopoint function. In the language of the halo model, the 1subhalo term refers to particles or mass elements within a same subhalo, while the 2subhalo term is due to those in distinct subhalos. The 1subhalo term is nothing else than the convolution of the subhalo density profile with itself
(18) 
The 2subhalo contribution contains identical terms which have the following form Sheth and Jain (2003)
(19) 
Using Eqs. (6) and (8), the convolution of the subhalo’s spatial distribution is
(20) 
where we have identified the twopoint subhalo correlation function , which encodes spatial correlation between pairs of distinct subhalos. Finally, the three last terms of Eq. (II.3) all have the same form and lead to a net contribution of . The connected twopoint correlation function of the substructure convergence field thus takes the form
(21)  
Noting that some of the integrals not involving in the second term exactly cancel the third term, we are left with
(22)  
The first two terms correspond to the 1subhalo and 2subhalo terms, respectively, while the last term, suppressed by an extra factor of , corresponds to the shot noise term, which only becomes important if the number of subhalos within the area of interest in the lens plane is small.
It is now straightforward to compute the convergence power spectrum by Fourier transforming Eq. (22). Using the following Fourier transform conventions:
(23)  
(24) 
the convergence power spectrum takes the form
(25) 
where is the wavevector, and where we have used the convolution theorem to perform the Fourier transform. We note that the independent part of the last term in Eq. (22) contributes an unobservable zeromode, which we dropped in the above. Here, is the Fourier transform of the subhalo twopoint correlation function . In the remainder of the paper we neglect the term in Eq. (II.3).
Up to this point, the only assumptions underpinning our calculation of the substructure convergence power spectrum are the statistical independence of each subhalo within a lens galaxy, and the fact that the subhalo internal properties do not depend on the subhalo position . We now introduce two simplifying assumptions:

We take the subhalo convergence profile to be circularly symmetric, implying that .

We assume that the subhalo twopoint correlation function is homogeneous and isotropic, hence leading to .
Here, . While subhalos are generally triaxial, projection effects and ensembleaveraging over all possible orientations and sizes of the subhalos’ ellipticity imply that the average convergence profile is close to circularly symmetric, hence our first assumption. Our second point amounts to assuming that the small area of the lens plane probed by strong lensing images is typical of other nearby lines of sight. With these assumptions, the Fourier transform of the subhalo convergence profile is
(26) 
where is the 0th order Bessel function.
The last step of the calculation is to reinstate the averages over subhalo mass and internal properties. We can write the total substructure convergence power spectrum as the sum of the 1subhalo and 2subhalo terms,
(27) 
where the 1subhalo term takes the form
(28) 
(the subscript has been dropped since it is now superfluous) and the 2subhalo term takes the form
(29) 
The amplitude of the 1subhalo term is approximately given by , where the quantity has been referred to as the “effective mass” in the lensing literature Refsdal and Stabell (1991); Neindorf (2003); Keeton (2009). This specific mass scale constitutes the primary dependence of the substructure power spectrum on the subhalo mass function, so we expect it to be one of the most constrained quantities with actual observations. The amplitude of the 1subhalo term can be approximated as . For a typical gravitational lens with Dalal and Kochanek (2002b), , and (given our choices for the source and lens redshift), we thus expect
(30) 
for scales larger than the typical size of a subhalo. On the other hand, the amplitude of the 2subhalo term is approximately , with very little dependence on the subhalo mass function. Given that typically and that can be important only on scales larger than the typical subhalo spatial separation, this term is generally subdominant compared to the 1subhalo term, except maybe on larger scales, depending on the size of .
Having derived the general expression for the lens planeaveraged substructure power spectrum, we can now apply it to realistic subhalo populations by specifying the probability distributions and the subhalo convergence profile . For definiteness, we make the following choices throughout the rest of this paper whenever we present numerical results: we assume a lens galaxy at redshift with virial mass and radius , kpc, and Einstein radius kpc. We take the source to be at .
Iii Truncated NavarroFrenkWhite subhalo population
iii.1 Characteristics of the subhalo population
In this section we compute the substructure power spectrum for a realistic population of smoothly truncated NavarroFrenkWhite subhalos. We are particularly interested in the strong lensing region, namely the region bounded more or less by the Einstein radius of the lens. Reference CyrRacine et al. (2016) performed a detailed analysis of the statistics of subhalo populations in strong lenses by looking at both the “local” (close to the Einstein radius of the host) and “distributed” (extending past the host virial radius) populations of subhalos and looking at their relative effects on lensing observables such as the lensing potential, deflection, shear and convergence. They found that the substructure contribution at a typical image position is largely dominated by the local subhalos.
The NFW density profile Navarro et al. (1997) has been found to provide a good fit to simulated CDM halos and is widely used to model the distribution of dark matter within galaxies and their satellites. This density profile (see Fig. 1) has an inner slope that goes as until it reaches the scale radius , where the slope steepens to . Formally, the NFW density profile leads to a divergent total subhalo mass. However, we expect tidal interactions to provide a finite truncation radius for a realistic subhalo orbiting within its host galaxy, hence leading to a finite subhalo mass. Here, we adopt the following truncated NFW profile (tNFW) Baltz et al. (2009) for our subhalos:
(31) 
which is also shown in Fig. 1. Here, is the threedimensional distance from the center of the subhalo and is the tidal radius. Observe that for , the density profile decays quickly as . Basically, our truncation scheme is meant to reflect that any dark matter particles outside are tidally stripped as the subhalo undergoes a full orbit within its host. The tidal radius thus evolves in time, generally getting smaller as the subhalo orbits within the tidal field of the host.
Projecting Eq. (31) along the line of sight leads to the following convergence profile for a tNFW subhalo Baltz et al. (2009)
(32) 
where
(33) 
(34) 
(35) 
The scale mass is related to the total subhalo mass via the relation Baltz et al. (2009)
(36) 
The parameter is similar to the concentration parameter, , which measures how concentrated the mass of a halo is since most of the mass is contained within . The tidal radius and virial radius are not necessarily the same however, so .
In the notation of Sec. II, the internal structure parameters for a truncated NFW subhalo are simply . Here, we adopt the following phenomenological relations between the internal structure parameters and the subhalo mass and position CyrRacine et al. (2016):
(37)  
(38) 
where we adopt below a fiducial value of (Vegetti and Vogelsberger, 2014; Duffy et al., 2008), and is a parameter that depends on the density profile of the host; for an isothermal profile , while for a subhalo outside the scale radius of an NFW host CyrRacine et al. (2016). The quantity is the threedimensional distance between the subhalo and the center of the host galaxy, and and are, respectively, the scale and truncation radii for a subhalo of mass at position . For a pivot mass , we adopt kpc (Vegetti and Vogelsberger, 2014), kpc, and kpc (Diemand et al., 2008; Springel et al., 2008).
In order to apply the result from the previous section, we need to know the distribution , which we assume can be written as
(39) 
We model the distribution for scale radii assuming that the scatter in the scale radiusmass relation Eq. (37) is normally distributed such that
(40) 
where is a Gaussian probability distribution with mean and standard deviation , and is the fractional scatter about the scale radiusmass relation given in Eq. (37). We take throughout the rest of this paper, but we note that this specific choice has very little impact on our results.
Noting that , where is the projection of along the line of sight and is the projection onto the lens plane, the distribution of tidal radii marginalized over can be written as
(41)  
where is the threedimensional distribution of subhalos within the lens galaxy and is a normalization factor equal to the projection integral,
(42) 
Under the assumption that the projected distribution of subhalos is uniform, the radial distribution of subhalos is simply equal to the inverse area of the strong lensing region, . The choice of to obtain this is not unique. However, in the limit that the strong lensing region is probing only a small projected area of the host lens galaxy, we can obtain a unique expression for even when the distribution of subhalos is nonuniform. As shown in Appendix A, the integral in Eq. (41) can be performed in this limit to yield
(43) 
In the following we model the host as being isothermal, for which as stated above. Note that has no dependence on the subhalo position within the host, consistent with the assumptions used in Sec. II.
Lastly, we express the mass probability distribution as a powerlaw function Springel et al. (2008)
(44) 
where and . This mass function is illustrated in Fig. 2 for different choices of . We note that the constant , which normalizes the subhalo mass function, and the average convergence are proportional to one another as per Eq. (II.2). Typical gravitational lenses have an average convergence in the range Dalal and Kochanek (2002b), so we normalize the subhalo mass function such that .
Although we do not require the convergence field to be Gaussian, nor do we assume it, we do limit the large nonGaussian contributions from the few most massive subhalos by setting an appropriate upper bound on the subhalo mass range included in our analysis. In practice, the maximum subhalo mass to include in the substructure convergence power spectrum calculation should be dictated by the data set used to measure it. Indeed, the spatial resolution, pixel size, and the signaltonoise ratio of the data specifies a subhalo mass sensitivity threshold below which a statistically significant direct detection of a subhalo is unlikely. For highquality spacebased optical data, this threshold could be as low as (Vegetti et al., 2014), while for interferometric data it could reach (Hezaveh et al., 2013). Here, we adopt a fiducial value of . The minimum subhalo mass we consider is . As we will show below, the specific choice of is largely inconsequential as long as .
iii.2 Power spectrum: 1subhalo term
We can now apply the formalism developed in Sec. II to a population of tNFW subhalos to study how the abundance, density profile, radial distribution, and subhalo sizes affect the the convergence power spectrum. In this case, Eq. (II.3) for the 1subhalo term becomes
(45) 
iii.2.1 Analytical discussion
For typical Poisson realizations of a population of spherically symmetric tNFW subhalos, we expect the behavior of the 1subhalo term to depend mostly on three quantities: a low power spectrum amplitude, a turnover scale corresponding approximately to the size of the largest subhalos, and an asymptotic high slope dictated by the small behavior of the subhalo density profile, which takes over for (defined below).
For small values, we expect the 1subhalo contribution to the power spectrum to plateau to a constant value since taking makes in Eq. (II.3), in which case and are independent. Another way to understand this low plateau is to realize that subhalos can be modeled as point masses, i.e. , on scales larger than the biggest subhalo’s truncation radius, hence leading to . With and our choice for the mass function parameters described above, we expect a low amplitude of kpc.
As is increased, the power spectrum begins probing the actual density profile of the subhalos, leading to a suppression of the power compared to the pure pointmass case. This turnover scale is determined by the truncated size of the largest subhalos, since this is the largest scale in the problem relevant to the 1subhalo term. We therefore expect that this turnover is going to occur near a scale that corresponds to the inverse of the tidal radius of the largest subhalo: .
As is further increased, the 1subhalo term probes the intermediate scales between the typical truncation and scale radii of the tNFW subhalo population. Finally, we expect the convergence power spectrum to asymptote to a powerlaw behavior at large where it is probing scales deep within the NFW scale radius. This power law can be determined by finding the small limit of the convergence profile given in Eq. (III.1),
(46) 
and taking the (2D) Fourier transform, which leads to
(47) 
This implies that for . We expect the power spectrum to reach this slope at a scale below that of the smallest scale radii in the population. It is therefore useful to define the wave number beyond which the convergence power spectrum is a simple power law determined by the inner density profile of the subhalos.
iii.2.2 Numerical results
Before ensembleaveraging over and , it is informative to consider the shape of the convergence power spectrum for specific values of and . Making the following choices:
(48) 
(49) 
the 1subhalo term takes the simple form
(50) 
Note that Eq. (49) is equivalent to having a constant ratio for , which is not generally the case. From our expressions for the scale and tidal radius we expect to lie in the range , depending on subhalo mass and position.
Figure 3 3 shows the power spectrum defined in Eq. (50). Panel (a) displays the features discussed in the preceding section, which have the expected behavior. The asymptotic low amplitude is kpc and matches the amplitude of the power spectrum of a population of point masses (black) with the same mass function. The truncation scale, which for kpc is kpc (dasheddotted gray), very closely matches the scale at which the power spectrum turns over, consistent with the fact that this scale corresponds to the sizes of the largest subhalos. Furthermore, past kpc (gray) the large behavior matches a power law (dashed red), which again matches our expectation since in this regime we are within the scale radius of even the smallest subhalos i.e., where the tNFW convergence goes as Eq. (46).
In the remaining panels we vary several parameters of relevance to the power spectrum. Panel (b) shows the effect of changing the density profiles of subhalos by changing . When we increase , we are keeping and fixed while increasing , which means that the subhalo size is increasing and subhalos are becoming less concentrated toward the center. This has the effect of decreasing power on small scales and decreasing .
Panels (c) and (d) both reflect changes in the subhalo mass function: the former shows the result of varying and , and the latter, the effect of making the power law shallower. Both changes affect the low amplitude as well as the distribution of power and slope on scales larger than ; to disentangle these two effects we keep the quantity fixed while changing the mass function, which makes the low amplitude remain the same. In this manner, we can isolate the effects of the subhalo mass function on the shape of the convergence power spectrum at high . In Panel (c) we see that decreasing by an order of magnitude adds power on small scales. Indeed, removing the largest subhalos and redistributing their mass among smaller subhalos causes an increase in , which adds power on small scales. Panel (c) also illustrates the impact of increasing from to . The resulting change to the convergence power spectrum is rather small, reflecting the fact that the more massive subhalos tend to dominate the behavior of the power spectrum. This also implies that the convergence power spectrum shows little sensitivity to the lowmass cutoff of the mass function. Finally, Panel (d) shows that, by making the power law shallower, we are reducing power on small scales. To understand this effect, we refer the reader to Fig. 2, where one can see that by making the slope shallower, we are decreasing the number of lowmass subhalos and are in fact increasing the number of subhalos more massive than the pivot mass. Note that despite the change in the shape of the power spectrum on intermediate scales, the spectra still match the power law of the fiducial case at .
Having gained some intuition into how different parameters in our model affect the power spectrum, we can move on to the more general case where we perform ensemble averages over the two intrinsic subhalo parameters: and . The 1subhalo power spectrum in this case is shown in Fig. 4. The fiducial model – shown in black in both panels – corresponds to the parameter values for and , given in Eqs. (40) and (43)), (isothermal lens) and .
In each panel we show the effect of changing one of these parameters. Panels (a) and (b) reflect changes in and , respectively. It is immediately obvious from Panel (a) that changing the index has little impact on the convergence power spectrum, beside from a slight redistribution of power at intermediate and small scales. This means that the power spectrum will have limited sensitivity to the host galaxy’s density profile; on the other hand, it also means that uncertainties on the density profile of the host will not prevent the power spectrum from being an effective tool to study subhalo populations.
Panel (b) of Fig. 4 demonstrates that the power law in the scale radiusmass relation can have a significant impact on the smallscale substructure convergence power spectrum. As we increase , the minimum scale radius decreases quickly, and so increases. In fact decreases by an order of magnitude as we change from 1/4 to 1/2. This has the effect of adding power on small scales, as discussed in Sec. III.2.1.
Another natural parameter to vary would be the scatter in the scale radiusmass relationship, . However, for a scatter of or less, the impact on the convergence power spectrum is much smaller than the change associated with varying the index , and we therefore do not show it here. We also note that for a scatter larger than , the approximate model presented in Eq. (40) likely breaks down at small subhalo masses, and should be replaced by a more realistic distribution of .
We find that the 1subhalo term for a population of tNFW halos is well fit by a function of the form
(51) 
where
(52) 
(53) 
(54) 
(55) 
(56) 
As shown, the parameters are determined by the truncation, the scale radius, the mass function, and the massconcentration relation. We note that this fit works best for values of , and starts deviating from the “true” curve for higher values of . In the above, we have defined
(57) 
The fitting function is shown as a dotted green line in Panel (a) of Fig. 3.
iii.3 Power spectrum: 2subhalo term
To find the total power spectrum we have to include the contribution of the 2subhalo term, given by Eq. (II.3). As explained in Ref. Chamberlain et al. (2015), the 2subhalo term receives contributions from two distinct effects. First, subhalos have, in general, a nonuniform spatial distribution ( from Eq. (6)) due to their interaction with the potential well of their host halo. This socalled “host” contribution simply reflects the fact that subhalos can be gravitationally bound to their host lens galaxy, hence leading to a local enhancement of the convergence’s twopoint function. Second, subhalos can form selfbound groups orbiting their host galaxy. Due to tidal interactions with the latter, however, these subhalo groups are not expected to survive for more than a few dynamical times, Chamberlain et al. (2015) and we thus foresee their contribution to be subdominant. So far, this contribution to has not been measured nor extracted from simulations, at least at the mass scale of interest (see Ref. Fang et al. (2016) for a measurement on cluster scales.). Due to this, we focus below on the host contribution, but the reader should keep in mind that the subhalo group contribution should be added in order to get a fully accurate estimate of the 2subhalo term.
As an illustrative example, we choose a radial distribution of subhalos that is cored and decays as for large ,
(58) 
where kpc correponds to the core size. The total power spectrum is shown in Fig. 5, together with the individual contribution of the 1 and 2subhalo terms. On large scales, for kpc, the 2subhalo term dominates, adding power and changing the low slope from a constant to a power law. On small scales, however, the 1subhalo term dominates (as expected), and the addition of the 2subhalo term leaves the power spectrum unchanged. Note that the oscillations at small come from having nonzero over a finite region in the lens plane.
Iv Truncated Cored subhalo population
In Sec. II we applied the convergence power spectrum formalism to a population of truncated NFW subhalos, since CDM halos in simulations seem to universally have NFW density profiles. We now apply the same methodology to a population of subhalos whose density profiles approximate what we expect SIDM subhalos to look like: cored at the center and with a large behavior similar to NFW. The idea is to gauge the extent to which the power spectrum differs for NFW and cored profiles, which could be indicative of the utility of this observable in discerning between CDM and a different dark matter scenario in which halos are predicted to have cores instead, like SIDM. As we have emphasized in preceding sections, there are essentially two types of ingredients that go into the convergence power spectrum: the statistical properties of the subhalo population and the internal subhalo parameters, which determine the surface mass density profile.
With respect to the first point, SIDM Nbody simulations have shown that, at least in the case of elastic scattering with cross section cmg, the spatial distribution and number density of subhalos are largely unchanged Vogelsberger et al. (2012); Rocha et al. (2013); Peter et al. (2013); Zavala et al. (2013). Indeed, we expect that the subhalo distribution on the lens plane will be largely intact with respect to the CDM case since the volume occupied from the outskirts of the lens galaxy to the edge of its central region, where dark matter selfinteractions can play a role, is many orders of magnitude larger than the volume occupied by the host’s core itself; in fact the latter makes up about of the total lineofsight volume. Furthermore, simulations find that there is essentially no change to the subhalo mass function for moderate dark matter selfinteraction cross sections (at least down to ; refer to Fig. 6 of Ref. Vogelsberger et al. (2012) to see both of these points).
With respect to the second point, there is a stark contrast between CDM and SIDM dark matter halos due to the appearance of a central core in the latter. A common cored density profile is the Burkert profile Burkert (1995),
(59) 
where is the core radius, and the scale mass is the mass within the core. Here we set , where is a constant that represents the size of the core as a fraction of the scale radius. Furthermore, we also add a smooth truncation term, resulting in a profile of the form
(60) 
where the total mass of the subhalo with this profile is given by
(61) 
We call this a truncated Burkert (tBurk) profile. Note that for a given , the intrinsic parameters for the tBurk subhalos are the same as for the tNFW ones: . This profile is shown in Fig. 6, where we show the tNFW profile and tBurk profile for . This choice for is motivated by the fact that Ref. Rocha et al. (2013) finds that for them, in Eq. (59) corresponds to the CDM value of