Flashlights: Properties of Highly Magnified Images Near Cluster Critical Curves in the Presence of Dark Matter Subhalos

Dark matter subhalos with extended profiles and density cores, and globular star clusters of mass 106–108 M ⊙ that live near the critical curves in galaxy cluster lenses can potentially be detected through their lensing magnification of stars in background galaxies. In this work, we study the effect such subhalos have on lensed images, and compare to the case of more well-studied microlensing by stars and black holes near critical curves. We find that the cluster density gradient and the extended mass distribution of subhalos are important in determining image properties. Both lead to an asymmetry between the image properties on the positive- and negative-parity sides of the cluster that is more pronounced than in the case of microlensing. For example, on the negative-parity side, subhalos with cores larger than about 50 pc do not generate any images with magnification above ∼100 outside of the immediate vicinity of the cluster critical curve. We discuss these factors using analytical and numerical analysis, and exploit them to identify observable signatures of subhalos: Subhalos create pixel-to-pixel flux variations of ≳0.1 mag on the positive-parity side of clusters. These pixels tend to cluster around (otherwise invisible) subhalos. Unlike in the case of microlensing, signatures of subhalo lensing can be found up to 1″ away from the critical curves of massive clusters.


INTRODUCTION
Recent work has shown that high magnification near galaxy cluster critical curves makes it possible to detect individual stars in background, strongly lensed galaxies at intermediate and high redshifts (Kelly et al. 2018;Chen et al. 2019;Kaurov et al. 2019;Chen et al. 2022;Kelly et al. 2022;Welch et al. 2022;Meena et al. 2022bMeena et al. , 2023)), realizing a prediction made 3 decades ago (Miralda-Escude 1991).A recent 192-orbit HST Flashlights project (PI: Kelly 15936; Kelly et al. 2022) was specifically designed to look for such high magnification events in the six Hubble Frontier Field galaxy clusters (Lotz et al. 2017).Lensing theory has evolved rapidly in the last several years to explain and take advantage of these new exciting events (Venumadhav et al. 2017;Dai & Pascale 2021;Diego 2022;Meena et al. 2022a).
Theory predicts that on the two sides of the cluster critical curve certain image properties are different.On the side further away from the cluster center, images have the same handedness as the unlensed (unobserv-able) sources, while on the side closer to the cluster center images are mirror-imaged compared to the source.The two sides are said to have positive and negative parities, respectively, and the corresponding images are called minima and saddle points, referring to the type of extrema they represent in the arrival time, or Fermat surface (Blandford & Narayan 1986;Schneider et al. 1992).While parity cannot be observed for images of point sources like stars, studies of microlensing by intracluster stars, and possibly black holes near cluster critical curve show that magnification properties of lensed images on the positive and negative parity sides of the cluster are different: for the former, the minimum magnification is within a factor of few of the smooth cluster's magnification at that location, and can be as high as thousands, whereas the magnification of the latter can be smaller by ∼ 10−100, and so these images can remain invisible for long stretches of time (see Figures 6 and 8 of Diego et al. 2018).This property is related to the lensing theorem that minima can never be demagnified, while there is no such restriction on saddles (Schneider et al. 1992).
High magnification of these events not only helps with the study of otherwise undetectable stars, but also allows for better resolution of the lensing mass distribution in the galaxy cluster, including intra-cluster stars, and dark matter subhalos.Lensing by subhalos near cluster critical curves is an emerging subject (Dai et al. 2018(Dai et al. , 2020)), and is also the subject of this paper.Dark matter subhalos share some similarities with stars and black holes as lenses, however, there are also important differences.We focus on the asymmetry between the positive and negative parity sides of the cluster, in the presence of subhalos, and use the asymmetry to devise tests to detect subhalos, and possibly determine the mass fraction and properties of subhalos, namely the core radius of their mass distribution.In the present paper we work with lensing simulations that solve the lens equation to find images of sources, leaving the application to observations to a later paper.

THE LENS MASS MODEL
The lens model we use for our simulations consists of the cluster macro-model, and a zoom-in portion of the cluster, near its tangential critical curve, populated with subhalos.

Cluster model
The cluster is smooth and elliptical, with a small core and isothermal density profile, ∝ r 2D −1 .In Appendix B we present results that use a somewhat shallower cluster density profile, ∝ r 2D −0.9 .Such slope values are typical in the regions where multiple images are found (e.g., Jauzac et al. 2019;Ghosh et al. 2022).Its critical curve is shown as the black curve in Figure 1.We will refer to it as the cluster CC in this paper.To model the effects of subhalos we use a square modeling window (gray), which encloses two merging images (red circles) of a source depicted by a blue star.The other images of that source are also shown as red circles, but we do not use these in the rest of the paper.We populate this square with subhalos, as described later, and leave the rest of the cluster subhalo-free.
We scale our cluster model to the observed SGAS J1226+2152 cluster, at z l = 0.43, used in Dai et al. (2020).The image pair near the cluster CC that the authors analysed originates from a source at z s = 2.93.For standard ΛCDM cosmology, that implies critical surface mass density for lensing of Σ crit = 0.41g cm −2 , and 1 is 5.61 kpc in the lens plane.The distance to the merging pair from the cluster center is 7.3 .Using this scaling, the size of the modeling frame (gray square) is 2.0 across, or 11.2 kpc in the lens plane.The total projected cluster mass in the modeling window is 1.243 × 10 11 M .We note that SGAS J1226+2152 is rather low mass, so when comparing our sky-projected distances quoted in arcsec, to that of some other cluster, one needs to rescale the distances.For example, MACS 1149 has an Einstein radius for sources at a redshift similar to the one used here, of ∼ 30 , a factor of ∼ 4.2 larger.Other massive cluster would have even larger sizes (see Table 5 of Johnson et al. 2014).
For the purposes of calculations, our code length unit in the lens plane is 0.002 , and to calculate image positions and magnifications we zoom in by a further factor of 10, to 2 × 10 −4 arcsec.In the analysis to get the observable properties, we use coarser pixels, 0.04 , that are closer to the size of an HST pixel.

Sources and their images
Our source is a portion of a galaxy, ∼ 2 × 0.5 , near the blue star symbol in Figure 1, such that it fills the modeling window in the lens plane.It consists of individual stars, represented by point sources, but our numerical calculations will hold for any source no larger than ∼ 1 pc. 1 The sources are randomly scattered in the source plane, but one could also select a subset of these that would follow a correlation function, for example, the one found for nearby star forming regions.Average surface density of stars in the source plane is 0.12 stars/pc 2 .This density is lower than that of all stars in an L galaxy like the Milky Way.However, only the intrinsically brightest stars, like giants, are interesting in the present context, and their number density is lower, especially if the search is done in a specific filter that would be sensitive to stars of a given type.
We solve the lens equation, i.e., forward lens few×10 6 stars per realization, but only about a third of them end up generating at least one image in the lens plane window.Image multiplicity ranges from 1 to 3-4 within the lens plane modelling window, but the typical multiplicity is 1 or 2, and the typical number of images is 2 × 10 6 .
In some of our calculations we are interested in lensing magnification, so these calculations assume the that all stars have the same intrinsic luminosity.Later in the paper, when considering observables, we use a luminosity function for stars, dn/dL ∝ L −2.5 (Kelly et al.The merging pair of images straddling the CC that we will model in detail is on the right side of the cluster.The box around it is the size of the modeling window we use, for example in Figure 2. The lower left corner of the box defines the origin of our coordinates.Our Cartesian coordinate system is oriented such that the image stretching is very nearly along the vertical y−axis, which aides in the analysis in Section 4.

2018)
. The luminosity normalization is not important because it merely rescales all fluxes by the same multiplicative factor, while our analysis deals with flux ratios, or lensing magnifications.

Subhalos in the lens plane
We populate the small gray square in Figure 1 with subhalos.One example of the total projected density is shown in Figure 2. The κ gradient due to the main cluster is clearly visible, while subhalos generate small deviations from smooth density gradient.
Subhalos were drawn from a mass function, dn/dm ∝ m −0.9 , used in Dai et al. (2020).They cover 2 decades in mass, 10 6 M → 10 8 M .With this, the average mass of subhalos is m = 2.46 × 10 7 M .Guided by previous work (Gilman et al. 2019;Dai et al. 2020), we use mass fractions in subhalos of about 1% and 0.5%, and in some cases 0.25%.With 1% (0.5%) the number of subhalos in the modeling window of 2 × 2 is ∼ 86 (∼ 43).Note that because of the finite number of subhalos, and the stochastic nature of random draws from the mass function, one cannot guarantee that any given realization of the lens plane window will have exactly 1% or 0.5% mass in subhalos.The actual mass fraction in subhalos in our Figure 2.An example of the lens plane modeling window (gray box in Figure 1) populated with subhalos whose mass fraction within this window is 1%.The color shows the κ density.The frame is 2 on the side.The white line is the critical curve.The center of the cluster is to the left of the frame.
runs tend to be somewhat different from these values, and are indicated in relevant figures.
If subhalos were added over the whole extent of the cluster, then the surface mass density of the smooth component of the cluster would need to be reduced by the corresponding amount, to keep the total density the same.However, we add subhalos only to the modeling window (small gray square in Figure 1); adding subhalo mass fraction of 1% translates into increasing the density, and hence the deflecting mass of the cluster by ∼ 0.01%, which is negligible.This is further supported by the fact that the distance of the critical curve from the cluster center is not affected by the addition of subhalos.
The density profile of subhalos is given by eq.A3 in Appendix A. For a = 0.05 it closely resembles the Navarro-Frenk-White profile for about 3 decades in radius, but falls off more sharply at large radii, imitating tidal truncation expected for subhalos living in a dense cluster environment.The simple analytical form of its relevant lensing properties aides greatly in reducing computational cost.
We take the core radius of the model to correspond to the radius where the 3D slope is approximately "isothermal", at r ≈ 0.1 in the natural units of the profile; see Figure 14.This core radius needs to correspond to a certain physical size that we estimate from published data.Consulting the LITTLE THINGS galaxy survey (Table 2 of Oh et al. 2015), the typical radii of 10 9 M galaxies is about 300 pc, with a lot of variation between individual galaxies.This value is in line with the results of Wolf et al. (2010).Our subhalos are less massive, so these values have to be scaled down.We chose the scaling motivated by the virial condition, where the cube of the virial radius scales linearly with the virial mass; we use r core ∝ m 1/3 , even though the core radius need not scale linearly with the virial radius.The average subhalo mass is ∼ 40 times less massive than 10 9 M , so core radii should be ∼ 3.5 times smaller, or ∼ 85 pc.
If our subhalos are to represent globular star clusters, instead of dark matter subhalos, then their core radii should be smaller, of order of few-10 parsecs (Faisst et al. 2022).To span the approximate range of possible core radii, we do three sets of calculations, with 22 pc, 44 pc, and 88 pc, but show only the 22 pc and 88 pc results in most of the plots.The mass range of our subhalos, 10 6 − 10 8 M , is reasonably appropriate for the higher mass end of globular clusters.
Globular clusters will have a similar gravitational lensing effect as subhalos, but globulars would contribute to the observable light, thereby altering the distribution of flux in the lens plane.The number density of globulars in the strong lensing regions of clusters is ∼ 1/ (Lee et al. 2022;Diego et al. 2023), so our modeling window would have ∼ 4 globulars, considerably fewer than dark matter subhalos (Gilman et al. 2017;Mao & Schneider 1998).Therefore we do not consider globulars as subhalos in this paper.We assume that all subhalos are completely dark.
For each run, few×10 6 stars from the extended source (Section 2.2) were forward lensed.All were multiply imaged by the cluster as whole, and some were further multiply imaged by subhalos.The fraction of sources that were multiply imaged by subhalos depends on the density profile of subhalos.Massive subhalos with density cusps or small core radii, can generate several additional images, but for the subhalo core radii chosen here, image multiplicity above 2 is rare.Lensed images that ended up outside of the lens plane modeling window were ignored.
One realization of the modeling window is not sufficient to draw statistically meaningful conclusions, because the typical number of subhalos in any 2 region is small, ∼ 20 − 80, and given that subhalos are drawn from a mass function and their positions are randomly chosen results in considerable scatter between realizations.Therefore we generate 16 realization of each subhalo scenario.Most figures in this paper that show statistical results use 16 realizations.

THE EFFECT OF SUBHALOS ON THE DISTRIBUTION AND NUMBER DENSITY OF HIGHLY MAGNIFIED IMAGES
In this section we present the effect that subhalos have on the image magnification and distribution in the lens plane.We will be using magnifications, not fluxes, which is equivalent to assuming that all sources have the same intrinsic luminosity.While magnifications and exact positions of individual images are not directly observable, understanding them will aid in devising observable tests for subhalos.
First, a few notes about the abbreviated names we use for certain commonly occurring concepts.We will refer to the positive and negative parity sides of the cluster CC simply as +ve and -ve sides of the cluster.By cluster 'sides' we mean the two portions of the cluster on either side of the CC of subhalo-less cluster.We use this same definition for both subhalo-less and subhaloed clusters.In other words, the cluster CC always refers to that of the subhalo-less cluster.CC of subhaloed cluster often has wiggles and loops, and it would be hard to differentiate between its two 'sides'.Subhalos themselves, especially those further from the cluster CC, often have their own CC loops, and we will specifically indicate that in the text.
The +ve and -ve sides of the cluster are not to be confused with +ve and -ve parities of images.While most of the image on the +ve cluster side have +ve parity, that is not always the case in the presence of subhalos.The same is true of the -ve cluster side and -ve parity images.
In the presence of subhalos, the images on the two sides of the cluster have different properties.Several papers, especially in the recent years have investigated the effect of microlenses, i.e., point lenses near cluster CC (Diego et al. 2018;Oguri et al. 2018;Kaurov et al. 2019), and a few have investigated subhalos near cluster CC (Dai et al. 2020).In this paper we focus on subhalos with finite central densities and hence core radii, that live in clusters whose density falls with distance from the center.We obtain some analytical approximate results, emphasize the role of cluster density gradient, and subhalo core radii, and discuss the similarities and differences between subhalo lensing and microlensing by stars.

Distribution of bright images in the lens plane
The two panels of Figure 3 show a portion of the lens plane near cluster CC (gray box in Figure 1).To enhance the visual impression of the effect of subhalos we show a superposition of 16 independent realizations, each with subhalo mass fraction of ≈ 0.5%.The gray  As expected, highly magnified images outline the cluster CC, but the presence of subhalos makes the CC deviate from that of the subhalo-less case.In some places the CC is so distorted that the string of highly magnified images is significantly away from the CC of subhalo-less cluster, by ∼ 0.1 , for a cluster of mass comparable to that of SGAS J1226+2152.The total length of the CC is somewhat increased due to the presence of subhalos, within the ∼ ±0.1 band due to small perturbations from subhalos, leading to a higher number of highly magnified images.
Lens plane image distributions like Figure 3 can be summarized by plotting the distribution of images with |µ| > 1000 as a function of their distance from the cluster CC, as we do in Figure 4.The distance on the horizontal axis is from the cluster CC; negative (positive) distances are on the -ve (+ve) side of cluster.The blue (green) histograms show cases with subhalo core radii r c = 22pc (88pc).The thin (thick) lines show cases with subhalo mass fraction of 0.42% (0.9%), while the gray histogram represents the case with no subhalos, and is a much narrower distribution than either the blue or green histograms.The +ve cluster side has more highly magnified images, with some lying further away from the CC: both blue and green histograms extend to further positive than negative distances.The asymmetry is especially pronounced in the green histogram.
This asymmetry, as well as the effect of subhalo core size on the properties of image (blue vs. green histograms in Figure 4) will be discussed in Section 4.
Figure 4.The distribution of highly magnified images, |µ| > 1000, on the -ve side (negative side of the x-axis), and +ve side (positive side of the x-axis) of the cluster CC.Two different subhalo core radii (blue, 22pc, and green, 88pc), and two different subhalo mass fractions (thin, 0.42%, and thick, 0.9%) are shown.Sixteen realizations of each case were added to generate these histograms.

Magnification distribution of sources and images
Because of the universal behaviour of magnification near caustics, our results should agree with the predictions that the area in the source plane, σ sp , where sources are magnified by > |µ| should scale as σ sp ∝ |µ| −2 .The differential distribution of source magnifications, N s (|µ|)d|µ|, in the lens plane should follow the same scaling, because σ sp ∝ |µ| 1 |µ| N s (|µ|)d|µ|, where the factor of |µ| on the right hand side comes from the fact that we are translating cumulative to differential counts, and the factor of 1/|µ| accounts for the area in the source vs. lens planes.
Figure 5 shows the distribution of source magnifications, summed over all the images in our 2 modeling window.The gray histogram represents the case with no subhalos.The blue (green) histograms are for cases with subhalos of r c = 22 pc (r c = 88 pc), and two values of subhalo mass fractions are shown.The high magnification tail of the distributions, above |µ| ≈ 100, is the same for all cases, regardless of subhalo properties.This is the prediction described in the previous paragraph.The faint end of the distribution contains saddle images due to subhalos.
The magnification distribution of images is plotted in Figure 6, with the left and right panels showing the +ve and -ve cluster sides, respectively.To isolate the effect of subhalo, the images in the immediate vicinity of the cluster CC, ±0.1 away from CC, were excluded.Because at that distance the typical magnification is |µ| ≈ 100, there is a slight drop in the image number brighter than |µ| ≈ 100.If image closer to the cluster CC were included, their great number would drown out the difference between the +ve and -ve sides of the cluster that exists slightly further away from CC.
The left and right panels of Figure 6 show images on the +ve and -ve sides of the cluster.Comparing these, it is again apparent that bright images, |µ| > 100, are found predominantly on the +ve side of the cluster.
The most striking feature of Figure 6 is the difference in the number of highly magnified images generated by subhalos with large core radii (88pc) on the +ve vs. - ve sides of the cluster: compare green distributions in the left and right panels.While small core radii (22pc) subhalos generate fewer images on the -ve side compared to the +ve side, by a factor of ∼ 2 − 3, the large core radii subhalos generate almost no images above some maximum magnification on the -ve side.In this figure, that maximum magnification is |µ| = 100 (because we removed a band of ±0.1 around cluster CC).In fact, on the -ve side, the number of bright images in the case with large-core subhalos and with no subhalos is about the same (gray and green histograms in the right panel).The reason for the absence of subhalo magnified images on the -ve side of the cluster in the case of large core radii is discussed in Section 4.4.
The image magnification distributions have two distinct bumps on the faint side of the peak.The images in the fainter bump are all positive parity images, all located inside the central, smaller critical curve near the centers of many subhalos, on both +ve and -ve cluster sides.The somewhat less faint bump is made up predominantly of negative parity images inside the larger, outer critical curve of subhalos, also on both sides of the cluster CC.These faint images are not the focus of our work.
From the image distribution we see that about 10 −5 of all images are highly magnified, with |µ| > 1000.The fraction is somewhat higher for subhalos with smaller core radii (22pc vs. 88 pc), and somewhat higher for subhalo mass fraction of ≈ 0.5% compared to ≈ 0.25%.But these differences are not more than a factor of 2.

Number of highly magnified images as a function of subhalos' CC length
Highly magnified images appear almost exclusively near critical curves.Therefore the number of images highly magnified by subhalos should depend on the total length of the subhalos' CC.This is suggested by Figure 3, where the highly magnified images tend to located on the subhalo CC (gray loops).Here we quantify this visual impression, and in Section 4 we discuss the underlying reasons for why subhalo CC are longer on +ve cluster side.
Figure 7 plots the log of the total length of all subhalos' CC (in code units) on the x-axis vs. the total number of images with |µ| > 100 on the y-axis.This is done separately for +ve (square symbols) and -ve (triangles) sides of the cluster.To capture the effect of just the subhalo CC, without including image magnified with the help of cluster CC, the band of width ±0.2 immediately surrounding the smooth cluster CC was excluded.The lower limit of image magnification, |µ| > 100, was chosen as a compromise: for lower |µ|, images are not sufficiently magnified, while for higher |µ|, the number To isolate the effect of subhalos, the band of width ±0.2 around the smooth cluster CC was excluded.Each point represents one realization, and there are between 8 and 16 realizations for most of the cases.If fewer points of a given type appear on the figure, that means either the subhalo length or the number of images with |µ| > 100 was zero, so these runs were omitted from the plot.We use 2 values for core radii, 22, 44 and 88 pc, and 3 value for subhalo mass fraction: ≈ 0.21%, ≈ 0.42%, and ≈ 0.9%, as indicated in the legend.In addition to the cases we present in the rest of the paper, here we also show cases with subhalo core radius = 44 pc, and subhalo mass fraction of 0.25%.The points roughly follow a diagonal.
of images drops rapidly (see Figure 6), resulting in poor statistics.
Figure 7 shows that the length of subhalos' CC is shorter on the -ve cluster side: for any given parameter set, triangles are mostly to the left of squares.The number of magnified images tracks the CC length well, with approximately the same linear relation on the +ve and -ve sides.(We construct a similar plot using |µ| > 300 images, which had typically, about 10 times fewer images, but the linear trend was still evident.) The length of the subhalo CC, which is closely related to the size of the caustic in the source plane, is the main (though probably not the only) factor that determines the number of highly magnified images.The factors that determine the length of subhalos' CC are discussed in Section 4. Because the CC length is larger on the +ve side, very bright images are found preferentially on the +ve side of the cluster.

General considerations
Since we are considering image properties near the CC, it is important to discuss the asymmetry between the +ve and -ve cluster sides, which affects image properties.This is addressed analytically in Venumadhav et al. (2017), Oguri et al. (2018), Diego et al. (2018) for microlensing, but these papers do not consider the situation where the macro density due to the cluster is different on the two sides of the cluster CC, and the effect of the subhalo core size.Here, we investigate these properties.Dai et al. (2020) and Kaurov et al. (2019) do simulations of cored subhalos, but do not explore the effect of these properties.
In the following subsections we discuss how image properties are affected by the cluster density gradient, the location of subhalo on the +ve and -ve sides of the cluster, and the core sizes of subhalos.The combined effect of these factors determines the length of subhalos' CC, and hence the abundance of highly magnified images, as presented in Section 3.3.We consider a simple example: two identical subhalos placed at equal distances on the opposite sides of the cluster CC.The resulting critical curves are shown in Figure 8.Though the subhalos are identical, the critical curves associated with them are different in two aspects: on the -ve side, the subhalo CC is smaller and forms an upright figure '8'.On the +ve side, the CC is larger and forms an '∞' sign.In the rest of this section we will describe how these and other properties arise.
Let the macro cluster potential be a singular isothermal sphere, SIS.Whether the cluster has a flat density core or not does not matter, since cluster CC is well outside it, where the density profile is falling as a power law.Its lensing potential is, where we ignore the normalization as it is not important for us.The corresponding dimensionless projected density and shear components are given by κ = 1 2r The total shear amplitude due to SIS is given by γ = (γ 2 1 + γ 2 2 ).SIS has κ = γ at all radii, and both the projected density and shear fall as 1/r.The two eigenvalues The critical curve is at radius r = 1, and the corresponding density is κ CC = 0.5.The shape of the iso-γ 1 and iso-γ 2 contours are shown in Figure 9, where (0, 0) corresponds to the cluster center.
The orientation of our Cartesian coordinate system is shown in Figure 1.The cluster CC stretches the image almost exactly vertically, along the y−axis.This makes the analysis in this Section easier, but the conclusions obtained here hold in general.

Why the subhalo CC's look like '8' and '∞'
The gray lines in the left panel of Figure 8 are tangential critical curves (see also Chang & Refsdal 1984;Oguri et al. 2018;Diego et al. 2018).In the absence of subhalos the line would have been an uninterrupted smooth approximate diagonal.The two subhalos create additional loops.The condition for the tangential critical curve in the presence of subhalos can be written as Here, subscripts m and s stand for the main cluster subhalo, respectively.
Figure 9 shows the pattern of γ's for the SIS lens.The center is at the center of the cluster, and our modeling window is on the positive x-axis, where y = 0 (see Figure 1), which means that γ 2m ≈ 0. Figure 9 can also be used to approximately estimate the γ values of subhalos.In that case the (0, 0) in that figure would correspond to the center of a subhalo.The subhalo density profile is not SIS, but the general four-leaf pattern with alternating positive and negative values is the same.That means that for subhalos γ 2s ≈ 0.
The density profile of any realistic subhalo in a cluster will not be a SIS.If its profile shape is similar to those from simulations, its central density slope will be somewhat shallower than r 3D −2 (or, r 2D −1 ), and considerably steeper than that at the outskirts, due to tidal truncation.For such a profile, κ s > γ s in the central region (in the limiting case of a very large flat density core, κ γ), but further out, κ s < γ s .(As a limiting case of centrally concentrated lens consider a point-mass lens where κ = 0 and so is always smaller than γ, outside of the actual point.)Let us consider such a centrally Figure 9. Contours of equal shear (γ = 0.1) of a SIS lens, eq.2; left: γ1; right: γ2.Positive and negative values are shown as solid and dashed lines, respectively.Location (0,0) is the lens center.Note that even for a non-SIS density profile, like that of a subhalo, this four-leaf structure with alternating positive and negative values is approximately the same.We use this structure in Section 4 to illustrate the values of γ1m, γ1s, and to show that γ2m ≈ 0, and γ2s ≈ 0, since our modeling frame is along the positive x-axis, where y = 0. concentrated subhalo.At the location of subhalo CC, which is some distance away from its center, its surface density is small compared to its shear, and we can assume that κ s ≈ 0. Tangential CC condition of eq. 3 is now reduced to Next, we use Figure 9 to estimate γ's for the SIS cluster lens.The center is at the center of the cluster, and our modeling window is on the positive x-axis (y = 0), and is small compared to the extent of the cluster, which means that γ 1m is negative everywhere in our modeling frame, and since we already established that γ 2m ≈ 0, we conclude that γ 1m ≈ −κ m .
The only unknown in eq. 4 is γ 1s .Because the other two terms in that equation refer to the cluster, their values will stay approximately constant near the subhalo CC.Therefore the shape of the subhalo CC will follow the shape of γ 1s ≈ const, i.e., one of the shapes in Figure 9.To determine which one, we need to figure out the sign of γ 1s for the two subhalos on the two sides of the cluster CC.
Let us first consider the case of γ 1s > 0. Equation 4 can be rewritten as 0 = 1 − 2κ m + γ 1s .This can be satisfied only if κ m > κ CC , where κ CC is at the location of the cluster CC, and is 0.5 for a SIS.Cluster density is larger than κ CC on the -ve side.The shape of γ 1s when its value is positive is an upright '8'.Therefore subhalo CC's will have that shape on the -ve cluster side.This is what we see in Figure 8. Now let us consider the case of γ 1s < 0. Equation 4 can be rewritten as 0 = 1 − 2κ m − |γ 1s |.This can be satisfied only if κ m < κ CC .Cluster density is lower than κ CC on the +ve side.The shape of γ 1s when its value is negative is an '∞' sign.Therefore subhalo CC's will have that shape on the +ve cluster side.This is what we see in Figure 8.

The effect of cluster density profile gradient
The density profile of the cluster plays an important role in the size of subhalo CC.Let us assume that the cluster density profile is a power law κ(r) = κ CC r −α , where r is the distance away from the cluster center in units of its CC distance, so at CC, r = 1, and the projected density at that location is κ CC .The density at a distance ∆r away from the cluster CC, i.e., at a distance r = 1+∆r from the cluster center can be obtained using Taylor expansion: where ∆r can be positive or negative.For a SIS which has α = 1, this becomes At very small distances from the CC, the [∆r] 2 term is negligible, and κ can be represented by a constant density slope.But for somewhat larger displacements from the CC, those relevant for subhalos, the [∆r] 2 term is important, and it makes the slope nonlinear.
On the +ve side, ∆r > 0, while on the -ve side, ∆r < 0, implying that |∆κ − | > |∆κ + | for the same value of |∆r| on both sides of the cluster CC.In other words, κ rises faster away from the CC on the -ve side, compared to its fall on the +ve side, consistent with the power law in the linear-linear representation.With this we can rewrite eq. 4 for the +ve side as where the last expression assumed that κ CC = 0.5.(We use −|γ 1s | to explicitly state that that term is negative.)Or, equivalently, For the -ve side, the corresponding equation is (10) In eq. 10 the γ 1s term has to compensate for two negative terms to satisfy the equation, while in eq. 9, the |γ 1s | has to compensate for a smaller value.Therefore, γ 1s,eq.10> |γ 1s,eq.9|.Larger values of shear are closer to the center of the subhalo, implying that on the -ve cluster side, subhalos' CC will be smaller in size.Note that this size asymmetry is due to the presence of the [∆r] 2 term, which has different signs in the two equations.In a cluster with a linear density profile, instead of a power law, this asymmetry in subhalo CC size would not exist.
The asymmetry is seen in the left panel of Figure 8.It is also present in the subhalos' signature on the +ve vs. -ve parity sides of the cluster in figure S4 (p.48) of Meneghetti et al. (2020), which shows caustics on both sides of the cluster CC.
It follows that steeper power law cluster density profiles will result in larger disparity between subhalo CC sizes on the two sides of cluster CC.On the other hand, if the cluster has a linear density slope, of the density is approximately the same on both sides, then subhalo CC sizes will be approximately same, which is seen in figure 4 (right panel) of Diego et al. (2018).
To quantify the effect of shallower power law cluster density profiles near the locations of CC, we also do our simulations for ∼ r 2D −0.9 profile, and present these results in Appendix B.

The effect of subhalo core size
The case of large subhalo core on the +ve vs. -ve sides can be addressed with a simple analytical argument based on the discussion earlier in this Section.In section 4.3 we considered a subhalo density profile with a small concentrated core, and saw that there existed solutions to eq. 9 and 10 for subhalo CC on either side of the cluster CC.In the case of a subhalo with a large flat density core, γ κ, and there may not be a solution to eq. 10, implying that subhalos may not even form their own CCs on the -ve side.On the +ve side, on the other hand, it would be easier to have solutions to eq. 9, for the same subhalo properties.
Absence of subhalo CC for large core radii, and on the -ve cluster side, means that there will not be any highly magnified images there.That is why there are no orange or green triangles (r c = 44 pc and 88 pc; -ve side) in Figure 7.
Even on the +ve side, the subhalo core size will affect the extent of its CC.To demonstrate that we perform a numerical experiment.We generate two new runs, each with one subhalo of the same mass, but different core radii, on the +ve cluster side, and at the same distance from the cluster CC.The same number of background point sources are used in both cases.Figure 10 shows the results: the subhalo in the top and bottom panels have core radii of 22 pc, and 88 pc, respectively.Curves of equal magnification are also shown, as well as the images of +ve (orange) and -ve (magenta) parity.The larger core radius results in subhalo's mass getting more spread out, and its CC and high magnification regions smaller, leading to fewer highly magnified images compared to a smaller core radius case.This is why in Figure 7 the blue squares (r c = 22 pc; +ve side) are to the upper right of the orange squares (r c = 44 pc; +ve side), which are in turn to the upper right of the green squares (r c = 88 pc; +ve side).

The effect of the shape of the arrival time surface
In the left panel of Figure 11 we show a smooth substructure-less cluster that generates images of a point source, of opposite parity on either side of its CC.The unsigned image magnifications are approximately the same.Black contours are those of the arrival time surface of two sources, whose images are the locations marked with red (+ve parity image), and blue (-ve parity image) circles.Since this is a subhalo-less cluster, +ve (-ve) parity image is on the +ve (-ve) side of the cluster.The unsigned magnifications of these are about equal, as expected, 29.86 and 27.87.
In the presence of substructures, either subhalos or microlenses that are concentrated enough to form two extra images, the +ve side will add a lemniscate contour in the time delay surface with a minimum and a saddle, or a lemacon contour with a maximum and a saddle (see figure 1 in Saha & Williams 2011).The -ve side will add a maximum and a saddle (see upper left of figure 2 in Saha & Williams 2011).Because minima are usually not demagnified below the magnification due to the smooth cluster, while saddles and maxima can and often are, the +ve side of the cluster will generally acquire more magnified images than the -ve side.
In the case of subhalos, which differ from microlenses by having an extended mass distribution and possibly a flat density core, an additional effect is in play, which makes minima on the +ve side even more magnified.The right panels of Figure 11 use the same sources, and same patches of the lens plane, but now introduce one subhalo at each of the locations marked by yellow crosses.The mass of each is 10 7 M , and the core radius = 88 pc.On the +ve side of the cluster a subhalo forms a lemniscate configuration with 3 images, as stated above.Since the contour levels are drawn at the same value for all 4 panels, their spacing can be used to judge the relative magnification of images.The shallowness of the arrival time surface on the +ve side indicates that the images will appear bright; in fact the total magnification of the 3 images is 93.01, compared to µ = 29.86 on the +ve cluster side image in the left panel of Figure 11.The shallowness of the arrival time can be explained as follows.
On the +ve side of the subhalo-less cluster the arrival time surface around an image looks like a bowl, and so hosts minima.Adding a nearby subhalo makes the profile of that bowl shallower (with or without creating extra images), resulting in higher magnification.
The situation on the -ve parity side of the cluster is different: the shape of the arrival time is that of a saddle, and so hosts saddle images.Adding a smooth bump due to a subhalo does not flatten out the arrival time shape, and hence does not increase the image magnification.This is seen in the spacing of the contours on the -ve side of the two panels.In fact, the image magnification in the presence of a subhalo (right side of the right panel) is lower (µ = −7.48)than that of the same image when subhalo is absent (µ = −27.87).
We conclude that on the +ve side of the cluster, subhalos, which have an extended mass distribution, and possibly a flat density core, will generally magnify images more than point mass microlenses would alone.The exception is when microlenses form more than 2 extra images on the -ve side, some of which will be minima (see the other 3 panels of figure 2 of (see the other 3 panels of figure 2 in Saha & Williams 2011).

OBSERVABLE SIGNATURES OF SUBHALOS
Our goal in this section is to identify tools to detect the presence of dark matter subhalos near cluster CC, and to estimate their properties.We use subhalos in the mass range 10 6 − 10 8 M only; more massive subhalos will likely produce a stronger signature.We describe two different types of subhalo signatures: their effect on the pixelated flux distribution (Section 5.1), and the effect of individual highly magnified images (Section 5.2).In this paper we work with simulations only, but plan to extend the method to observations of clusters, specifically from the ongoing HST Flashlights project.

Pixel-level signatures
First, we devise a metric suitable for working with flux collected in detector pixels, with an applied PSF, which we model as a Gaussian with a width of σ = 0.04 .(A similar problem has recently been analysed by Bayer et al. (2023) for subhalos in galaxy-scale lenses.)We bin the flux from images into pixels, and then apply PSF smearing of width 0.04 , comparable to that of optical HST PSF's.2

Local flux peaks and holes
In Figure 3 we saw that away from the cluster CC, highly magnified images are strongly associated with subhalos' CC, on the +ve cluster side.The amplitude of association depends on the subhalos' density profile and the cluster density profile.Because lensing magnifies point sources and at the same time dilutes their number density, the lens plane number density of highly magnified images is decreased, by the same factor as the magnification.Therefore, near subhalos on the +ve side, which generate a lot of highly magnified images with a low number density, we also expect to find considerable fluctuations in pixel fluxes between different lens plane locations.
To quantify that, we first introduce a metric to quantify local fluctuations, and then describe how we measure their distribution in the lens plane.Flux variations, i.e., the heights of local peaks and the depths of local holes, are calculated as, where f in is the flux in a 0.04 × 0.04 pixel, and f out is the average flux in the surrounding pixels, up to p pixels away.For p = 1 (p = 3), the number of surrounding pixels is, (2p + 1) 2 − 1, i.e., 8 (48).After experimenting with a few values of p we settled on p = 3.A value of h = 0.02 corresponds to a flux ratio (f in /f out ) of about 1.05 of the central pixel compared to the average of its neighbors, or −0.05 magnitudes.Peaks have h > 0, while holes have h < 0. We measure the distribution of peaks and holes as follows.At a distance d from every subhalo we calculate the average value of |h| , weighted by the mass of the subhalo, The sum is over all pixels, and all subhalos, but leaving out subhalos in a band around CC of ±0.1 .We then normalize these by the corresponding values of subhaloless clusters, to get a correlation function-like quantity, In subhalo-less clusters, instead of subhalo positions we use randomly generated points in the lens plane.In Figure 12 we plot eq. 13.As the previous plots, solid lines represent the +ve side of the cluster.These show a positive correlation with subhalos, which means that subhalos are surrounded by local flux peaks and holes of larger amplitude than peaks and holes further away from subhalos.Dashed lines represent the -ve side of the  |h| (eq. 11) in the presence of subhalos, compared to the case where subhalos are absent, as a function of the distance between it and subhalos (eq.12 and 13).The average is over all the subhalos outside the ±0.1 band around smooth cluster CC.Sixteen individual realizations were used for each line.The legend indicates model parameters.There is a correlation between peaks/holes, and subhalos on the +ve cluster side (solid lines), but not on the -ve side (dashed lines).
cluster; locations of peaks and holes on the -ve side of the cluster do not correlate with subhalo locations.
The fact that subhalos are surrounded by regions of higher fluctuations in pixel flux is reminiscent of the principle used to measure distance to ellipticals using surface brightness fluctuation method (Tonry & Schneider 1988).In those studies, galaxies are assumed to have the same surface brightness, and the reason for some galaxies having fewer, but brighter stars, instead of more numerous fainter stars per sky area is the distance to the galaxies.In the case of lensing, it is the magnification of images coupled with the decrease in their number density on the sky.In both cases the underlying principle is the conservation of surface brightness.
Since the correlations quantify the average number density profile of peaks and holes around stacked subhalos, we do not expect w(d) in Figure 12 to be different for cases of low vs. high subhalo mass fraction, but we might expect a difference between the cases of small and large subhalo core radii.A small difference is visible, with the smaller core subhalos generating somewhat stronger correlations.
The main conclusion from Figure 12 is that for subhalo mass fraction of ∼ 0.5 − 1% (and presumably higher) the amplitude of flux variations between pixels in a radius 0.2 around subhalos is significantly higher than at larger distances.Furthermore, such fluctuation will be found predominantly on the +ve cluster side.That means flux variations may signal the presence of dark subhalos on the +ve side; on the -ve side subhalos do not generate higher than typical flux variations.

Detecting pixels affected by subhalos
The correlation function, w(d) described above, is one possible way to detect subhalos.Another observable test for subhalos would take advantage of the fact that in the presence of subhalos, +ve and -ve sides of the cluster have different degrees of flux variations between pixels.When subhalos are absent, the differences between the two side of the cluster are small.
In the left panel of Figure 13 we plot the fraction of 0.04 pixels in the lens plane that have h value indicated on the horizontal axis.The holes and peaks are on the left and right sides of each panel, respectively.The model parameters are shown in the left panel.A band of ±0.1 surrounding cluster CC was omitted to isolate the effect of subhalos from that of the cluster CC.The solid and dashed lines represent the +ve and -ve sides of the cluster CC, respectively.The gray line shows the subhalo-less case.At low |h| values all distributions are very similar.At |h| 0.03 − 0.04, the solid blue and green lines are above the gray lines, meaning there are more pixel-to-pixel variations in flux in the presence of subhalos.More importantly, the solid lines tend to be above the corresponding dashed lines, implying that more of the pixel-to-pixel variations in flux are on the +ve cluster side.This is the type of signature that is in principle observable.|h| = 0.04 corresponds to 0.1 magnitudes difference in flux.If all sources have the same intrinsic luminosity, as the left panel assumes, flux peaks and holes are probably only barely measurable given the typical noise level of observations.
The right panel of Figure 13 is the same as the left, except that the sources were drawn from a luminosity function of the form, dn/dL ∝ L −2.5 .Here, the distribution is significantly wider than on the left, implying that larger values of |h| are possible; in fact, for |h| > 0.05 (or, 0.125 magnitudes), there should be very few (fraction < few×10 −4 ) holes and peaks on the -ve cluster side, or in a subhalo-less case.On the +ve side of a cluster with subhalo, the fraction of pixels with these value of |h| show be several times higher.Given these numbers, determining if subhalos are present or absent may be possible using one region of the cluster of size 1few arcseconds, but extracting subhalo parameters, for example, subhalo mass fraction or core sizes will require more than one region in a galaxy cluster, and possibly more than one cluster.
To sum up, a signature of subhalos is that the number of peaks and holes of progressively higher amplitude gets larger on +ve side of the cluster compared to the -ve side.This asymmetry should be looked for at some distance away from the cluster CC, ∼ ±0.1 , so as not to confuse it with that due to microlensing by stars near the cluster CC.
We note that pixel-level detection relies on the fact that point sources, and hence lensed images, have a finite projected number density on the sky.If the source density were very high, or if the source were a truly continuous, though not necessarily uniform surface brightness distribution, these tests would not work.They rely on discreteness and stochasticity of image distribution.

Individual highly magnified images
When subhalos magnify images by a factor of a few 100 to 1000, their flux will contribute to the total flux of the pixel they are in, as discussed above.But if an image is magnified by more than a few 1000 it can be detected on its own.This has been the case with in-dividual star detections, (e.g., Kelly et al. 2018).Very high magnifications are possible by subhalos alone, or with the help of microlensing.

Magnification by subhalos alone
Figure 6 shows that the high magnification tail of image distribution in the lens plane scales as |µ| −2 .Extending that relation to higher magnifications suggests that fewer than 1 in 10 6 images will attain |µ| > 10 4 , and be detectable as brightened images of individual stars.Since the number of images per ∼ 1 − 2 sized region is of order of a million, there is a reasonable chance of detecting subhalo-brightened images of individual stars.If these happen to lie close to the cluster CC, it will be hard to ascertain that they were brightened by a subhalo, since microlensing is more likely to be responsible.Subhalos are less likely to produce rapid time variability, compared to microlensing.

Magnification by subhalos and microlensing
If an image that was moderately magnified by a subhalo also gets microlensed, then the combined magnification can be high enough to make it detectable as an individual highly magnified star.
To estimate the statistical effect of including microlensing by stars and other point lenses in the clus-ter, one would convolve the image luminosity function with the microlensing probability distribution.Convolving steep image flux distribution with shallower magnification probability distributions, like that of microlensing, will result in a shallower observed luminosity function.Consulting Figure 6, which plots image luminosity function, we see that microlensing will have the largest impact on the images on the -ve cluster side and subhalos with larger core radii, shown as the thick dotted green line, because their image flux distribution has a sharp drop off at intermediate |µ|'s.Therefore convolving with microlensing PDF will extend these distributions to higher magnifications.This will increase the number of highly magnified images on the -ve side, but their number should still remain small compared to the ones the +ve side.
The PDF of image magnifications from microlensing differ on the +ve and -ve sides, as is shown in figure 14 of Diego et al. (2018): on the -ve side, where image demagnifications are usually possible, the PDF has a flat plateau for |µ| up to ∼ 300, and a high magnification tail scaling as |µ| −3 .

Where to look for subhalos
The calculations in this paper were tailored to a cluster of mass similar to that of SGAS J1226+2152, which was studied by Dai et al. (2020) who found variations in flux on either side of the cluster CC.Other cluster candidates include Hamilton's object (Griffiths et al. 2021), and Lyman-α detection in RXJ 1347.5-1145(Richard et al. 2021).
In general, any images straddling the cluster CC are good targets to look for subhalo lensing, especially not in the immediate vicinity of the cluster CC, because as we saw in Figures 3 and 4, subhalo magnified images can appear further away from the cluster CC, where they are less likely to be confused with microlensingbrightened events.Microlensed stars are also more likely to occur in close pairs such as in Kelly et al. (2018); Chen et al. (2019), and are more likely to give rise to time variability, compared to subhalos.

Subhalos vs. non-standard forms of dark matter
While dark matter is generally assumed to consist of massive non-relativistic particles, it is possible that its nature is very different from that.Wave, or fuzzy dark matter, composed of ultralight particles, can introduce both positive and negative mass fluctuations about the mean density throughout the cluster as a result of quantum interference (Schive et al. 2014;Hui et al. 2017;Laroche et al. 2022;Amruth et al. 2023).These pervasive fluctuations have been shown to produce significant perturbations to the brightness and positions of lensed images in galaxy scale lenses (Amruth et al. 2023), while also being able to reproduce the observations of specific lensing systems that smooth lens models typically struggle to explain.These results highlight the need for some form of mass substructure in the lens.However, the effect of wave dark matter is expected to be smaller on cluster compared to galaxy scales and remains to be explored in detail, with preliminary results presented in Kelly et al. (2022).CDM subhalos, on the other hand, lead to an asymmetry between +ve and -ve sides of clusters.The presence or absence of asymmetry may be a test for the nature of dark matter.

Clusters with shallower density profiles
So far we worked with the cluster lens whose 3D density profile is that of a singular isothermal sphere.However, shallower profiles are also possible at radii where cluster CCs are formed.In Appendix B we consider a situation very similar to the one above, but with projected density ∝ r −0.9 .Shallower cluster profiles produce notable differences: highly magnified images extend further away from the cluster CC (Figure 15), and should therefore be looked for further from the cluster CC.The flux peaks and holes are stronger correlated with subhalos (Figure 16, right panel), and the amplitude of flux peaks and holes, |h|, are larger (Figure 17), and will therefore be more easily detected than in clusters with steeper density profiles.

SUMMARY AND CONCLUSIONS
In this paper we investigated the effect of dark matter subhalos on lensing near the tangential critical curves of galaxy clusters.We started by first calculating the properties-distribution and magnification-of lensed images highly magnified by subhalos.We then used analytical and numerical experiments to examine and explain the various effects of subhalo lensing.In the last portion of the paper we explored if subhalos have detectable signatures.We divide these into two categories: pixel-level signatures, where subhalo magnifications contribute toward variations in flux between neighboring detector pixels, and images that are so highly magnified that they can be detected on their own as individual stars in source galaxies.In the latter case high magnification can be due to subhalos alone, or aided by microlensing.Our work is preliminary; the actual analysis will need to take into account noise and its properties, which can be non-trivial, arising not only from noise in empty fields, but from intracluster light as well.
Just like microlensing by intracluster stars, subhalos can magnify stars in background galaxies by factors of 100's to 1000's.In both cases, high magnifications happen preferentially on the positive parity side of the cluster.However, there are some important differences between subhalo lensing and microlensing: • Because subhalos are generally significantly more massive than microlenses means that their influence extends further away from the cluster critical curve (Figures 3,4 and 15).If magnified images are found far from the cluster CC, they are more likely to have been magnified by subhalos rather than microlensing (Meena et al. 2022b).
• If the cluster density profile is a power law (instead of linear), subhalos of a given mass will produce longer critical curves around themselves on the positive parity side, compared to the negative side.This can lead to a significant asymmetry on the two sides of the cluster: there will be more highly magnified images on the positive parity side (Figures 3,4 and 15).The asymmetry is stronger for steeper cluster density slopes, and is more pronounced for subhalos than microlenses because the former's effect extends further away from the cluster critical curve (Section 4.3).
• The number of highly magnified images is proportional to the length of the subhalo critical curves, with approximately the same linear relation applicable to both positive and negative parity sides of the cluster (Figures 7 and 16, left panel).
• The fact that subhalos have an extended mass distribution, possibly with a flat density core (Section 4.4) increases the asymmetry between the image properties on the two parity sides of the cluster, making it more likely that bright images will appear on the positive side (Section 4).Because of their density cores, the asymmetry is stronger compared to the case with microlenses, and can help differentiate between microlens-and subhalomagnified images as well as between small and large subhalo core radii.
• Subhalos with larger core radii generate shorter critical curves around themselves (Figure 10), resulting in fewer highly magnified images.For core radii 50 pc subhalos generate no critical curves at all on the negative cluster parity side.
• Another consequence of subhalo mass being significantly larger than that of microlenses is that their lensing effects are not transient or time variable, unless subhalo lensing is augmented by stellar microlensing, or if the subhalo is low mass, 10 6 M , and the velocity vector of the subhalo relative to the source star is well aligned with the the image stretching due to the main cluster.
LLRW, PLK, and TT would like to acknowledge HST programs GO-15936 and GO-16278, and SNAP program GO-16729  One can show that ∇ψ(r) = α(r), and that the second derivatives of ψ can be combined through the lensing Poisson equation to give κ(r): 1 2 (ψ ,xx + ψ ,yy ) = κ.The 2D density profiles for a range of a values are shown in Figure 14.A range of shapes can be obtained by varying a. Smaller a result in smaller flattish density cores and gentler transition to steeper outer profiles.On the other hand, larger a have extended flat density cores, and a sharper transition to outer steep slope.The outer profile always goes as r −3 in projection, ensuring that mass converges quickly.
For comparison, we plot NFW profile in brown.It has been rescaled vertically and in radius to show the close resemblance to one of our functions: 3rd from the top, with a = 0.05.The small vertical arrow at the bottom of the plot shows the NFW's scale radius r s .NFW declines slower at larger radii, as r −2 in projection.

B. SUBHALOS IN A CLUSTER WITH A SHALLOWER DENSITY PROFILE
In the main paper we considered a cluster whose density profile slope is that of an isothermal sphere, i.e., ∝ r 2D −1 in projection.One of our findings is that the cluster density slope is an important factor in determining the properties of highly magnified images formed near its critical curve.For that reason, in this Appendix we consider a somewhat shallower projected density profile, that is also consistent with observations of clusters, ∝ r 2D −0.9 .The plots presented here (Figures 15, left and right panels, 16, left and right panels, and Figure 17) are similar to the corresponding ones in the main paper (Figures 3,4,12,and 13).
The conversions relating linear and angular scales in the source and lens plane are the same as in the main paper, but the number density of points sources in the source plane is 0.03 stars/pc 2 here.As in the main paper, we are concerned only with intrinsically luminous stars, which are rate, motivating the low number density.
The distribution of images around the cluster CC here is wider than for the steeper cluster profile; in fact, the right panel of Figure 15 shows that the distribution does not peak at the cluster CC; the images are uniformly distributed within ∼ ±0.15 .This flattening of the distribution is the direct result of the shallower cluster density profile, which makes the region where (1 − κ − γ) ≈ 0 into a wider band around the CC.It is similar to the case of microlenses that reach optical depth ∼ 1 near cluster CC.The same effect is also seen in models with wave dark matter where the key parameter is not cluster slope or microlensing optical depth, but the mass of the boson.The asymmetry seen in the steeper cluster density case, where bright images were preferentially found on the +ve cluster side, is also present here, but is seen further away from the cluster CC, at 0.15 .The spatial correlation of local flux peaks and holes with subhalo locations is stronger here (right panel of Figure 16) compared to the case of steeper cluster profile.That means that flux peaks and holes are more likely to indicate the locations of subhalos.The amplitude of local flux peaks and holes (i.e., values of |h|) are also larger here (Figure 17).

Figure 1 .
Figure 1.Cluster-scale critical curve.The source is the blue star symbol and its images are shown with red circles.The merging pair of images straddling the CC that we will model in detail is on the right side of the cluster.The box around it is the size of the modeling window we use, for example in Figure2.The lower left corner of the box defines the origin of our coordinates.Our Cartesian coordinate system is oriented such that the image stretching is very nearly along the vertical y−axis, which aides in the analysis in Section 4.

Figure 3 .
Figure 3.A 2 box in the lens plane, showing a piece of cluster CC in the presence of subhalos with mass fraction ≈ 0.5%.The figures contain a superposition of 16 different realizations, to enhance the visualization of the effect of subhalos.Subhalos themselves are not marked to avoid overcrowding.Critical curves are shown as gray.The +ve and -ve parity sides of the cluster are indicated with '+' and '-' brown signs.Left: Subhalo core radii are rc = 22 pc.The projected density in the cluster is ρ ∝ r2D −1 .number of subhalos per realization is ∼ 43.Images with +ve and -ve parity with |µ| > 1000 are red and blue, respectively.Right: Subhalo core radii are rc = 44 pc.The projected density in the cluster is ρ ∝ r2D −0.9 .The number of subhalos per realization is ∼ 48.Images with +ve and -ve parity with |µ| > 100 are red and blue, respectively.Note that aside from the band immediately around the cluster CC, highly magnified images in both panels are found preferentially on the +ve parity side.See Figure4for the histogram of image distribution as a function of distance from the CC for the case presented in the left panel.

Figure 5 .
Figure5.Source magnification distribution.Source magnification is the sum of unsigned image magnifications that end up in our 2 lens plane modeling window.The extra saddles created by subhalos are mostly responsible for the low magnification side of the distribution.The two green (blue) histograms represent subhalos with core radii rc = 88 pc (rc = 22 pc).The thin (thick) lines show subhalo mass fraction of 0.21% (0.42%).The gray line (mostly hidden behind other lines) depicts the case with no subhalos.

Figure 6 .
Figure6.Distribution of unsigned image magnifications in the image plane.The color and line types are indicated in the plot.To avoid the high magnification images associated with the cluster CC itself, we left out a band of width ±0.1 around the cluster CC.Left panel: Positive parity cluster side.Right panel: Negative parity cluster side.Subhalos with larger core radii (green histograms) do not produce any very bright images, so the gray (no subhalos) and green histograms nearly overlap.

Figure 7 .
Figure7.Log of the number of magnified images with |µ| > 100 vs. log of the subhalo CC length (in code units).To isolate the effect of subhalos, the band of width ±0.2 around the smooth cluster CC was excluded.Each point represents one realization, and there are between 8 and 16 realizations for most of the cases.If fewer points of a given type appear on the figure, that means either the subhalo length or the number of images with |µ| > 100 was zero, so these runs were omitted from the plot.We use 2 values for core radii, 22, 44 and 88 pc, and 3 value for subhalo mass fraction: ≈ 0.21%, ≈ 0.42%, and ≈ 0.9%, as indicated in the legend.In addition to the cases we present in the rest of the paper, here we also show cases with subhalo core radius = 44 pc, and subhalo mass fraction of 0.25%.The points roughly follow a diagonal.

Figure 8 .
Figure 8.An illustration of subhalos' CC on the +ve and -ve sides of the cluster.(These runs were not used for the actual calculations.)Left panel: Image, or lens plane with two identical subhalos equidistant from the cluster CC.As discussed in Section 4 lensing rotates the orientation of the subhalos by ∼ π/2 on the opposite sides of the cluster CC, and size is smaller on the -ve side due to density gradient of the main cluster.Right panel: The portion of the source plane from which the images in the left panel would originate.The 6 different colors show regions giving rise to 6 different image multiplicities.Note that only those images that would appear in the image plane square of the left panel are counted.This figure is for illustration purposes only; these subhalos are more massive than those used in the main calculations.Subhalos used in calculations rarely generate 5 or 6 image configurations.

Figure 10 .
Figure 10.Zoom-in regions around two individual isolated subhalos, in two separate controlled runs.Both are on the +ve side of cluster CC.Both subhalos have the same total mass, but different core radii.Top: Subhalo core radius = 22 pc; Bottom: Subhalo core radius = 88 pc.The critical curves are the gray dotted lines; images of positive and negative parity are shown as orange and magenta points respectively.The density of source in the source plane is the same in both cases.The 3 sets of black curves show unsigned magnifications of 300 (thick), 100, and 30 (thin).The dilution of images in regions of high magnification is apparent in both cases.

Figure 11 .
Figure 11.Arrival time surface around the locations of individual subhalos, indicated by yellow crosses.Positive (negative) parity images are red and blue circles, respectively.The left and right sides within the left and right panels show the +ve and -ve sides of the cluster, as marked on the plots.Left panels: Subhalo mass set = 0, to show the location of images for a subhalo-less cluster.The magnifications of the two images are 29.86 (red circle), and -27.87 (blue circle).Right panels: Each subhalo has mass of 10 7 M .The magnification of the -ve parity image on the -ve side (blue circle) is -7.48, while the total unsigned magnification of the 3 images of the lemniscate on the +ve side is 93.01.To facilitate comparison, contour levels on all 4 panels are the same values, and are spaced by about an hour in observer frame.

Figure 12 .
Figure 12.Average value of the pixels' local flux peak height or hole depth,|h| (eq.11)  in the presence of subhalos, compared to the case where subhalos are absent, as a function of the distance between it and subhalos (eq.12 and 13).The average is over all the subhalos outside the ±0.1 band around smooth cluster CC.Sixteen individual realizations were used for each line.The legend indicates model parameters.There is a correlation between peaks/holes, and subhalos on the +ve cluster side (solid lines), but not on the -ve side (dashed lines).

Figure 13 .
Figure13.Observable signature of subhalos: the fraction of pixels with values of h, eq.11.Positive and negative values of h are for local peaks and holes, respectively.Solid and dashed distributions represent +ve and -ve sides of cluster CC, respectively.A band of ±0.1 around the cluster CC was omitted.See the legend for model parameters.Left panel: assuming that all sources have the same luminosity; Right panel sources have a luminosity function dn/dL ∝ L −2.5 .Each line is an average over 16 realizations of the 2 modeling window.
. PLK acknowledges NSF AST-1908823 grant.JMD acknowledges the support of projects PGC2018-101814-B-100 and MDM-2017-0765.AKM and AZ acknowledge support by Grant No. 2020750 from the United States-Israel Binational Science Foundation (BSF) and Grant No. 2109066 from the United States National Science Foundation (NSF), and by the Ministry of Science & Technology, Israel.

Figure 14 .
Figure14.The function we use to represent subhalos is shown as blue curves.The shape parameter a determines the inner profile shape.In this paper we use a = 0.05, which closely resembles NFW (brown curve), except for large radii where our function falls off quicker, approximating tidal truncation.The 2D slope of d ln ρ2D/dr = −1, which corresponds to an "isothermal" slope in 3D, is shown as a dashed gray line.

Figure 15 .
Figure 15.Left: Similar to Figure 3, but for a shallower cluster density slope.Right: Similar to Figure 4, but for a shallower cluster density slope.Note that the horizontal axis is wider than in the former figure.