Dark Matter under the Microscope: Constraining Compact Dark Matter with Caustic Crossing Events

(Kelly et al. 2018, hereafter K18) present the first observations of a single high-redshift star in a background, lensed spiral galaxy at redshift z = 1.49 (Smith et al. 2009; Zitrin & Broadhurst 2009) being magnified by a factor of several thousand by a galaxy cluster MACSJ1149.5+2223 (hereafter MACS1149) at z = 0.544 (Ebeling et al. 2007).

This event was discovered serendipitously while monitoring a lensed supernova (SN) behind the cluster (SN Refsdal; Kelly et al. 2015, 2016; Rodney et al. 2016). The light curve of the star shows at least one prominent peak in the spring of 2016 that lasted ∼2 months. A first event, named Icarus or LS1/Lev 2016A by K18, produced a peak in the light curve that lasted several weeks; after the peak, the flux returned to its original value. This event is interpreted as a crossing of a bright background star through a microcaustic produced by one of the stars (or star remnant) in the intracluster medium. At a position separated by 026 from this initial peak, a second peak (named Iapyx, or LS1/Lev 2016B by K18) appeared between 1 and 2 months after the first event faded, and lasted less than 3 months. No object was observed at this second position in the previous 10 years or in the months after it vanished. This second event is also interpreted as a microlensing event of the same background star (and a different microlens in the intracluster medium). However, in this case one possible interpretation is that the low-magnification region around a microlens was hiding the background star for >10 years with occasional brief periods of high magnification (see K18 for other possible interpretations, including binary stars). A potential third event discussed by K18, Perdix or Ls1/Lev 2017A, is found 01 away from the second event. If confirmed, this third event could be produced by the same network of microcaustics, but possibly involving a different background star.

As discussed by K18, the picture described above is consistent with the expected behavior of a background star traveling at typical relative velocities of ∼1000 km s⁻¹ and a lens plane populated with a density of stars that is compatible with the observed intracluster light (ICL) at the position of the two events. Microlensing events are expected to be produced by the stars responsible for the ICL (and their remnants). As shown in earlier work, the light curve of an object being lensed by a field of microlenses may contain high- and low-magnification periods. This behavior has been predicted in previous papers (see, for instance, Chang & Refsdal 1984; Kayser et al. 1986; Paczyński 1986). In particular, Chang & Refsdal (1979, 1984) were the first (to our knowledge) to recognize that counterimages may have very low fluxes (disappearing below the detection limit of a given instrument) for some periods of time, in agreement with the observed behavior of the Iapyx event (counterimage of the Icarus event of K18).

A massive galaxy cluster acts as a cosmic telescope that enlarges the images of background galaxies. However, near a critical curve (CC), small changes in the deflection field result in large changes in the magnification. These small changes in the deflection field can be produced by small masses in the range of a stellar mass or below. In this situation, as we will show, the galaxy cluster may act also as a cosmic microscope since it effectively enlarges any imperfection in the deflection field near the CC caused by microlenses. Microlensing near a cluster CC has the interesting feature that the individual micro-CCs around the microlens (and corresponding microcaustics) get enlarged by a factor that is larger the closer they are to the main CC (see the discussion of this effect in Section 2). This allows, in principle, probing small-mass microlenses as we approach the cluster CC.

Earlier work has explored the behavior of counterimages during caustic-crossing events in smooth potentials (from galaxies to clusters). Miralda-Escude (1991) considers, as in this work, the case of a single star crossing a caustic from a smooth lens model. He estimates the maximum magnitude of a lensed background star at the time of caustic crossing, as well as the rate of events based on the surface brightness of a background galaxy (this case is also discussed by Chang & Refsdal 1979, 1984 and Schneider & Weiss 1986). The combined effect of overlapping caustics from an ensemble of microlenses embedded in a stronger gravitational field has been also studied in detail (Gott 1981; Young 1981; Chang & Refsdal 1984; Kayser et al. 1986; Paczyński 1986), in particular in the context of quasar (hereafter QSO) microlensing (Chang & Refsdal 1979; Irwin et al. 1989; Witt et al. 1995; Metcalf & Madau 2001). Kayser et al. (1986) and Paczyński (1986) show how a large number of microlenses embedded in a deep potential can redistribute the magnification, producing complex light curves of a background source. For certain configurations (see, e.g., Figures 9 and 10 of Kayser et al. 1986), the magnification splits into compact regions of large and low magnification. As shown in these papers, a source traveling across this field may disappear suddenly when entering one of the low-magnification regions, only to reappear at some time later as a bright source.

This type of behavior resembles the observed flux in the Icarus and Iapyx events. However, when the microlenses are very close to the CC (a fraction of an arcsecond), the magnification pattern exhibits features that have not yet been studied in detail. Paczyński (1986) investigated the general case of high optical depth of microlenses embedded in a galaxy or cluster potential, but he ignored the effect of shear and focuses on areas in the lens plane that are not close to the main CC. Kayser et al. (1986) included the shear term from the large deflector (cluster or galaxy) in their calculations, but again did not study the particular case of short distances to the main CC. As noted by Paczyński (1986), this regime is computationally very expensive (owing to the very large magnifications involved that require the mapping of a small field in the source plane into a very large field in the image plane), and could not be studied in detail in those early papers.

Some authors have focused their attention on the high-magnification regime (see Wambsganss 1990; Schechter & Wambsganss 2002, and references therein) in the context of QSO microlensing, but these high magnifications are still modest (a few tens at most) compared with the more extreme values (several hundred to several thousand) considered in this paper and do not reveal some of the properties of the lensed images that are accentuated with extreme magnification (see Section 3.2). The smaller magnifications found in QSO microlensing are partially due to the larger intrinsic size of the background source. As we will show later, the maximum magnification attained by a background source scales as the inverse of the square root of its radius. For QSOs, the radius is related to the half-light radius of the accretion disk. These disks are typically of the order of 10 light days, when observed in the optical, and about an order of magnitude smaller when observed in X-rays (Chartas et al. 2009; Dai et al. 2010; Jiménez-Vicente et al. 2012) for typical supermassive black holes. This radius is known to scale with the mass of the black hole (Morgan et al. 2010; Jiménez-Vicente et al. 2015). When compared with the radius of a giant luminous star, the accretion disks around QSOs are approximately a factor 10³–10⁴ times larger. Consequently, the maximum magnification attained by a lensed giant star can be up to two orders of magnitude larger than the corresponding one for QSOs. This is an important advantage that comes with the added bonus that the smaller stellar radii translate into shorter-lived events which are easier to monitor (days as opposed to years). Although QSO microlensing is not directly comparable to the work presented in this paper, there are many similarities. Earlier papers focusing on the interpretation of QSO microlensing contain useful insights that are applicable to this work when the magnifications are significantly higher. Schechter & Wambsganss (2002) present interesting similarities with some of the results given here, in particular when discussing the statistics of the magnification around microminima and microsaddle points.

In this paper, we explore for the first time the regime of very short distances to the main CC (or, equivalently, very high magnification), motivated by the observation of the two (or possibly three) intriguing events discussed by K18.²¹

In the case of a galaxy cluster, its larger size translates into a greater magnification of a background object. Also, if a significant fraction of dark matter (DM) is made of compact objects like primordial black holes (PBHs), galaxy clusters are ideal to study microlensing by these objects since it is possible to find CCs (with high optical depth for microlensing) relatively far away from member galaxies and reduce the impact of stars (or remnants) in these galaxies that could produce similar microlensing events. The case of PBHs is interesting since they are (still) a valid candidate for DM (or at least a fraction of it) in some mass regimes (see, for instance, Carr et al. 2010, 2016a; Clesse & García-Bellido 2015). The fraction of DM that can be in the form of PBHs has been constrained for different PBH masses. The possibility that PBHs constitute a sizable fraction of the DM is interesting and has been studied extensively, although PBHs are excluded as the primary component of DM in virtually all mass ranges. Bird et al. (2016) proposed that at around 30 M_⊙ there is still a range of masses that have not been convincingly ruled out (see also Sasaki et al. 2016; Clesse & García-Bellido 2017, for a related result). Interestingly, if a significant fraction of DM is in the form of PBHs with M ≈ 30 M_⊙, events like the collision of two black holes with these masses would be more common, facilitating the interpretation of the first LIGO detection (Abbott et al. 2016). This interpretation, however, is not supported by the second LIGO event with significantly smaller masses. On the other hand, more recently a new LIGO event as well as a LIGO/Virgo event imply detections of massive pairs of BHs (M_BH ≈ 20–30 M_⊙), implying a higher than expected abundance of BHs with M_BH ≈ 30 M_⊙) (Abbott et al. 2017; The LIGO Scientific Collaboration et al. 2017).

Analyses of multiply imaged QSOs have found that the observed microlensing signal is incompatible with the hypothesis that ∼30 M_⊙ PBHs make up most of the DM (see Mediavilla et al. 2017, and references therein). The same work concludes that the fraction of mass in the form of microlenses can still be as high as 20% of the total mass, but with the most likely mass of microlenses being below 1 M_⊙. If confirmed, this key work leaves little room for the hypothesis that PBH with ∼30 M_⊙ can make a significant fraction of the DM (∼10%) unless extended mass functions (instead of the monochromatic or bimodal models considered by Mediavilla et al. 2017) can have a significant impact on the results, or the size of accretion disks around QSOs are an order of magnitude larger than what has been considered so far (the latter point being an important source of uncertainty in this and other work). Moreover, we should note that in Mediavilla et al. (2017), the limit of high optical depth (for microlensing) does not seem to be explored and, as we will show later, in this regime the fluctuations in flux are smaller owing to the constant presence of multiple overlapping microcaustics that tend to average out the observed integrated flux. It would not be surprising to have constraints from QSO microlensing that differ from (or even contradict) those derived from microlensing of background stars (this work). If one finds that tensions between these regimes exist, some of the assumptions made in each regime will have to be reviewed. Constraints from microlensing in our local environment (the Magellanic Clouds) are weaker, and recent work has shown that uncertainties in these constraints can be as high as one order of magnitude (Green 2017).

Carr et al. (2017) review the constraints for PBHs using more realistic extended mass functions and conclude that one could allow for as much as 10% of the DM in the form of PBHs in the mass range M_BH ≈ 25–100 M_⊙ (although this work does not include the results of Mediavilla et al. 2017). This limit of 10% will be adopted in this paper as an upper limit for the fraction of PBHs in this mass range. At smaller masses, constraints on the fraction of PBHs allow for a modest fraction of DM below M ≈ 1 M_⊙ (see, however, Inomata et al. 2017; Kühnel & Freese 2017, for the mass range M_BH ≈ 10⁻¹⁰–10⁻⁸ M_⊙). These constraints tighten at very low masses. A lower limit for the PBHs of M ≈ 10¹¹ kg can already be established from theoretical grounds and observations of γ-rays (see, e.g., Kim et al. 1999). Below this mass, no PBHs are expected to exist as they should have evaporated by now (down to the Planck mass). This limit can be increased a little from detailed observations of the γ-ray background, since sufficient PBHs with masses near the above limit would be a strong source of γ-rays in our vicinity, which is not observed. Continuous monitoring of background galaxies intersecting a cluster CC provides an excellent data set for constraining the abundances of PBHs based on their lensing signature. Combining these data with models of the full stellar population in the lensing plane can address many of the systematic biases inherent in past measurements.

In this paper, we explore a different technique to constrain the fraction of compact DM, paying particular attention to the mass range relevant for the three most significant LIGO events. We show how microlensing events by relatively small masses can take place thousands of years before (or thousands of years after) a bright star in the background galaxy crosses the position of a cluster's main CC. Hence, the probability of observing a microcaustic crossing event is considerably increased when compared with earlier work that only considered the crossing of the main cluster caustic. As mentioned earlier, as the background star approaches the main cluster CC, the sensitivity to detect progressively smaller microlenses grows, offering a unique opportunity to probe masses that could not be tested otherwise. This provides an exciting opportunity to set limits on the fraction of DM in the form of compact objects in low-mass regimes that are difficult to study otherwise.

This paper is organized as follows. In Section 2 we describe the basic properties of the magnification near a CC. Section 3 presents results based on numerical simulations with a focus on the structure of caustics in the source plane. In Section 4 we explore in detail the disruption of the CC in the image plane when microlenses populate the lens plane. We predict in Section 5 the behavior of the observed flux (light curve) of a background star traveling through a field of microcaustics. In Section 6, we predict how events like Icarus will disappear (or first appear) once the last (or first) microcaustic is crossed. Section 7 considers the prospects for constraining compact DM with this type of observation. Some of our results are discussed in Section 8, and we conclude in Section 9.

This paper is very much related to K18. While it presents the theoretical (lensing) and numerical (simulations) background of K18, the reader is directed to that paper for a detailed discussion of the particular Icarus and Iapyx events, including their interpretation. In this paper we refer to the Icarus and Iapyx events when appropriate or relevant for the discussion. Throughout the paper we assume a cosmological model with Ω_M = 0.3, Ω_Λ = 0.7, and H₀ = 70 km s⁻¹ Mpc⁻¹. For this model, 1'' = 6.45 kpc at the distance of the cluster MACS1149 (z = 0.544) and 1'' = 8.4 kpc at the distance of the background source (z = 1.49).

Besides CC, several terms will be used often in this paper. We refer to the CCs and caustics around microlenses as micro-CCs and microcaustics, respectively. The cluster CC and caustic that would form if there were no microlenses are the main CC and main caustic, respectively. Macro-images are the counterimages that would have formed if there were no microlenses and the lensing potential were sourced only by the cluster. A macro-image in a region filled with microlenses usually breaks up into smaller portions that we refer to as micro-images and also as bits. Then, around a CC we expect to find two macro-images, each composed of several smaller micro-images or bits. When the optical depth of microlenses is relatively small (${{\rm{\Sigma }}}_{\mathrm{microlens}}/{{\rm{\Sigma }}}_{\mathrm{crit}}\lt 0.01$) and the macro-images form very close to the main CC, the resulting group of micro-images is usually stretched along a straight line, following the direction of the cluster deflection field. Because of this geometry, we refer to this group of micro-images as a train of micro-images, or simply as a train. At low optical depth of microlenses, a background source will form typically two trains (or macro-images), one on each side of the CC. At higher optical depth, a single background source can form more than two trains, and each train can contain multiple smaller micro-images. In this sense, we can think of Icarus and Iapyx as unresolved macro-images which consist of even smaller bits or micro-images. The surface mass density of microlenses, Σ, is used in two contexts. In its broader sense we simply use Σ. When Σ takes the value of 7 M_⊙ pc⁻² (as found by K18 at the position of Icarus/Iapyx.²² ), we refer to it as Σ_o. Sometimes we express Σ in units of Σ_o and use f = Σ/Σ_o. Toward the end of this work we use another variable, F, that should not be confused with f. We use F to refer to the fraction of the total mass that is in the form of compact objects (whether this is made of stars from the ICL, PBHs, or both). By construction, F is always smaller than 1 while f can be larger than 1.

A gravitational lens is characterized by the lens equation

$$\begin{eqnarray}&&\beta =\theta -\alpha (M,\theta ),\end{eqnarray} \tag{ 1 }$$

where β is the position of the background source, θ is the observed position in the sky of the lensed image, and α(M, θ) is the deflection angle produced by the lens with mass M. The dependence of α(M, θ) on the position θ results in Equation (1) being nonlinear. Consequently, for a given position β, it is sometimes possible to find multiple solutions to Equation (1) with each solution representing a different counterimage. Each counterimage is magnified by a factor μ. Since lensing preserves the surface brightness of the background source, μ can be defined as the ratio between the observed size (i.e., area) of the counterimage, dΩ_θ, and the intrinsic size of the background source, dΩ_β. For a given lens model the deflection field α(M, θ) is known, and the magnification can be computed in a given position from the derivatives of the deflection field. The inverse of the magnification is defined as

$$\begin{eqnarray}&&{\mu }^{-1}={a}_{1}{a}_{2}=(1-\kappa -\gamma )(1-\kappa +\gamma )={(1-\kappa )}^{2}-{\gamma }^{2},\end{eqnarray} \tag{ 2 }$$

where κ and γ are the convergence and shear (respectively), and are combinations of the derivatives of the deflection field. We introduce the inverse of the magnifications, a₁ and a₂, that will be used later in this work. At a tangential CC, a₁ = 0. On each side of the CC, a₁ takes positive and negative values (parity). The sign of a₁ gives the parity of the image, so images with negative parity have a₁ < 0 and images with positive parity have a₁ > 0.²³

Counterimages that form near a CC can be magnified by very large factors. At the CC, the magnification diverges and dβ/dθ = 0. We can take advantage of this property to Taylor expand the lens equation around the CC,

$$\begin{eqnarray}&&\beta ={\beta }_{o}+\displaystyle \frac{d\beta }{d\theta }(\theta -{\theta }_{o})+\displaystyle \frac{1}{2}\displaystyle \frac{{d}^{2}\beta }{d{\theta }^{2}}{(\theta -{\theta }_{o})}^{2}+....\end{eqnarray} \tag{ 3 }$$

We choose β_o and θ_o as the origins of the reference systems in the source and image plane, respectively, and redefine β =β − β_o and θ = θ − θ_o. The second term cancels out at the position of the CC, leaving to second order

$$\begin{eqnarray}&&\beta ={\theta }^{2}/{\rm{\Theta }},\end{eqnarray} \tag{ 4 }$$

where we have defined the constant ${{\rm{\Theta }}}^{-1}=(1/2){d}^{2}\beta /d{\theta }^{2}$. At the position of the CC (θ = 0) we satisfy the condition μ⁻¹ = 0, and to first order ${\mu }^{-1}=d\beta /d\theta \propto \theta$. Hence, in the image plane we obtain for the total magnification (i.e., the magnification of the two images on each side of the CC)

$$\begin{eqnarray}&&\mu =\displaystyle \frac{{\mu }_{o}}{\theta }\end{eqnarray} \tag{ 5 }$$

near the CC, where μ_o is a constant that depends on the lens mass and geometry. This condition is satisfied for most lenses up to θ ≈ 1'' (see Figure 1). The asymptotic behavior when θ ≫ 0 is μ = 1 in the external side of the CC (a₁ > 0) and μ = 0 at the position of the lens for a point-source lens.

Zoom In Zoom Out Reset image size

The magnification in Equation (5) is expressed in the image plane. In terms of the position in the source plane, we can use Equation (4) to obtain

$$\begin{eqnarray}&&\mu =\displaystyle \frac{{\mu }_{o}/\sqrt{{\rm{\Theta }}}}{\sqrt{\beta }}.\end{eqnarray} \tag{ 6 }$$

The maximum magnification is obtained when the source touches the caustic—that is, when the distance from the center of the source to the caustic equals the radius of the source, R. Then we obtain

$$\begin{eqnarray}&&{\mu }_{\max }=\displaystyle \frac{{\mu }_{o}/\sqrt{{\rm{\Theta }}}}{2\sqrt{R}}.\end{eqnarray} \tag{ 7 }$$

Equations (4)–(6) are very useful for characterizing the properties of the counterimages near a CC. The values of μ_o and Θ can be obtained for a given lens and at a given position after fitting several positions near the CC. For a cluster like MACS1149 at the position of the Icarus event (and a background source at z = 1.49), μ_o ≈ 150'' and Θ ≈ 68'' for the model of D16. These values may change by as much as a factor of ∼2 for alternative models that still predict the CC in the correct position depending on the slope of the potential at the position of the CC, but we will adopt them below as realistic examples (see Section 8.1). For these values of μ_o and Θ, if the background star is a giant star with radius R = 100 R_⊙, the maximum magnification can be as high as μ_max ≈ 10⁶ at the CC near the Icarus position. If a background source is moving in the source plane with a constant velocity ${v}_{p}=d\beta /{dt}$ in a direction perpendicular to the caustic, the apparent observed velocity of the counterimages in the image plane is

$$\begin{eqnarray}&&{v}_{\mathrm{obs}}=\displaystyle \frac{d\theta }{{dt}}=\displaystyle \frac{{\rm{\Theta }}}{2{\mu }_{o}}{v}_{p}\mu ,\end{eqnarray} \tag{ 8 }$$

where we have used Equation (4) to relate θ with β and replaced θ with μ using Equation (5) after doing the derivative. Hence, at total magnifications μ ≳ 1000, the counterimages would appear to move at superluminal speeds (for this particular configuration). This has an interesting implication, since counterimages that move with larger apparent velocities can cover a larger region in the lens plane, probing more substructure as the source moves across the lens plane. Conversely, this can be seen as the microcaustics being more cramped in the source plane as we approach the main caustic.

If DM is clumped (like in some wave DM models), or if a significant fraction of it is made of compact objects such as PBHs, the probability that in a fixed period of time a given counterimage passes behind a clump of DM or a PBH will be higher near a CC, where counterimages probe the lens plane at a faster rate. The most interesting scenario to constrain the fraction of DM in compact objects is near tangential CCs. Radial CCs are normally close to the center of the cluster, where the ICL or the stars from the brightest cluster galaxy (BCG) can overwhelm the possible signature from compact DM.

In this work, we will focus on the case of tangential CCs, where the exploitation of crossing events may be most fruitful (see, however, M. Chen et al. 2018, in preparation, where a microlensing event near a radial CC in the same cluster, MACS1149, is used to unveil a supermassive BH ∼10 kpc from the center of the BCG).

Tangential CCs form when 1 − κ − γ = 0 (while $| 1-\kappa \,+\gamma | \gt 0$). If the lens plane is populated by small microlenses, they will contribute to the convergence (or surface mass density) with a small factor ${\bar{\kappa }}_{* }(\theta )\ll 1$, where ${\bar{\kappa }}_{* }(\theta )$ is the mean surface density of a point-like star with mass M within a radius θ; that is, ${\bar{\kappa }}_{* }(\theta ){{\rm{\Sigma }}}_{\mathrm{crit}}=M/(\pi {\theta }^{2})$. Near the CC, the condition for diverging magnification becomes $1-\kappa -\gamma ={\mu }_{t}^{-1}={\bar{\kappa }}_{* }(\theta )$, where $\mu \,={\mu }_{t}{\mu }_{r}={({a}_{1}{a}_{2})}^{-1}$. Hence, we can conclude that

$$\begin{eqnarray}&&{\theta }_{{\rm{E}}}=\sqrt{\displaystyle \frac{M{\mu }_{t}}{\pi {{\rm{\Sigma }}}_{\mathrm{crit}}}}.\end{eqnarray} \tag{ 9 }$$

The above result has profound implications. A microlens at a position near a CC, where the magnification is μ_t ≈ 1000, will behave (to first order) like an isolated microlens, but will be a thousand times more massive (see Figure 2). As shown in Section 8.2, in the last moments before a star crosses the cluster caustic, the magnification can become of order 10⁶, allowing the detection of substructures with masses comparable to a Jupiter mass (see also Equation (23) in Paczyński 1986, where an expression similar to Equation (9) is introduced as the dimensionless radius).

Zoom In Zoom Out Reset image size

Despite being based on some approximations, like neglecting higher-order terms, the expressions above seem to match remarkably well the results derived from detailed numerical calculations. Figure 3 shows the change in effective radius (defined as the perimeter divided by 2π) for a microlens that is isolated (no external field; bottom curve) and for a microlens that is embedded in a lensing potential with ${\mu }_{t}\approx 100$ (top curves). All curves grow with radius as $\sqrt{M}$, but in the case of the microlens in a lensing potential, the amplitude is increased by a factor $\sim \sqrt{{\mu }_{t}}$ as predicted by Equation (9).

Zoom In Zoom Out Reset image size

The magnification around a microlens in a field with external shear and convergence has previously been studied in detail (see, e.g., Schechter & Wambsganss 2002). In Appendix A, we present a brief and simplified description of a single microlens in an external field at high magnification.

The results presented in the previous section (see also Appendix A) give us useful insights into the behavior of the magnification around a microlens near a CC. However, in most realistic scenarios the CC region will be populated by a number of microlenses having, in general, different masses. In order to explore this more realistic regime, we resort to numerical simulations where the magnification is computed from simulated data.²⁴

We will assume the background source is a luminous giant star, which have radii ranging between ∼100 and ∼1000 R_⊙. We adopt the value R_star = 1000 R_⊙. The results presented in this paper are virtually the same for smaller stars, except for the maximum magnification reached when a microcaustic is being crossed (after the star touches a given microcaustic). In this case, the maximum peak magnification would grow by a factor $\sqrt{1000/{R}_{\mathrm{star}}}$ (see Equation (6)) if the microlens is at sufficiently large distances from nearby microlenses. The deflection field in the simulations contains a smooth component from the large-scale cluster potential and small-scale fluctuations from the microlenses. Normally, extremely luminous (background) stars can be found in star-forming regions where the density of stars may be relatively high. In this work we ignore the effect of neighboring stars and consider only the effect over one of those stars. If several background stars were moving on the web of caustics, each star would produce a series of peaks as they cross microcaustics. In realistic scenarios, only the brightest stars will produce peaks that can be measured, with the remaining stars contributing to a stochastic background of small fluctuations in the light curve. For instance, K18 argue that in order to explain the Icarus event, the background star needs to be extremely luminous and hence very rare.

For the smooth-scale potential we assume a realistic lens model, in particular the lens model of D16 for the cluster MACS1149 in the region of the Icarus and Iapyx events; see K18, where the same model was also used to interpret the observations. The model produces a CC (for a background source at z = 1.55) that falls in between the positions of Icarus and Iapyx. Following K18, in this work we assume that the CC is exactly between Icarus and Iapyx as predicted by various models (Richard et al. 2014; Oguri 2015; D16; Kawamata et al. 2016). However, the reader should note that other interpretations are also possible where, for instance, the CC is closer to Icarus than to Iapyx, in which case the two events would be produced by two different background stars (see K18 for a more detailed discussion of this and other alternative interpretations). The hypothesis that there are multiple bright stars in the background moving between microcaustics would also be supported if the third event of K18 is confirmed as an additional microlensing event. In this case, since this new image does not appear aligned with the direction where counterimages of the same star are expected to form, a second background star would be needed.

Alternatively, when a fast computation is needed over a large area (for instance, in Figure 1), we use the analytical model from Blandford & Kochanek (1987) for the smooth component, where we tune the lens parameters to produce a magnification pattern similar to the model above. In Section 3.2, some of the results computed in the source plane come from a small region at very high resolution. For this particular case, we use a simplified model for the macromodel that is given by just two parameters—the surface mass density (κ_smooth) and shear (γ_smooth). These two parameters are considered constant, which is a valid assumption given the small area being simulated. Tuning κ_smooth and γ_smooth to the desired values allows us to quickly produce simulations with a variety of magnifications from the macrolens model.

The microlenses are assumed to be point masses and are randomly distributed. Unless otherwise noted, the masses of the microlenses are drawn from a Spera et al. (2015) initial–final mass function where the only surviving stars in the intracluster medium are less massive than 1.5 M_⊙ (above this mass, the remnants of more massive stars are also included in the simulation). The same model is also discussed by K18 together with other alternative models. The mass function is normalized to match the inferred stellar surface mass density at the Icarus position (as estimated by Morishita et al. 2017).

We place stars in a region (or extended region, hereafter) which is slightly bigger than the final simulation region (or target region, hereafter). This is done in order to minimize edge effects. The extended region contains the target region plus buffer zones around it, of 0.2 mas width each extending in the vertical direction. (The left and right edges of the simulation are not used for the computation of the light curves, so we do not add an extra buffer on these two edges.) This buffer zone is sufficiently large to account for the individual effect of the largest microlenses that could be found beyond (but near) the edge of the target simulated region. The target region is a band of width 1 mas and length 10 mas, aligned in the direction where counterimages form and it is contained in the extended region of width 1.4 mas and length 10 mas. The total number of microlenses included in this extended region is 18,686, and they are placed randomly within the extended region. The total mass of the microlenses in the extended region is ∼4000 M_⊙. When PBHs are included in the simulation, we use the same distribution of stars and add the effect of randomly placed PBHs. The number of PBHs is determined by the fraction of total mass that is in the form of PBHs. This number scales as N_PBH ≈ 30 F_PBH per mas², where F_PBH is the fraction (in percent) of mass in the form of 30 M_⊙ PBHs. This results in 420 F_PBH PBHs (of 30 M_⊙ each) in the extended region. Finally, we subtract from this extended region the contribution to the deflection field from a smooth mass distribution with the same surface mass density as the microlenses (stars plus remnants, or stars plus remnants plus PBHs), so the total surface mass density remains constant.

The simulations are made at a resolution of 1 μas in the image plane. As mentioned earlier, the target region is a band of width 1000 pixels and length 10,000 pixels in the direction where counterimages are expected to form (i.e., at an angle α_c ≈ −40°). The length of the simulated box (∼10 mas/cos(α_c)) maps into a corresponding length in the source plane of ∼100 μas. This is enough to follow a moving background source at z ≈ 1.5 with v ≈ 1000 km s⁻¹ over ∼1000 yr. This resolution is sufficient to resolve the elongated arcs that form when a background star crosses a microcaustic, if the background star is at least a few tens of solar radii.

We assume that the source is moving perpendicular to the main caustic of the cluster. The simulated light curves have a weak dependence on this angle, since the microcaustics are stretched by a large factor in the direction of the main caustic of the cluster. Only if the background sources are moving in a direction very close to the main axis of the microcaustics would the simulated light curves be significantly different—but this is unlikely since the probability of moving in this narrow range of angles is small. If the source is not moving perpendicular to the main axis of the caustics but at a different angle, the simulated light curves would still be very similar to those presented in this work, except that the time it takes for a source to cross a microcaustic would be stretched by a factor cos(α_sc)⁻¹, where α_sc is the angle between the main axis of the microcaustic and the direction of motion of the source. The reader can find videos at https://cosmicspectator.org/2017/06/30/dark-matter-under-the-microscope/ (movies 5 and 6) extracted from the simulations and showing the effect of the motion of a source as it travels through a web of caustics that is moving parallel to, or at an angle with, the main axis of the microcaustics. A source moving parallel to the caustics may be interpreted as moving in a region with a small surface density of microlenses (see also Section 3.2). Also, if the velocity is very small, it may be erroneously interpreted in a similar way.

When the lensed images form farther from a micro-CC, the dimension of the micro-images is typically smaller than the pixel size of the simulation. In this case (but also when the micro-image is resolved), the real dimension of the micro-image is computed at the subpixel level (making use of approximations that allow resolving scales much smaller than the simulation pixel). This is achieved by interpolating the deflection field so any position in the source plane can be mapped into the corresponding interpolated position in the image plane, effectively achieving infinite resolution. Simple, fast interpolations are sufficient because the deflection field is extremely smooth. The smoothness of the deflection field is guaranteed since it is simply the superposition of the deflection field from the cluster and that from the microlenses. The former is orders of magnitude larger than the latter, so the small perturbations from the microlenses do not break the smooth condition needed for the simple interpolations. The only place where the simple interpolation may break down is when one is looking for counterimages very close to the microlens, since in this case the deflection field diverges. Luckily, these positions correspond to the lowest magnification regions, so those counterimages can never be observed.

The magnification is computed as the ratio of the total area in the image plane divided by the area of the background star in the source plane (i.e., $\pi {R}_{\mathrm{star}}^{2}$), and we neglect limb-darkening effects (this would add a small correction during a caustic crossing event that is most important in its last moments). We also neglect interference effects, since both the background stars and the microlenses are sufficiently large.

3.1. Multiple Microlenses with the Same Mass

The large magnifying power of a galaxy cluster near its CC can allow for detailed study of both the background objects and the substructure in the lens plane itself. A point-like microlens with mass M (like a star or a BH) in the lens plane will behave like a lens with effective mass μM (see Equation (9)). In the simple scenario where the microlens is isolated (i.e., no other microlenses nearby), the magnification μ of a cluster at a distance less than 01 from the CC can easily reach values above μ = 1000 for a point-like background source. At this magnification, a background compact bright object such as a giant star will be boosted by ∼7.5 mag. This boost factor may be sufficient to make luminous stars at z > 1 detectable with deep observations. In this situation, a microlens with mass M = 10⁻² M_⊙ in or near the line of sight to the background star and close to the cluster CC will behave (in terms of its lensing effect) like a microlens with mass M > 1 M_⊙. This makes it possible to detect the microlens in the light curve of the background source.

If no microlenses are present in the lens plane, on small scales (less than 01) the cluster CC can be approximated by a straight line (see the left panel in Figure 4), and the magnification grows as the inverse of the distance to the CC (Figure 1). The corresponding caustic is equally well described by a straight line, but the magnification grows as the inverse of the square root of the distance to the caustic. Hence, a large magnification of several thousand requires an incredibly small separation between the background star and the caustic, making this type of configuration very rare, and observing a caustic crossing very unlikely (see the middle panel in Figure 4, where we show the incredibly narrow region in the source plane that maps on the image plane in the left panel).

Zoom In Zoom Out Reset image size

A cluster caustic crossing event is expected to be very short lived (several hours or a few days, depending on relative velocity and star radius) and involves very large magnifications when the lens plane contains no microlenses. The dependence on the square root of the separation between the star and the caustic means that the precise moment of the caustic crossing event can be predicted, since the observed flux evolves as $1/\sqrt{t-{t}_{{\rm{o}}}}$, where t is time and t_o is the time of crossing. When microlenses are included in the lens plane, the situation can be very different since microlenses can significantly disrupt the CC (see Figure 5). The disruption is most prominent near the CC and decreases with the distance to it (see Figure 4). In the source plane, the corresponding caustic region gets expanded by a factor that depends on the number of microlenses and their masses. A larger caustic region means that observing a microlensing event becomes more likely, since a moving background star may intersect multiple microcaustics in a given period of time (as opposed to intersecting just the main CC of the cluster). One of these microcaustics can also be intersected many years before the source crosses the position of the main caustic. This translates into a dramatic increase in the probability of seeing caustic events before (but also after) crossing the position of the main caustic.

Zoom In Zoom Out Reset image size

3.2. Source Plane Interpretation

The plots in the previous section show the magnification pattern in the image plane. However, the magnification in the image plane for individual micro-images is normally not observed (unless the number density of microlenses is very small). To better understand how the magnification of the multiple images works, it is more useful to represent it in the source plane. The mapping between the image and lens planes is done through the lens equation.

In this subsection we present a small portion of the source plane computed at higher resolution than the simulations used in the bulk of this work. The higher resolution is attained by interpolating the deflection field from the main simulation to smaller scales. Figure 6 shows the image plane, source plane, and a cross-section of the magnification in the source plane for both sides with negative and positive parities. The image planes display the characteristic hourglass shapes discussed in previous sections. The source plane shows the familiar overlapping and stretched diamond shapes. The side having negative parity clearly reveals the gap with low magnification in the central region of the diamond-shape caustics. This is a familiar result found in earlier work (see, for instance, Chang & Refsdal 1984; Schechter & Wambsganss 2002). This gap of low magnification results in periods of low flux in the observed light curves, on the side with negative parity. The bottom plot in each panel shows the cross section along the diagonal line in the source plane.

Zoom In Zoom Out Reset image size

The circles in the source plane represent a source with a radius of ∼1.5 μarcsec. That source would simultaneously see the caustics from the negative and positive parity sides, but we have separated them here for clarity. The corresponding amplified image is shown as an ellipse (with small distortions) in the image plane. A source of this size would produce only one counterimage on the side with negative parity, and another one on the side with positive parity (i.e., only two macro-images and no micro-images) since the size of the source is significantly larger than the characteristic scale of the microcaustics. Also, the counterimage on the side with negative parity could not be hidden by microlenses with masses similar to those in the simulation, since portions of the source would always overlap with regions of high magnification in the source plane.

Note how a source that is significantly larger than the width of the microcaustics produces counterimages that are notably less magnified on the side with negative parity. This phenomenon would not take place if there where no microlenses, since in that case the magnification within the circular region would be very similar for both parities. This can be understood if we integrate the total magnification within the circle in the source plane. The gap between the caustics results in a smaller total magnification in the enclosed area. Consequently, relatively small regions in the source plane, of angular size a few times the typical size of the caustics (like the circular regions in Figure 6 or R ≈ 1.5 μas ≈ 0.01 pc at z ≈ 1.5), could show a ratio in the flux between the two counterimages (positive/negative parity) of ∼1.3 (although the exact value depends on multiple factors like source size, mass of microlenses, distance to the macro-CC, etc.). A similar property has been exploited in the context of QSO microlensing, as for instance by Mediavilla et al. (2017). An explanation of this phenomenon is given by Schechter & Wambsganss (2002), where the authors demonstrate with a simple toy model how macrominima and macrosaddle points can have their magnifications affected significantly in the presence of microlenses.

Finally, movies 7 and 8 at https://cosmicspectator.org/2017/06/30/dark-matter-under-the-microscope/ show how the observed magnification depends on the point where the trajectory of the background star intersects the microcaustics.

In order to get a better view of the source plane in the different scenarios, we make an ensemble of alternative simulations at even higher resolution where we vary both the magnification of the macromodel and the surface mass density of microlenses. We adopt the same redshift for the cluster lens and background source as in the Icarus and Iapyx events. For this particular set of simulations, we adopt a simplified model where the surface mass density from the macromodel (κ_m) and the shear (γ_m) are fixed, instead of adopting the model from D16 used in the main simulations. This is a valid approximation since the simulated region is very small. By varying κ_m and γ_m, we can easily simulate a given region of the lens plane with the desired magnification from the macromodel. For simplicity, we also assume that the component of the shear in the vertical direction is zero (that is, ${\gamma }_{m}=\sqrt{{\gamma }_{1}^{2}+{\gamma }_{2}^{2}}={\gamma }_{1}$), so the deflection field has its main component in the horizontal direction. Since both κ_m and γ_m can be expressed in terms of the derivatives of the deflection field (α), these relations can be reversed, and we can also express the derivatives of the potential as a function of κ_m and γ_m,

$$\begin{eqnarray}&&{\alpha }_{x}^{x}=\displaystyle \frac{1}{2}({\kappa }_{m}+{\gamma }_{m}),\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&{\alpha }_{y}^{y}=\displaystyle \frac{1}{2}({\kappa }_{m}-{\gamma }_{m}),\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&{\alpha }_{x}^{y}={\alpha }_{y}^{x}={\gamma }_{2}=0,\end{eqnarray} \tag{ 12 }$$

where ${\alpha }_{i}^{j}$ is the derivative of the i component of α with respect the coordinate j.

With these equations we can describe the deflection field (except for an irrelevant constant) of the macromodel. The deflection field from a population of microlenses is added linearly. When considering microlenses, instead of κ_m one needs to use ${\kappa }_{m}\,\mbox{--}\,{\kappa }_{* }$, with κ_* the convergence from the surface mass density of microlenses. This guarantees that the total convergence (cluster plus microlenses) is equal to the target κ_m.

The specific values of κ_m and γ_m are determined by the value of the magnification to be simulated, μ_m = μ_t × μ_r. One can easily find that for the side with negative parity ${\gamma }_{m}=({\mu }_{r}^{-1}+{\mu }_{t}^{-1})/2$ and ${\kappa }_{m}=1-{\gamma }_{m}+{\mu }_{r}^{-1}$, while for the side with positive parity we have the condition ${\gamma }_{m}=({\mu }_{t}^{-1}-{\mu }_{r}^{-1})/2$ and ${\kappa }_{m}=1\,-{\gamma }_{m}-{\mu }_{r}^{-1}$. We vary μ_r, μ_t, and κ_*, producing a set of simulations for the cases with positive and negative parities. The simulations consider a large circle of radius 0.465 mas where we place the microlenses randomly. The pixel scale is 31 nas and we compute the total deflection field in a narrow horizontal band of width 744 μas and height 93 μas. By construction, this area maps into an area a factor μ = μ_r × μ_t times smaller in the source plane. Without loss of generality, we fix μ_t = 1.5, so in the source plane and to first order, the simulated region maps into an area of dimension 93/1.5 × 744 μ/1.5 arcsec². As discussed below, at high optical depth, the effective magnification is smaller than that for the macromodel, so the simulated region can be larger than that inferred from the values above.

By inverse ray tracing, we compute the magnification in the source plane after interpolating the original deflection field to achieve effective resolutions of ∼3 nas per pixel (or ∼1000 solar radii at z = 1.5). A small area of the source plane is shown in Figure 7 for each simulation. The results in this plot are divided into two groups. On the right side we show the source plane at fixed κ_* but varying magnification. On the left we show the simulated source plane at fixed magnification but varying κ_*. The magnification considered for this example is moderate (μ = 30), but it serves our purpose as it better shows the structure of the caustics in the source plane. For larger magnifications the behavior would be qualitatively similar to what is shown on the left block of Figure 7. Since we have the sides with opposite parity projecting back into the same region in the source plane, at a given pixel in the source plane one would get a bundle of rays (from the inverse ray tracing method) coming from the side with positive parity and a different bundle coming from the side with negative parity. The mapping between the image plane and the source plane around a CC can be visualized as a sheet of paper being folded by its middle point (CC). In the source plane, the fold represents the caustic, with the two halves of the sheet forming two overlapping planes. A source will project into these two planes; when unfolded (i.e., the image plane), the sheet of paper will show two images which are symmetric with respect to the folding line.

Zoom In Zoom Out Reset image size

To better illustrate the differences between the sides with positive and negative parities, we show the source plane for each of the two planes within it, described at the end of the previous paragraph, and also compute the statistics of each plane separately. This makes sense when comparing with observations, since the statistics of the observed lensed image depends on the parity, as we show later, and has been demonstrated in earlier work (e.g., Schechter & Wambsganss 2002). The top row shows the plane with positive parity while the bottom row shows the plane with negative parity, where the characteristic microsaddle points with low magnification can be appreciated clearly. The two columns with κ_*/κ = 0.05 and κ_*/κ = 0.15 in the left block contain a stellar component consistent with the upper limit in K18 (that is, with 19 M_⊙ pc⁻² or κ_*/κ = 0.012 and a Salpeter spectrum in the low-mass regime) plus microlenses with 30 M_⊙, mimicking a monochromatic population of PBHs. As mentioned above, the value 19 M_⊙ pc⁻² is motivated by the updated estimate in K18. This is almost three times more mass than the value of Σ_o = 7 M_⊙ pc⁻² used in the rest of this work, but the careful reader will notice that κ_*/κ is instead a factor of four larger than the value of 3% corresponding to Σ_o. This is due to the fact that for this set of high-resolution simulations with 19 M_⊙ pc⁻², the value of κ considered is 0.66 instead of the value from the model in D16 (κ = 0.9) used in the rest of this work. κ ≈ 0.66 is the value required by μ once μ_t is fixed to 1.5 and μ_r ≫ 1 as described above. The left panel with κ_*/κ = 0.004 has approximately four times fewer microlenses than the model with κ_*/κ = 0.012 and no PBHs (of 30 M_⊙). The last column with κ_*/κ = 0.37 represents a population of PBHs but with a power-law spectrum (for the mass function) similar to that used to simulate the stellar component from the ICL.

On the right side of Figure 7 we show similar plots, but this time the surface mass density of microlenses is fixed (to a value consistent with the upper limit on κ_* in Kelly et al. 2018; κ_*/κ = 0.012 and no PBHs) and we vary the magnification. The case with μ = 30 can be compared directly with the cases presented in the left block.

Figure 7 makes the degeneracy between κ_* and μ evident. When κ_* is sufficiently high, the source plane saturates with overlapping caustics. The same is also true for moderate values of κ_* but at large magnifications. In this case the overlapping of the caustics is produced by the high magnification. In order to quantify the differences, we compute the probability distribution function (PDF) of the magnification in the source plane from these simulations.

The result is shown in Figure 8. The left plot shows the PDF for the case with fixed magnification and varying κ_*. For small values of κ_*, the PDF shows a clear peak near the magnification of the macromodel, μ_m. The side with positive parity peaks at slightly smaller values than μ_m while the side with negative parity peaks at values very close to μ_m. At high magnification, both sides behave very similarly, with the tail of the PDF falling like μ⁻³, typical of isolated lenses. At low magnifications, the differences are significant between the two parities. When κ_* is increased, we observe that the differences between the two parities grow as well. The peaks of the PDF separate more, with the one from positive parity clearly below and that from the negative parity clearly above μ_m. At high magnifications, there seems to be an excess of probability on the side with negative parity with respect to the side with positive parity. Also, the probabilities of having high and low magnification increase as κ_*.

Zoom In Zoom Out Reset image size

As we increase κ_*, the peaks in the PDF disappear, and for sufficiently large values of κ_*, the PDFs of both parity sides start to resemble each other. This is the saturation regime at which the notion of sides with positive and negative parity loses its meaning (the sign of a₁ = 1 − κ − γ can adopt positive and negative values on both sides of the main CC). The right plot in Figure 8 shows the case of fixed κ_* and varying magnification. The model with μ_m = 30 (black curves) is the same as the black curve model in the left plot. We observe a similar trend, but now the excess of probability at high magnification in the side with negative parity is more evident (specially at μ = 150 and μ = 600). Also, for μ_m = 2400, the PDF of the sides with negative and positive parity are almost similar and deviate from μ⁻³ at high magnification. Instead, the PDF resembles a log-normal distribution, which typically appears in multiplicative processes. This is attained by the combined effect of multiple overlapping caustics and the high magnification of the macromodel. Finally, we note that, similarly to what happens in the left panel, when we increase the magnification of the macromodel, the average magnification in the simulated region deviates from the one we would have obtained from the macromodel (i.e., with no microlenses). In particular, we find that in the presence of microlenses with κ_*/κ = 0.012 and in the side with positive (negative) parity, the averages of the magnifications are 29.9 (30.2), 150.4 (144.5), 564.8 (566.3), and 1455.7 (1345.3) for macromodel magnifications of 30, 150, 600, and 2400, respectively.

From results like those shown in Figure 8, one can extract important properties of the magification, but they do not contain all the information. The magnification is highly non-Gaussian, as shown in this figure, and hence the PDF alone gives an incomplete picture of the problem. For instance, the PDF plots shown in Figure 8 do not account for the correlations that are evident in Figure 7. Higher-order statistics like the correlation function or power spectrum are useful discriminators in this type of situation.

In this section we study in more detail the effects of microlenses at the position of, or very close to, the CC. For the microlenses we adopt as a reference the Spera et al. (2015) model normalized to ${{\rm{\Sigma }}}_{o}\approx 7\,{M}_{\odot }$ pc⁻² (similar to the surface mass density in surviving stars inferred at the position of Icarus).

When microlenses are present in the vicinity of the CC, the infinitesimally narrow CC widens, with overlapping critical lines that form a complex network (see, e.g., Figure 5). This network extends up to a maximum range that depends on the total number and masses of the microlenses. In the case of Figure 5, this network extends well beyond the displayed field of view. The change in the network when the amount of microlenses is varied is made more evident in Figure 9. The magnification pattern gets shifted around, with regions of high magnification becoming regions of low magnification, and vice versa. When the amount of microlenses is small (i.e., small fraction f = Σ/Σ_o), the main CC becomes sharper by trading high magnification for lower magnification with the surrounding area in the lens plane. As more microlenses are added, the disruption becomes more serious, and at some point between f = 0.003 and f = 0.01 in Figure 9 the CC itself transforms into a network of micro-CCs. The extension of this network around the main CC marks the region where microlensing events are more likely to be observed. An interesting consequence is that the typical magnification one would expect is changed by the addition of microlenses. In Figure 10 we show the median of the distribution of magnifications (of individual micro-images) computed at different distances from the position of the main CC. For each distance, the distribution of magnifications and its median are computed in an area of 1.9 × 0.5 mas². Adding microlenses thus results in a reduction of the typical magnification of micro-images that one would have obtained without them. This median magnification, however, cannot be normally observed, since the lensed image forms (typically) an unresolved train of micro-images and what we observe is the total flux of all micro-images (an exception being at low optical depth, where the the total flux is usually given by one micro-image). We show later, however, that if the lens plane is populated by massive microlenses (a few tens of solar masses), the separation between micro-images can reach a few milliarcseconds, opening the door to future high-resolution observations of the individual micro-images.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

This change in the magnification is also evident in Figure 5, where we display a small region around the main CC in the case of the smooth lens compared with the magnification pattern when microlenses are added (diagonal band). The figure also shows a lensed background object (with ∼0.01 pc radius), or train of micro-images, at the moment of maximum magnification. The lensed image breaks up into multiple smaller components. For smaller background sources (such as a large star), the lensed image would break up into even more smaller pieces.

When f is sufficiently small, the effect of the microlenses is small and the magnification behaves like in the smooth lens model case, except when we approach the CC. At short distances from the CC, even small microlenses can have a significant impact on the magnification pattern. As f grows, the range at which the CC gets disrupted grows as well. For values of f ≈ 0.001 the disruption is still significant up to scales of a few milliarcseconds. In this situation, if macro-images are being formed on both sides of the main CC at a distance of a few milliarcseconds, a telescope like the Hubble Space Telescope (HST) observing the unresolved macro-images would start to see not only a change in flux over time but also a change in the observed position of the peak, since the observed images would appear to be jumping back and forth between the two sides with opposite parity. For values of f ≈ 1 it is impossible to determine the exact location of the main CC and the magnification pattern is completely disrupted over a scale of hundreds of milliarcseconds. If f ≫ 1 the disruption can extend to scales of an arcsecond and these types of microlensing events would be much more common. The fact that no similar microlens event has ever been reported before Icarus/Iapyx is a simple indirect indication that the optical depth of microlenses cannot be much higher than that from the stellar component (with κ ≈ 10⁻³–10⁻²).

The results presented in the previous section show examples of the magnification pattern near the main CC and in a narrow region of the source plane when microlenses are present. However, observations can only sample the magnification in the source plane (unless micro-images are resolved). Through observations we can measure the total flux as a function of time—that is, the light curve.

To simulate the light curves, we place a background star in the source plane moving with a relative velocity of 1000 km s⁻¹ toward the main caustic. The results presented in this section can easily be rescaled (stretched or compressed) to any velocity (v) by the factor v/1000. A relative velocity of 1000 km s⁻¹ is a reasonable assumption given the redshifts of the lens and source. Figure 11 shows a small segment of the simulated light curves for Icarus (blue) and Iapyx (red). We assume the Spera model (left panel) and the Spera+PBH(10%) model (right panel). In the case when microlenses are only ICL stars (left panel), the light curves for Icarus and Iapyx can be very different when microlenses populate both sides of the main CC. Macro-images on the Iapyx side (a₁ < 0) can disappear for periods of 10 years or more.

Zoom In Zoom Out Reset image size

This is a consequence of the low-magnification regions present between the semi-diamond-shape caustics discussed in Section 3.2, on this side of the main CC.

As shown by Schechter & Wambsganss (2002), macrosaddle points (like those on the side with a₁ < 0) can be fainter because they can lack microminima, while macrominima must have at least one microminimum.

Another interesting difference is the amplitude of the peaks, which can be higher on the Iapyx side (see Section 3.2, where the same effect is observed in the PDF). This tradeoff between low-magnification periods and higher peaks conserves the total integrated flux when integrating over long periods of time. We find that for surface mass densities of microlenses comparable to Σ_o, the average of a light curve converges (both on the side with positive and negative parities) toward the value of the macromodel when averaging the light curve over a few hundred years. The right panel of Figure 11 shows the corresponding light curve when 10% of the mass in the lens plane is substituted by PBHs with a mass of 30 M_⊙ each (this case also includes the microlenses from the ICL). In this case, the light curves on both sides of the main CC are more similar, and we do not observe periods of low magnification on the side with a₁ < 0. We also note that when the fraction of PBHs is high, the clustering of PBHs introduces large-scale temporal and spatial correlations in the magnification pattern that can result in long periods of relatively low or high magnification, as shown in the right panel of Figure 11. As we show below in Section 7, when the optical depth of microlenses is sufficiently high (for instance, when 10% of all mass is in microlenses), we are in the saturation regime and the properties of the light curves must be similar. We also notice that the peaks on both sides are also smaller as a consequence of the reduction in Einstein radius. The associated Einstein radius of the microlenses no longer scales like ${(M{\mu }_{t})}^{1/2};$ more precisely, the effective μ_t is smaller when we reach the saturation regime. This is an interesting result since it implies that an event like Icarus would require, on average, a brighter star if a significant fraction of the DM in the lens plane is made of PBHs. At the same time, hiding Iapyx for at least 10 years seems unlikely in this scenario and would require the presence of a second (and fainter) star in the source plane to explain Iapyx.

This argument is in agreement with the probabilities estimated by K18, which show that a scenario where a sizable fraction (more than a few percent) of the DM is in the form of compact lenses is less agreeable with the data than a scenario where the microlenses are just the ICL stars. Also, K18 found that simulations suggest that the existing data have sensitivity to the mass function of stars and remnants (given the assumptions made about the magnification and stellar mass density). However, the observed light curve presented by K18 does not have as many data points as one would desire to achieve good constraining power. A clear discrimination between different models (other than a preference for models with a low fraction of PBHs or microlenses forming binary systems as discussed by K18) is not yet possible. More data are needed (more peaks, better cadence, smaller photometric errors) in order to clearly distinguish between the different possible scenarios. We hope that regular monitoring of this cluster will provide such data in the near future.

If we imagine we can monitor Icarus and Iapyx for ∼10,000 yr, we would witness the moment where both events merge and disappear as the background star crosses the last microcaustic (see Section 6). If the deflection field were perfectly smooth (without microlenses), the light curves of Icarus and Iapyx would be featureless and the flux would increase as $1/\sqrt{t-{t}_{o}}$, as the background star approaches the main caustic, and would reach a magnification of ∼10⁶ in the last moments before disappearing. Conversely, if the lens plane near the CC is populated with microlenses from the ICL stars, PBHs, or both, the light curves will be rich in features like those shown in Figure 11. In this scenario, we will observe hundreds or thousands of peaks, each having a magnification between ∼10³ and ∼10⁴. The total flux integrated over the 10,000 yr will be the same, but the light curves will be very different.

It is also interesting to track the fluxes of the individual micro-images. In Figure 12 we compare a fragment of the simulated light curve (for Icarus) with the magnification in the lens plane at the position where the brightest micro-image is formed. During the valley periods, both magnifications coincide because there is only one dominant micro-image (there are other counterimages but with significantly smaller fluxes), so the total flux of the macro-image is basically given by the flux of the dominant micro-image. During a peak, new bright micro-images form around the micro-CC. The total flux is then typically larger by a factor 2 or 3 than the flux of the brightest micro-image. The result shown in Figure 12 demonstrates how, at moderate optical depth, during most of the time the train of micro-images is very compact in size (only one dominant micro-image during the more frequent valley periods). However, during a peak event, it may be possible that the separation between different portions of the train may be large enough so the overall size of the train of micro-images could be distinguished from that of a point source. This will be explored in more detail later in this paper.

Zoom In Zoom Out Reset image size

5.1. Estimating Rates

Based on our current understanding of the Icarus/Iapyx events, we should expect to see more peaks in the future at the position of these events. Assuming the source is heading toward the main caustic, for transverse velocities on the order of 1000 km s⁻¹, it might take ∼10,000 yr for the source to cross the main caustic and become unobservable. Conversely, if the background source is moving away from the main caustic, we are witnessing the lensed source ∼10,000 yr after it first crossed the position of the main caustic. Based on our model for the underlying microlenses and the cluster potential, we can predict how the last moments (or first moments) of the multiply lensed images will take (or took) place. For simplicity we assume that the source is heading toward the main caustic from now on. The results presented in this section can be reversed to describe the situation where the background source is moving in the opposite direction—that is, away from the main caustic.

Figure 13 shows the light curve moments before and after the background source crosses the main caustic. When approaching the last moments, the macro-images are broken into multiple groups or trains. Each train comprises a bright micro-image surrounded by smaller ones, usually aligned in the direction of the cluster deflection field. The last moments in the light curve show how the individual trains of micro-images disappear at different moments. When one of the trains vanishes the total flux drops by a large fraction. The last surviving micro-images have small fluxes and roam close to the largest microlenses (that is, in regions with low magnification). In principle, every single microlens (assuming they are truly point sources) would have a small micro-image near its center since the deflection angle of a point source diverges at the position of the microlens (there is a nice visual demonstration of this effect by Lewis et al. 1993 and Witt 1993). These small micro-images can survive a much longer time than the trains. However, these low-magnification images have such small fluxes that they can be neglected in terms of their contribution to the observed flux.

Zoom In Zoom Out Reset image size

Figure 13 also shows the maximum extension, or size, of the set formed by all micro-images. We compute this size based on the micro-images that have magnification larger than 50. At time ∼50 yr in the left plot, there is only one surviving micro-image with magnification above 50 and the apparent size of the train drops to zero (a second micro-image appears and vanishes quickly at times ∼58 and ∼65 yr). The size evolves as $\sqrt{| t-{t}_{o}| }$, which is inversely proportional to the apparent velocity of the macro-images. Although not shown in the plot, there is a long tail past 100 yr where micro-images with small magnification survive for a long period close to the largest microlenses (with μ ≈ 1 or below). This tail and the evolution of the size can be better appreciated when the contribution from microlenses is smaller.

In the right panel of Figure 13 we show a similar plot, but this time the masses of all the microlenses are divided by a factor 30. Even though this may not represent a realistic scenario, it does mimic a situation where the optical depth of the microlenses is a factor 30 times smaller than in the left panel of Figure 13. (The deflection field scales as the mass; hence, rescaling the masses is a much faster way of producing multiple realizations with different lensing optical depths). A more realistic simulation corresponding to a factor 30 times smaller surface mass density would contain fewer (but heavier) stars and the right panel in Figure 13 would have fewer (but more pronounced) peaks.

At low optical depth (right panel of Figure 13), the light curve in the last moments starts to resemble the expected behavior of the classic caustic crossing event with a smooth model (see, for instance, Miralda-Escude 1991), with the typical $1/\sqrt{| t-{t}_{o}| }$ change in flux (red solid line). In this case the relative separation between the two trains of micro-images follows closely the $\sqrt{| t-{t}_{o}| }$ law and could be extrapolated to earlier times to predict the crossing time from real observations or constrain the relative velocity between the background source and the caustic. A background source close to the time of caustic crossing will produce two macro-images separated by a distance that is proportional to $\sqrt{| t-{t}_{o}| }$. The proportionality constant involves the relative velocity between the source and the caustic. Future high-resolution observations (for instance, with the Extremely Large Telescope and a maximum expected resolution of ∼1 mas) may be able to discriminate between point sources and extended sources with sizes larger than a few milliarcseconds. Monitoring the change in the size of the unresolved pair of macro-images can be used to constrain the relative velocity between the source and caustic. Similarly, the rise in observed flux in the last moments exhibits a similar dependency with the relative velocity between the source and the caustic. From Figure 13, the period of flux decline after the peak is sensitive to the number and type of microlenses. Nevertheless, in the presence of microlenses, a caustic-crossing event needs to be monitored regularly for several decades to produce light curves that capture the rise and decline of the flux as the caustic network is being crossed. This will require an effort similar to the production of light curves for QSOs and other variable objects.

A popular candidate for DM is PBHs (see, e.g., Carr et al. 2016a, for a recent discussion). They are formed in the first moments of the universe and, other than through gravitational interactions, they do not play any significant role in nucleosynthesis or baryonic physics after inflation (with the exception of very small PBHs that can evaporate quickly and inject energy into the universe, or the very massive ones that can transform mass into energy through their associated accretion disks). Other than in these extremely low- and high-mass regimes, PBHs would be very elusive, interacting only through gravity with the surrounding matter. A wide range of masses has been excluded for PBHs as a significant contributor to the DM. PBHs below a certain mass (∼10¹¹ kg) would have evaporated by now (Hawking 1974). PBHs with slightly higher masses would be strong sources of γ-rays.

Current γ-ray observations have ruled out this mass range (Barnacka et al. 2012; Carr et al. 2016b). At the other extreme, a large number of massive PBHs (M > 100 M_⊙) would disrupt globular clusters or impact the baryonic physics in the early universe. At intermediate masses, there are constraints of varying strengths from pulsar timing or microlensing (Kashiyama & Seto 2012; Niikura et al. 2017). A window where the constraints on the fraction of PBHs still need to be improved is $M\approx 30\,{M}_{\odot }$. Bird et al. (2016) made the interesting suggestion that the first LIGO event (that involved two BHs with masses in this range) could have been the coalescence of two PBHs in this mass regime (see also Sasaki et al. 2016; Clesse & García-Bellido 2017). If true, this would offer a natural explanation for the extreme BH masses measured by LIGO (Abbott et al. 2017; The LIGO Scientific Collaboration et al. 2017).

In this section we explore the possibility of constraining the fraction of PBHs with a caustic crossing event. If PBHs make up a significant fraction of the DM, they would disrupt the caustic. As seen in Figure 9, the disruption is proportional to the optical depth of microlenses. This scaling can be easily derived from basic principles (for simplicity, we assume all microlenses have similar masses). The optical depth, τ, can be estimated as the number of microlenses, N, times the area of their associated Einstein ring ($\pi {\theta }_{{\rm{E}}}^{2}$) per unit area, A, as

$$\begin{eqnarray}&&\tau =\displaystyle \frac{N\pi {\theta }_{{\rm{E}}}^{2}}{A}.\end{eqnarray} \tag{ 13 }$$

If we introduce the surface mass density of microlenses with mass M, Σ = NM/A, and make use of Equations (9) and (5), we get

$$\begin{eqnarray}&&\tau =\displaystyle \frac{{\rm{\Sigma }}}{M}\displaystyle \frac{4\pi {GM}}{{c}^{2}}\displaystyle \frac{{D}_{{ds}}}{{D}_{d}{D}_{s}}{\mu }_{t}=(4.2\times {10}^{-4}){\rm{\Sigma }}({M}_{\odot }{\mathrm{pc}}^{-2})\displaystyle \frac{{a}_{2}{\mu }_{o}}{\theta },\end{eqnarray} \tag{ 14 }$$

where we remind the reader that a₂ is the inverse of the magnification in the direction parallel to the main CC.

As expected, τ scales with Σ and is independent of the masses of the microlenses. If enough PBHs are present near the main CC, at some point their CCs start to overlap (or, similarly, the microcaustics overlap in the source plane). This corresponds to τ = 1. From Equation (14), we can infer the width of the saturation region around the main CC. Since Equation (5) gives the total magnification and each train of micro-images carries half the magnification, after accounting for this factor of 2 the width of the saturation region is simply

$$\begin{eqnarray}&&{\rm{\Delta }}{\rm{\Theta }}=(4.2\times {10}^{-4}){\rm{\Sigma }}({M}_{\odot }\,{\mathrm{pc}}^{-2}){a}_{2}{\mu }_{o}.\end{eqnarray} \tag{ 15 }$$

This simple approximation works remarkably well when tested with the simulated data (see Figure 9, where the length of the yellow lines is derived from Equation (15) for the model of D16, a₂μ_o ≈ 0.2 × 150'' = 30''). Despite the assumption of similar masses made to derive Equation (14), we find that when the optical depth is computed exactly for a realistic distribution of masses (by integrating the areas within the Einstein ring for each microlens), the disagreement is approximately only a factor of 2 or 3. In particular, Equation (14) predicts τ = 0.7 at the position of Icarus (assuming the model of D16 with a₂μ_o ≈ 30'') and for ${\rm{\Sigma }}\approx 7\,{M}_{\odot }\,{\mathrm{pc}}^{-2}$, while the integration of a realization of the Spera model between $M=0.01\,{M}_{\odot }$ and the maximum mass in the realization (M ≈ 70 M_⊙) renders τ = 0.3. For alternative models like those of Woosley et al. (2002) and Fryer et al. (2012), we obtain similar values of 0.26 and 0.28, respectively. From Equation (15) and for the model of D16, we obtain ΔΘ ≈ 01.

Interestingly, this distance is smaller than (but comparable to) the separation found between Icarus and Iapyx (ΔΘ ≈ 026). This suggests that the adopted Σ may actually be close to its real value. Much higher values of Σ would result in a larger disruption of the main CC and would make it very difficult to observe a nearly constant flux for a period of ∼10 yr at the position of Icarus. Conversely, significantly smaller values of Σ would translate into very small probabilities of observing a peak like that witnessed for Icarus in 2016 May.

We can define the effective surface mass density of microlenses ${{\rm{\Sigma }}}_{\mathrm{eff}}={\rm{\Sigma }}{a}_{2}{\mu }_{o}/\theta$ (see Equation (14)), which is inversely proportional to the distance to the CC. Hence, at large separations from the CC, both Σ_eff and τ tend to zero and the observed light curve should be featureless (i.e., no microlensing events and just a slowly varying flux). On the other hand, a given area in the image plane at ΔΘ/2 from the main CC maps into a larger area in the source plane than the same area in the image plane would if it were at a smaller distance from the main CC (i.e., smaller θ). Thus, one would naively expect the probability of seeing microlensing events to be maximized at the point where we reach the saturation regime. Future observations of similar events should be found at comparable distances, or that scale with the estimated surface mass density of microlenses as shown in Equation (15).

This argument can be used to indirectly infer that the fraction of DM in the form of compact objects (like PBHs) must be low (independently of their masses). The value of Σ ≈ 7 M_⊙ pc⁻² corresponds to a fraction of the total mass of ∼0.3%. The relatively smooth behavior of the Icarus light curve (1–3 minor peaks in the last ∼10 yr) suggests that the fraction of DM in microlenses cannot be significantly higher than 3% (or the mass of the PBH is orders of magnitude smaller than 1 M_⊙). If the fraction of DM were much higher, the width of the saturation region would extend to farther distances. According to Equation (14), a 3% fraction would saturate a region of ∼1'' around the CC. Current observatories have the resolution and sensitivity to witness such events but none has been reported prior to Icarus/Iapyx. If we assume that the spiral arm extending all the way to the CC in MACS1149 contains a nearly constant surface density of stars (as suggested by its nearly constant surface brightness), during the last ∼10 yr we would have expected to see additional peaks farther away from the CC. The fact that both Icarus and Iapyx (but also Perdix, also called LS1/Lev 2017A; see K18 for more details) appear at distances ∼01 from the CC points to a low fraction of DM (<1%) in the form of microlenses. This conclusion is, in principle, independent of the microlens mass. However, this is not entirely true at sufficiently low microlens masses since in this case the microlenses start to behave as a smooth distribution of DM (their associated Einstein rings would be too small compared with the dimension of the macro-images). Even if individual events cannot be resolved, deep observations around the CC can reveal fluctuations that are not observed at larger distances from the CC. The strength and extension of these fluctuations can be used to constrain Σ_eff and hence the fraction of compact objects.

In the low-mass microlens regime, the sensitivity to the mass of the microlens depends on two variables: (i) the surface mass density of massive microlenses, Σ, and (ii) the underlying magnification. If massive microlenses (like, for instance, stars and remnants from the ICL) have a large Σ, small microlenses play a negligible role (the magnification pattern is entirely dominated by the massive ones). Only when the contribution from massive microlenses is small, can the signature of small microlenses (like PBHs with planet-like masses) be detected. As shown in Section 2, the Einstein ring associated with a small microlens can be highly amplified when it lies very close to a cluster CC. In the source plane, larger microcaustics from massive microlenses would overlap on top of the smaller microcaustics, overwhelming them. A clear illustration of this effect was shown in the right panel of Figure 4, where M = 0.01 M_⊙ microlenses can disrupt the main CC in a way that could be quantified during a caustic crossing event. The ideal scenario for constraining compact DM in the low-mass regime is to monitor a background galaxy that lies behind a caustic that is itself far away from any member galaxy or with minimal ICL. This will minimize the effect of more massive microlenses (like stars) offering a clean view of the structure of DM in the small-mass regime. The same argument can be applied to PBHs with higher masses. If a caustic crossing event is found in a region where the contribution from the ICL is known to be negligible and the event shows signatures of microlensing, PBHs could be constrained more easily than in a region where the ICL produces microlensing events.

Giant arcs are ideal candidates, since the exact location of the CC can be estimated with great accuracy owing to symmetry principles. Knowing the exact location of the CC is very important to break the degeneracy between the microlens mass and the magnification (see Equation (9)). In particular, elongated thin blue arcs tend to form in regions where κ ≈ γ ≈ 0.5, since then a₁ = 1 − κ − γ ≈ 0 and a₂ = 1 − κ + γ ≈ 1. The first (a₁ ≈ 0) condition produces arcs with very large tangential magnification while the second condition (a₂ ≈ 1) is needed to produce very thin arcs (small radial magnification). Thin giant arcs are normally found in elliptical potentials in regions relatively far away from the central BCG and where the ICL is moderate or small, reducing the negative impact of microlenses from the ICL that could be mistaken for PBHs. Blue arcs are also more likely to contain younger giant stars that can be very luminous and make ideal background sources.

From the observational point of view, having data with a high cadence (one data point every week or two weeks and daily when the flux aproaches the maximum) is important for constraining the shape of the peaks. The width of the peak depends most strongly on the ratio R/v, where R is the radius of the background source and v the relative velocity between the source and the caustic. Deep observations of the event are also useful for improving the photometry and ruling out (or confirming) smaller peaks that could be hidden between the most prominent ones. Also, if the background source consists of a binary system, different events would have a similar pattern with two consecutive peaks.

7.1. PDF of the Light Curves

The PDF of the observed magnification (extracted from light curves spanning ∼400 yr) is shown in Figure 14. Clearly, a compensating effect is taking place when F is higher. The PDFs look remarkably similar independent of the value of F. The situation is very similar to the result presented in Figure 8. When F is higher, the PDF shown in Figure 14 shifts toward lower magnifications. However, in the observed light curve, moderate magnifications (μ ≈ 10³) have similar probabilities. This is a consequence of the trains breaking up into smaller bits as F grows. Perhaps the most interesting difference is the trend observed at very low magnifications (μ ≈ 10) and very high magnifications (μ ≈ 10⁴). When F is high, it is more difficult to hide an entire macro-image from the observer (i.e., the total flux of the macro-image is very low) on the side where a₁ < 0. This is a consequence of the overlapping caustics in the source plane, as we have seen earlier. For high magnifications, Figure 14 also shows a deficit of bright peaks when F is higher. When F is higher, the side with negative parity (a₁ < 0) and the side with positive parity (a₁ > 0) behave similarly, as shown in Figure 14. The fact that the observed light curves from Icarus and Iapyx look very distinctive suggests indirectly that F must be ≲1%.

Zoom In Zoom Out Reset image size

7.2. Power Spectrum

The distribution of peaks and valleys depends on the optical depth of the microlenses. A correlation function, or power spectrum, of the light curves can highlight hidden features in the light curve such as correlations that are not easily identified in the PDF. Correlations between peaks are expected, for instance, when a single microlens produces a double peak. The separation between the peaks is, in this case, generally proportional to the square root of the mass of the microlens. We compute the power spectrum of the light curves and show the result in Figure 15 for the models with only ICL stars (Spera) and for models where we add to the Spera model a fraction of the total mass in the form of compact DM (or PBHs).

Zoom In Zoom Out Reset image size

The power spectrum seems to be more sensitive to the optical depth of microlenses on the side where a₁ < 0. This behavior is also found in the PDF of the light curves. Similarly to Schechter & Wambsganss (2002), we find that fluctuations are accentuated when the fraction of microlenses relative to the smooth component is smaller. Interestingly, on the side with a₁ < 0 and for the cases with higher F (PBH = 3.3% and PBH = 10%), the power spectrum behaves like pink noise (or 1/f noise) on large scales, a behavior which is found in many situations in nature. One of these situations is relevant for this work and refers to self-organized criticality²⁵ systems introduced by Bak et al. (1987; see Jensen 1998; Christensen & Moloney 2005, for more detailed descriptions of these systems).

The 1/f nature of the power spectrum on large timescales (more than a few months) indicates that there is no characteristic timescale on the side where a₁ < 0.

A similar argument could have been used for the PDF of the magnification or the light curves, which also exhibit power-law behaviors typical of systems where self-organization is involved. In this case, it is the magnification itself that is redistributed across the image plane, owing to the constraint that the total magnification (or observed flux) remains constant when integrating over long periods of time independent of the substructure in the lens plane. When F is large (PBH = 10%), the properties of the power spectrum resemble those of time series that fall in the categories of stationary fractals (Rodríguez 2014), which can be linked to renewal processes (a generalization of a Poissonian process) in which the time intervals between events are not correlated. This is the expected behavior for the separation between events when the distribution of projected microlenses is random.

The energy spectral density (ESD) can be used to measure the total amount of signal contained in the power spectrum,

$$\begin{eqnarray}&&E=\int | \mu (t){| }^{2}{dt},\end{eqnarray} \tag{ 16 }$$

which by Parseval's theorem is equivalent to the integral of the square of the power spectrum. Figure 16 shows the ESD for the Icarus and Iapyx simulated light curves and for the models shown in Figure 15. A clearer trend is observed on the side where a₁ < 0, suggesting that ESD of the fluctuations on this side of the main CC may be more capable of discriminating between different optical depths of microlenses.

Zoom In Zoom Out Reset image size

To produce accurate power spectra that can discriminate among different models it will be necessary to monitor these events for decades, similarly to what is being done for QSOs and other variable objects. With smaller samples, the two-point correlation function offers an interesting alternative.

Earlier work suggests that, from a statistical point of view, it is very difficult to distinguish the mass distribution of the microlenses for a fixed surface mass density and shear (Wyithe & Turner 2001). In contrast, Schechter et al. (2004) show how the PDF of magnifications is sensitive to the mass distribution of microlenses. Our results partially support this. Some differences exist for similar optical depths, but the differences are relatively small. Interestingly, these differences seem to be more accentuated on the inner part of the main CC (a₁ < 0).

One of the most striking features in the light curves, and which is also sensitive to the optical depth of microlenses, is the "hiding" periods when the flux of a macro-image falls below the detection limit. These low-magnification periods observed in our simulated light curves are similar to those found by Kayser et al. (1986) or Paczyński (1986) when the dimensionless surface mass density of stars, σ = Σ_s/(1 − Σ_c), is negative (here Σ_s is the surface mass density of stars and Σ_c is the surface mass density of the smooth distribution of matter). The counterintuitive results of Schechter & Wambsganss (2002) can be also interpreted in light of our own results. They find that "contrary to naive expectation, diluting the stellar component of the lensing galaxy in a highly magnified system with smoothly distributed DM increases rather than decreases the microlensing fluctuations caused by the remaining stars." Our simulations show that increasing the fraction of mass in compact objects results in lowering the typical magnification and consequently the fluctuations. Also, overlapping caustics smooth out large fluctuations.

Although in our simulations, σ is always positive (Σ_s > 0 and Σ_c < 1), making the masses of the microlenses negative in the simulation changes the magnification pattern of the a₁ > 0 side into a pattern similar to the a₁ < 0 side and vice versa. That is, it would be equivalent to the cases discussed by Kayser et al. (1986) and Paczyński (1986) with σ < 0. When comparing the histograms of intensities (or magnifications) in Paczyński (1986) with our results, we find a similar trend. Negative values of σ (or negative optical depth following Paczyński's argument) reproduce our results for a₁ < 0.

The work of Paczyński (1986) also presents another interesting result relative to the extension of the train of micro-images, or the halo in Paczyński's terms. As σ increases, the size of the halo of micro-images grows as well. In our case, we observe the same phenomenon with the trains extending over larger distances along the direction where micro-images form (i.e., in the direction of the deflection field). When Σ_s is sufficiently large, micro-images start to form in a more complex pattern, breaking the straight line found when Σ_s is low and adopting a configuration that more resembles the halo described by Paczyński (1986).

Perhaps the most interesting practical application of caustic crossing events is their ability to constrain the fraction of compact DM. For a given surface mass density of microlenses, Σ, the optical depth scales linearly with Σ. One can use this fact to set limits on Σ by observing the slow rise (or decline) in flux of a magnified background star as it approaches (or departs from) the cluster caustic. If Σ is sufficiently large and the position of the caustic or CC can be determined with relative precision, microlensing events should become ubiquitous thousands of years before (or after) the star crosses the cluster caustic. This can be visualized easily in the source plane for the case of high Σ.

The rate of observed events like Icarus at the position of background arcs known to cross caustics gives us indirect information about the level of disruption of cluster caustics and/or about the luminosity function of the background object. This is easier to see in the ideal case where the caustic is not disrupted. If the mass–luminosity relation follows analytical models (L ∝ M³ in normal stars where gas pressure and gravity are balanced, and L ∝ M in heavy stars where radiation pressure overwhelms gas pressure), we should expect a classic IMF that falls like ${dN}/{dM}\propto {M}^{-2.3}$ to translate into luminosity functions ${dN}/{dL}\propto {L}^{\alpha }$ with −2.3 < α < −1.1 depending on the mass (or luminosity) of the star. As discussed by K18, owing to the 1/μ² probability of being magnified by a factor larger than μ, if ${dN}/{dL}\propto {L}^{-2}$, the smaller probability of being magnified by a larger factor μ gets exactly compensated by the larger abundance of lower-luminosity stars, making all unlensed luminosities have the same probability of being observed. A stellar luminosity function, dN/dL, that is steeper than L⁻² results in lower-luminosity stars being more likely to be observed above a certain flux limit in caustic crossing events than higher-luminosity stars. In contrast, if dN/dL is shallower than L⁻², brighter stars will be more likely to be observed than the less luminous ones. When the cluster caustic is disrupted by microlenses, the argument above remains the same, since the probability of being magnified by a factor larger than μ still retains its fundamental form at high magnifications, μ⁻² (see Figure 14), except that now the normalization of the probability is smaller and very rare extremely luminous stars may be required to produce the observed flux.

8.1. Uncertainties in the Lens Model

Despite the high quality and excellent agreement between the different lens models that have been published for clusters like MACS1149, there are still significant uncertainties in the lens models that limit the capability of using CC crossing events as probes of DM (or substructure in general). For the work presented here, one of the most relevant systematics is the uncertainty in the magnification near the CC. Meneghetti et al. (2017) and Priewe et al. (2017) show how discrepancies of order 50% in the magnification are typically found between state-of-the-art lens models that are otherwise considered accurate. Here we compare the predictions made by two of these models derived under very different assumptions.

The first one is the free-form model from D16 derived using the WSLAP+ free-form code (Diego et al. 2005, 2007). The second is from Kawamata et al. (2016) derived using the parametric code GLAFIC (Oguri 2010). The two models correctly predict the position of the CC between Icarus and Iapyx, and have successfully anticipated the position and time of reappearance of SN Refsdal. However, as shown in Figure 17, the two models disagree by a factor ∼2 in the predicted magnification at Icarus and Iapyx. (A third model, by Richard et al. 2014, predicts a magnification similar to the D16 model with μ_o = 150 ± 10). The reason for this disagreement can be found in the small differences in the deflection field near the critical region. The condition for a CC to occur is that a certain combination of derivatives of the deflection field cancel out. In particular, the inverse of the magnification is given by

$$\begin{eqnarray}&&{\mu }^{-1}=1-{\alpha }_{x}^{x}-{\alpha }_{y}^{y}+{\alpha }_{x}^{x}{\alpha }_{y}^{y}-{\alpha }_{x}^{y}{\alpha }_{y}^{x},\end{eqnarray} \tag{ 17 }$$

where α_x and α_y are the deflection field in the x and y directions and ${\alpha }_{i}^{j}$ is the partial derivative of α_i with respect to the coordinate j ($i,j=x,y$). From the expression above, it is clear that small changes in the deflection field near the point where the inverse of the magnification is zero can have a significant impact on the delicate balance between the derivatives in Equation (17). The magnification is then more sensitive to the small differences in α near the CCs, where in relative terms these differences are larger.

Zoom In Zoom Out Reset image size

Small differences in the inferred deflection field between different lens models are expected. For instance, each method assumes different profiles for the halos around member galaxies or for the cluster halo(s), which results in potentials that are more or less shallow.²⁶ Also, the inferred positions of the background sources used to constrain the lens models differ from model to model. These unknown positions are treated as variables that need to be determined together with the mass distribution (or potential) of the galaxy cluster. By comparing the models of D16 and Kawamata et al. (2016), we find that CC positions, κ, and γ generally agree to within percent levels (except in the surroundings of small member galaxies considered in one model but not in the other). However, at the position of Icarus and Iapyx, the deflection angle disagrees by a factor of ∼40%. In absolute terms, the difference between deflection angles is small (only ∼3''), but since both deflection angles are relatively small in this part of the lens plane, it translates into a large relative difference. The ∼3'' difference in the deflection field originates in the different predicted positions of the background source used to constrain this part of the lens plane. The small differences in the derivatives of the deflection field around the CC translate into a factor of ∼2 in the magnification.

The uncertainty in the magnification has direct implications on our ability to constrain the masses of the microlenses. As shown in Section 2, the effective lensing mass of the microlens scales as $\sqrt{\mu }$. A factor of 2 uncertainty in the magnification then translates into a factor of $\sqrt{2}$ uncertainty in the mass of the microlenses. Future observations may need to rely on independent calibrators of the magnification. One such calibrator could be an SN Ia in a background galaxy and that happens to be near the CC. For instance, in Rodney et al. (2015), SN HFF14Tom is used to compare the observed magnification with the one predicted by several lens models, finding in general good consistency but also a small systematic bias in many of the lens models that tend to overpredict the observed magnification by ∼10%–20%.

8.2. Exploring the Small-mass Compact DM Regime

The best constraints on the compact DM when Σ is small should be obtained in a portion of the CC where the contribution from ICL stars (and remnants) is negligible. If ICL stars are populating the lens plane at the position of the CC, they will limit the ability to constrain the amount of compact DM. The ideal situation is to monitor a blue arc (to maximize the number of luminous stars) known to cross a CC and located far from the ICL or member galaxies. Ideally, such an arc contains a very bright star that is, on average, visible when magnified with the macromodel magnification in order to construct its light curve to probe DM (in contrast, it will be difficult to distinguish if different events correspond to different stars). In such situations, the possibilities of detecting small-size microlenses increase significantly. Moreover, as discussed earlier in the paper, monitoring arcs behind a galaxy cluster CC has the added advantage of probing a wide area in the lens plane for a fixed amount of time. This is a consequence of the increase in density of microcaustics at high magnifications. The reader can find an interesting discussion on this subject in Oguri et al. (2018).

If the contribution from ICL stars is very small, the light curve should contain very few peaks and the main caustic would be minimally disrupted, allowing the approach to the caustic to have maximal magnification. As shown earlier in Figure 13, at optical depths of microlensing which are a factor of 30 smaller than the value found at Icarus/Iapyx (i.e., Σ = Σ_o/30), the fluctuations in the light curve may still be significant, and the microlenses from the ICL may overwhelm the signature of an unknown population of microlenses with low mass (for instance, PBHs with masses smaller than 1 M_⊙). However, the results shown in the right panel of Figure 13 may not be the most realistic representation of a population of ICL microlenses with a small optical depth from the ICL. Estimating the precise amount of ICL that minimally disturbs the CC (needed to explore the low-mass regime of a hidden population of microlenses) is not trivial, but we can make some approximations.

First, and for the sake of simplicity, we assume that all the stars in the ICL have 1 M_⊙. (This is not a bad approximation, since we want to estimate the number density of objects which, for most IMF models, is dominated by stars with masses close to 1 M_⊙). Then, we make use of the scaling for the effective Einstein radius of the microlens, ${\theta }_{e}\sqrt{{\mu }_{t}}$ (Q = 1) (see Appendix A). Since from Equation (5), ${\mu }_{t}=\mu /{\mu }_{r}={\mu }_{o}/({\mu }_{r}\theta )$, we can find the value of θ for which the effective Einstein radius of the microlens is equal to θ, the distance to the CC. We can require that the separation between microlenses be at least twice this angle in order for the CC not to be completely disrupted.

For the values in the D16 model, we then infer that the surface mass density of stars must be $\sim {10}^{6}{M}_{\odot }/{\mathrm{arcsec}}^{2}$, or roughly two orders of magnitude smaller than the value found at the position of Icarus and Iapyx. Even though this may seem like a big difference, we have to realize that Icarus and Iapyx are found relatively close to the cluster center (∼7'' from the BCG), where the contribution from the ICL is still important. A simple fit to the total observed light around the BCG of MACS1149 (in the F160W band) shows how the unresolved flux around the BCG falls as $1/{d}^{2}$ (up to ∼30'' from the BCG), where d is the distance to the BCG. This is also the expected behavior if the ICL traces the DM, which from simulations (but also observations) falls as ρ(r) ∝ r⁻³ in the outer regions of a cluster (with $\rho (r)$ the three-dimensional mass density). Hence, in order for this unresolved flux to fall two orders of magnitude with respect to the ICL found at 7'', one should move a factor of ∼10 (that is, ∼1') from the BCG (less if the ICL falls faster than the assumed $1/{d}^{2}$ at larger radii, since there are fewer galaxies that can lose stars to the ICL). The CC of clusters at z ≈ 1.5 does not reach these distances, but at higher redshifts (z > 3) it can extend up to ∼1' along the main axis of the cluster. However, at these redshifts, the flux of a background star is a factor ∼5.5 times smaller than at z = 1.5 (or ∼1.8 mag fainter for z = 3). One would need extremely luminous stars at high redshift and a powerful telescope capable of observing their redshifted emission. Luckily, such combinations will exist in the near future. The James Webb Space Telescope (JWST) may be able to observe the first Population III stars through caustic crossing events (see Windhorst et al. 2018, where this scenario is studied in detail). Population III stars crossing a caustic far from the projected center of the cluster (to avoid the contamination from ICL microlenses) can then be used to constrain the fraction of compact DM to unprecedented levels.

This scenario offers the possibility of constraining PBHs in the less explored low-mass regime (M < 0.1 M_⊙). If we consider the ideal situation of a smooth deflection field with negligible contribution from the ICL stars, a microlens with small mass, M, close enough to the CC will behave as a larger microlens with mass μM, and will produce a peak before the maximum flux is reached. If the background star is large, this peak may merge with the primary peak, making the identification of the microlens harder (see Figure 18). However, if the background star is sufficiently small (few tens of solar radii), the peaks may be distinguished and the relative height of the peak due to the microlens can be used to constrain the mass of the microlens. Figure 18 shows how a microlens with a Jupiter mass at z = 0.55 can produce a change in flux of almost half a magnitude for a period of ∼1 day.

Zoom In Zoom Out Reset image size

Finally, as discussed in the previous subsection, the ideal target area would contain a magnification calibrator such as an SN Ia near a CC, which can be used to directly estimate the magnification in that part of the lens plane and reduce the uncertainty in the lens model magnification. SN Refsdal is relatively close to Icarus, but it is also too close to a medium-size member galaxy which overwhelms the cluster contribution.

This work studies the particular case of microlenses very close to a cluster CC. As one approaches the CC, the magnification changes rapidly, affecting the way microlenses disrupt the magnification pattern. The main results are summarized below.

1. Superluminal substructure probe. We have shown how the observed velocity of the macro-images is proportional to the magnification. Consequently, the apparent transverse motion of the observed macro-images becomes superluminal when the distance to the CC is sufficiently small and the magnification is sufficiently large. Under these circumstances, moving macro-images probe more substructure in the lens plane than a classic microlensing event.

2. Distinctive light curves. As found in previous work, the light curves on each side of the CC look quite different (at low optical depths and magnifications smaller than 1000), with periods of very low flux on the side where a₁ < 0 not present on the a₁ > 0 side. This deficit in flux on the side with a₁ < 0 is compensated by brighter peaks when a micro-CC is intersected by the background star. At high optical depth (τ > 1) or extreme magnifications (from the macromodel) of several thousands, there is little distinction between the positive and negative parity sides.

3. Reduced maximum magnification. Macro-images of a small background source, like a bright star, crossing a cluster CC can be magnified by factors of ∼10⁶, provided there is no substructure that disrupts the CC. When microlenses populate the lens plane, the maximum magnification is reduced, but it can still reach values of several thousands.

4. Flux borrowing. Microlenses can borrow photons from the main caustic thousands of years before the background star crosses the position of the main caustic, making the observation of these events much more likely than previously thought. This produces peaks in the light curve with magnification factors of several thousands. If microlenses are ubiquitous, multiple peaks are expected hundreds or thousands of years before the background star crosses the main caustic. The number of peaks (and their intensity) depends on the density of microlenses and their masses. The constraint on flux conservation implies that, during the crossing of the main caustic, the magnification is orders of magnitude smaller than for the smooth DM model when microlenses are disrupting the main caustic. Moderate disruption of cluster caustics increases the probability of detecting microlensing events of fainter background stars, while significant disruption produces microlensing events much farther away from the main CC that can be observed only when the background star is extremely luminous.

5. Saturation and 1/f noise. When microlenses populate the lens plane, there is a region around the main CC where the optical depth exceeds the critical value of 1. The width of this region depends on the surface mass density of microlenses and the properties of the cluster deflection field. Within this region, the microcaustics overlap in the source plane and the magnification pattern is changed substantially. The power spectrum of the light curves shows features on larger timescales that resemble the ubiquitous 1/f noise present in self-organized criticality or renewal processes studied in other fields.

6. Constraints on PBHs. The condition required for the interior macro-image (a₁ < 0) to disappear during an extended period is that there exist low-magnification areas in the source plane. This effect was also predicted in earlier work, and it was shown to be important for low contributions from microlenses (see, e.g., Schechter & Wambsganss 2002). When enough massive microlenses are added in the lens plane, for instance if PBHs make up a significant fraction of the DM (larger than a few percent), it becomes harder to satisfy this condition since the microcaustics start to overlap on top of the low-magnification regions. Thus, the lack of extended periods where one of the macro-images vanishes can be used to set a limit on the fraction of DM in the form of microlenses (like PBHs). So far, the data provided by Icarus seem too sparse to place strong constraints on the fraction of DM that can be in a compact form (although they appear to favor a low fraction of DM as PBHs). Continuous monitoring of both Icarus and Iapyx in the near future will allow us to derive such strong constraints after the identification of additional peaks and the precise modeling of these peaks together with the extended periods between peaks. Interestingly, K18 report a possible nearby third event (Perdix, or LS1/Lev 2017A) that could correspond to a different background star being lensed by the same web of caustics.

7. Complementary to QSO microlensing. An obvious advantage of this type of observation when compared with QSO microlensing is that, owing to the smaller intrinsic size of the background source (a few solar radii as opposed to ∼0.01 pc (∼0.001 pc) for optical (X-ray) accretion disk) and the larger apparent motion of the macro-images, the timescale of a microlensing event is much shorter (by 3 or 4 orders of magnitude) and the magnifications during maxima are larger (by 1.5 or 2 orders of magnitude). Light curves spanning several years can intersect several microlenses and can be used to determine a census of microlenses. The constraints derived in this way on the population of microlenses are then affected mostly by local substructure and less by projection effects acting on larger scales.

8. Cosmic microscopes. When the optical depth of microlenses is very small, a caustic crossing event can be used to probe very small masses near the CC. A microlens at cosmic distance with a mass similar to that of Jupiter can produce features in the light curve that can be observed with current technology provided the radius of the background star is sufficiently small. This opens the exciting possibility of probing compact DM in the low-mass regime, unreachable by other means, where at extreme magnifications of order 10⁶, a Jupiter-like microlens would behave as a microlens with a mass ∼1000 M_⊙. Monitoring (with JWST) of high-redshift galaxies rich in luminous Pop III stars that are intersecting a CC far away from the ICL would be valuable, since these are ideal targets for this kind of study (Windhorst et al. 2018).

The authors wish to thank Juan M. Lopez, Miguel Angel Rodriguez, Gary Bernstein, Evencio Mediavilla, and the anonymous referee for inspiring discussions and useful comments. J.M.D. acknowledges the support of projects AYA2015-64508-P (MINECO/FEDER, UE), AYA2012-39475-C02-01, and the consolider project CSD2010-00064 funded by the Ministerio de Economia y Competitividad. J.M.D. acknowledges the hospitality of the Physics Department at the University of Pennsylvania for hosting him during the preparation of this work. T.J.B. gratefully acknowledges the Visiting Research Professor Scheme at the University of Hong Kong. M.O. is supported in part by World Premier International Research Center Initiative (WPI Initiative), MEXT, Japan, and JSPS KAKENHI Grant Numbers 26800093 and 15H05892. M.J. acknowledges support by the Science and Technology Facilities Council (grant number ST/L00075X/1). Support for S.R. through HST programs GO-13386, GO-14199, and GO-14208 was provided by NASA through grants from the Space Telescope Science Institute (STScI), which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. A.V.F. and P.L.K. were supported by NASA/HST grants GO-14041, GO-14199, GO-14208, GO-14528, GO-14872, and GO-14922 from STScI; they are also grateful for financial assistance from the Christopher R. Redlich Fund, the TABASGO Foundation, and the Miller Institute for Basic Research in Science (U.C. Berkeley).

This appendix presents the basic formalism for a single microlens near a CC. We consider the simple scenario where the deflection field from the cluster is oriented in the vertical direction (y axis in the equations below). In this configuration, the CC from the cluster will be a horizontal line, the tangential magnification will be that in the vertical direction, and the radial magnification will be in the horizontal direction. A reader familiar with lensing models may find this definition of tangential and radial magnifications counterintuitive, since we are assigning the radial magnification to a direction that is tangential with respect to the source of the deflection field. This, however, is the case for some tangential arcs that stretch in directions perpendicular to the CC, like for instance the arc containing the Icarus event discussed by K18.

If we place a microlens (point source) with Einstein radius r_e at a position (x_o, y_o), the deflection field around the point source is ${\alpha }_{x}={{xr}}_{e}^{2}/{r}^{2}$ and ${\alpha }_{y}={{yr}}_{e}^{2}/{r}^{2}$, where ${r}^{2}={x}^{2}+{y}^{2}$ and locations (x, y) have their origin at $({x}_{o},{y}_{o})$. The mapping between the image plane and the source plane is no longer given by Equation (4) (in the main text); instead, it can be described by

$$\begin{eqnarray}&&\beta =\left(\begin{array}{c}(x+{x}_{o})/{\mu }_{r}\\ {(y+{y}_{o})}^{2}/{r}_{{\rm{E}}}\\ \end{array}\right)-{r}_{e}^{2}\left(\begin{array}{c}x/{r}^{2}\\ y/{r}^{2}\end{array}\right),\end{eqnarray} \tag{ 18 }$$

where the first term accounts for the deflection field from the cluster (${r}_{{\rm{E}}}={\rm{\Theta }}$ and y = θ in Equation (4)) and the second term corresponds to the deflection field from the point source. We have assumed that the deflection field of the cluster is such that it stretches the images in the y direction by some large factor ${\mu }_{t}=1/{a}_{1}$ and in the x direction by a smaller factor ${\mu }_{r}=1/{a}_{2}$ (i.e., μ = μ_tμ_r). The constant μ_r is the small eigenvalue (for a tangential CC) of the magnification (i.e., ${a}_{2}={\mu }_{r}^{-1}\,=1-\kappa +\gamma$). Similarly, ${a}_{1}={\mu }_{t}^{-1}=1-\kappa -\gamma$. The inverse of the magnification, μ⁻¹, can be computed as the determinant of the Jacobian between the source-plane positions and the image-plane positions:

$$\begin{eqnarray}\displaystyle \frac{d\beta }{d\theta }=\left(\begin{array}{cc}{a}_{2}+({x}^{2}-{y}^{2})\displaystyle \frac{{r}_{e}^{2}}{{r}^{4}} & 2{xy}\displaystyle \frac{{r}_{e}^{2}}{{r}^{4}}\\ 2{xy}\displaystyle \frac{{r}_{e}^{2}}{{r}^{4}} & 2\displaystyle \frac{(y+{y}_{o})}{{r}_{{\rm{E}}}}-({x}^{2}-{y}^{2})\displaystyle \frac{{r}_{e}^{2}}{{r}^{4}}\end{array}\right).\end{eqnarray} \tag{ 19 }$$

Computing the determinant of the expression above results in

$$\begin{eqnarray}&&{\mu }^{-1}={a}_{1}{a}_{2}-({a}_{2}-{a}_{1})\displaystyle \frac{{x}^{2}-{y}^{2}}{{x}^{2}+{y}^{2}}\displaystyle \frac{{r}_{e}^{2}}{{r}^{2}}-\displaystyle \frac{{r}_{e}^{4}}{{r}^{4}},\end{eqnarray} \tag{ 20 }$$

where ${a}_{1}=2(y+{y}_{o})/{r}_{{\rm{E}}}$. When there is no microlens (r_e = 0), μ⁻¹ is equal to the original ${\mu }^{-1}={\mu }_{t}^{-1}{\mu }_{r}^{-1}={a}_{1}{a}_{2}$. The term $({x}^{2}-{y}^{2})/({x}^{2}+{y}^{2})$ is equal to $\cos (2\delta )=Q$, where $\tan (\delta )\,=y/x$. Q changes sign in different quadrants as shown in Figure 19.

Zoom In Zoom Out Reset image size

Approximate solutions to Equation (20) can be found after a few simple assumptions are made. Analytical solutions, even though approximate, are very useful for unveiling scalings with the microlens mass or distance to the CC of the cluster that can be exploited later in statistical analyses. To find an analytical solution of Equation (20), we assume that ${a}_{1}=2(y+{y}_{o})/{r}_{{\rm{E}}}\approx 2{y}_{o}/{r}_{{\rm{E}}}=\mathrm{constant}$. This approximation is valid when y_o ≫ y, which corresponds to a distance from the cluster CC much larger than the Einstein radius of the microlens. A simple solution (CCs) can be found after setting this equation to zero, multiplying by r⁴, and making the variable change $x={r}^{2}/{r}_{e}^{2}$. This results in the following general solution,

$$\begin{eqnarray}&&{r}^{2}={r}_{e}^{2}\displaystyle \frac{({a}_{2}-{a}_{1})Q}{2{a}_{1}{a}_{2}}\left[1\pm \sqrt{1+\displaystyle \frac{4{a}_{1}{a}_{2}}{{({a}_{2}-{a}_{1})}^{2}{Q}^{2}}}\ \right].\end{eqnarray} \tag{ 21 }$$

Depending on the values of a₁, a₂, and Q, the quantity r² can be positive or negative. If r² is negative there is no real solution. On each side of the CC, ${a}_{1}=1-\kappa -\gamma$ changes sign (and equals zero on the CC). The solutions of Equation (20) can be divided into two sets, one having a₁ > 0 and the second with a₁ < 0. Below we discuss these solutions based on simple approximations.

For simplicity we assume that ${\mu }_{r}$ is constant and positive. Near a tangential CC, μ_r (or similarly a₂) changes very slowly around a given point, so this is a good approximation for our purposes. Also, in our case (microlenses near a tangential CC), ${\mu }_{t}\gg {\mu }_{r}$ (or a₂ ≫ a₁) and the term (a₂−a₁) is always positive. For convenience we define

$$\begin{eqnarray}&&F=1\pm \sqrt{1+4{a}_{1}{a}_{2}/{(({a}_{2}-{a}_{1})Q)}^{2}},\end{eqnarray} \tag{ 22 }$$

and in particular,

$$\begin{eqnarray}&&{F}_{+}=1+\sqrt{1+\displaystyle \frac{4{a}_{1}{a}_{2}}{{({a}_{2}-{a}_{1})}^{2}{Q}^{2}}}\end{eqnarray} \tag{ 23 }$$

and

$$\begin{eqnarray}&&{F}_{-}=1-\sqrt{1+\displaystyle \frac{4{a}_{1}{a}_{2}}{{({a}_{2}-{a}_{1})}^{2}{Q}^{2}}}.\end{eqnarray} \tag{ 24 }$$

Since ${a}_{2}\gg {a}_{1}$, we can approximate ${F}_{+}\approx 2$ and ${F}_{-}\,\approx -2{a}_{1}/({a}_{2}{Q}^{2})$. The term $Q=\cos (2\delta )$ can be either positive or negative depending on the quadrant: Q < 0 in π/4 < δ < 3π/4 (hereafter quadrant 1) and 5π/4 < δ < 7π/4 (hereafter quadrant 3), and Q > 0 in −π/4 < δ < π/4 (hereafter quadrant 4) and 3π/4 < δ < 5π/4 (hereafter quadrant 2). For clarity we mark these quadrants in Figure 19.

A.1. Solutions for a₁ > 0

A.2. Solutions for a₁ < 0

This section includes various animations of a background star moving between the network of caustics. Figures 21 and 20 show how counterimages are being formed and destroyed as a very large background star moves with 1000 km s⁻¹ relative velocity toward the main caustic. The radius of the star is very large (70,000 R_⊙) in order to better visualize the counterimages that otherwise would be unresolved at this resolution. This large radius is comparable to the X-ray emitting region of a QSO so the movies are useful also to see how these QSO regions would form and disappear if they could be resolved. Figure 20 is for a case where the microlenses are the stars and remnants from the intracluster medium. Figure 21 adds to these microlenses PBH with 30 M_⊙. The surface mass density of the PBH makes up to 1% of the total surface mass density of the cluster. In both cases, the top panel shows counterimages in the side with positive parity and the bottom panel is for counterimages in the side with negative parity. Note how in the case where PBHs are present counterimages can form relatively far away from each other. Future observations with milliarcsecond or sub-milliarcsecond resolution can potentially resolve the different locations and set strong limits on the fraction of PBHs.