Exotic Image Formation in Strong Gravitational Lensing by Clusters of Galaxies – II: Uncertainties

Galaxy clusters as strong lenses offer a unique probe to the matter distribution as lensing does not discriminate between the visible and dark matter distribution (e.g., Blandford & Narayan, 1992; Kneib & Natarajan, 2011). At the same time, by magnifying distant sources (which would have otherwise remained unobserved), it also opens up the possibility of studying the high-redshift Universe. Observation of these cluster lenses led us to the detection of multiply imaged supernovae (Kelly et al., 2015, 2016), which can help us in constraining the Hubble constant. On the other hand, the observation of highly magnified stars (Kelly et al., 2018) can help us in probing the mass-function of intra-cluster medium (Diego et al., 2018) and observing the pop III stars (Windhorst et al., 2018). Due to the high magnification provided by lensing, the detection of high-redshift ($z>7$) galaxies (e.g., Bradley, et al., 2008; Coe, et al., 2013; Watson, et al., 2015) can help in determination of the galaxy luminosity function (e.g., Atek et al., 2018).

In order to use galaxy clusters as a probe, we need to model their mass distributions. However, as the data from the observation is limited, one cannot model the mass distribution of these clusters with arbitrary precision and resolution. As a result, to reconstruct the lens mass distribution, different groups start with different sets of prior. For examples, light trace mass (LTM, Broadhurst et al., 2005; Zitrin et al., 2009) parametric mass reconstruction assumes that the mass distribution in the cluster lenses follows the light distribution. On the other hand, the non-parametric (free-form, Liesenborgs et al., 2007) mass reconstruction methods do not rely on any preliminary information related to the mass model and only take into account the strong and weak (Liesenborgs et al., 2020) lensing data. Hybrid mass reconstruction methods take input from both parametric and non-parametric approaches (Sendra et al., 2014). Due to the finite amount of data and different set of priors, it is possible that these different methods give different results when applied on the same cluster lens (e.g., Smith et al., 2009; Zitrin & Broadhurst, 2009). Hence, it is very important to compare these different techniques in case of simulated as well as real lenses to improve these methods and for more robust predictions (e.g., Meneghetti et al., 2017; Priewe et al., 2017; Rodney et al., 2015; Remolina González, Sharon, & Mahler, 2018; Raney et al., 2020). Each of these methods results into the best-fit mass model parameters and uncertainties associated with them. One avenue to observe the effect of these uncertainties is the cluster lens magnification maps. For the best-fit lens mass model, one will have a certain area in the lens plane that gives higher magnification than a threshold value. However, once we account for the statistical uncertainties, the area with a magnification greater than the threshold area is not a unique number; instead, it will be denoted by a range. As a result, various observational estimates (like the number of highly magnified galaxies) are also subjected to these uncertainties.

For a source to be highly magnified, the source must lie close to the caustic in the source plane and the corresponding images are formed near the corresponding critical line in the lens plane. The magnification factor by which the source is magnified also depends on the source size: smaller the source, larger the magnification factor, and vice versa. This is why one can have a magnification factor of thousands for lensed stellar sources (Kelly et al., 2018), whereas the maximum magnification factor for galaxy sources can only be in hundreds. These caustics and critical lines are the singularities of the lens mapping and are of two different types: stable and unstable. Fold and cusp are stable singularities as they are present for all source redshift, whereas beak-to-beak, umbilics, and swallowtails are unstable singularities as they only occur for specific source redshifts for a given lens system. So far, the number of known image formation near unstable singularities is very small: one near hyperbolic umbilic (Orban de Xivry & Marshall, 2009; Zitrin et al., 2010), handful near swallowtail (e.g., Abdelsalam, Saha, & Williams, 1998; Suyu & Halkola, 2010), and none near elliptic umbilic. With the upcoming observing facilities like: Euclid: (Laureijs, 2009), James Webb Space Telescope: (JWST, Gardner, et al., 2006), Nancy Grace Roman Space Telescope: (WFIRST, Akeson, et al., 2019), Vera Rubin Observatory: (LSST, Ivezić, et al., 2019), the number of strongly lensed system is expected to increase by more than an order of magnitude, and it is likely that we will encounter image formations near these unstable singularities. Hence, it is timely to explore the properties of these unstable singularities in detail.

In our previous work, we have described the method to locate these unstable singularities for a given lens model (Meena & Bagla, 2020, hereafter MB20) and applied this method in the case of actual cluster lenses (Meena & Bagla, 2021, hereafter MB21). The final output of this method is a singularity map containing all the point singularities and $A_{3}$-lines. Since these $A_{3}$-lines correspond to cusps in the source plane, they mark the high-magnification regions in the lens where image formation near cusps will occur (three or more lensed images lying near to each other). The point singularities depend on the second and higher-order derivatives of the lens potential. Hence, they are very sensitive to the presence of small/intermediate-scale structures and the variations in the lens parameters. As a result, the singularity maps corresponding to the parametric and non-parametric mass models can be very different for the same cluster lens As pointed out in MB21, parametric mass models give significantly larger number of point singularities due to the small scale structures compared to non-parametric mass models.

In our current work, we study the effect of the statistical uncertainties associated with the reconstructed lens mass model parameters on the singularity maps and the higher-order singularity cross-section. In order to study the effect of mass model uncertainties on the singularity maps, we have considered two simulated galaxy cluster lenses, Irtysh I and II from Ghosh, Williams, & Liesenborgs (2020, hereafter GWL20) and all of the six galaxy cluster lenses from the Hubble Frontier Fields (HFF) survey. The lens mass reconstruction of all of these cluster lenses has been done using grale¹¹1https://research.edm.uhasselt.be/jori/grale2/ (Liesenborgs et al., 2007). We choose the grale mass models as the corresponding singularity maps are simplest and provide a lower limit on the point singularity cross-section (MB21). Further, as the grale mass model does not make any assumption about density profiles of substructure, the reconstruction is independent of the nature of dark matter. Due to the simple nature, we can construct singularity maps for a larger region of the lens plane without introducing spurious point singularities. Apart from that, different grale runs for a cluster lens give independent mass maps. Hence, there is no correlation between different mass maps reconstructed for a lens using grale that may not be true for the parametric mass models. For simulated clusters, Irtysh I and II, we construct the singularity maps for the original mass models, for the individual runs (with 150, 500, 1000 lensed images) and for the final best-fit mass model. For HFF clusters, we construct the singularity maps for the individual runs and for the best-fit mass model obtained by the averaging of the individual runs.

As we know that the $A_{3}$-lines in the singularity maps correspond to the cusps in the source plane. Hence, a singularity map also gives the information about the number of cusps formed in the source plane at a given source redshift. By drawing the critical lines for a given source redshift overlaid on the singularity map, one can calculate the number of cusps in the source plane by counting the points where a critical line and corresponding $A_{3}$-line cut each other. As the grale best-fit mass models provide the simplest singularity maps, they also give the lower limit on the cusp cross-section. Hence, it is worthwhile to calculate the lower limit of the cusp cross-section in the HFF clusters using the grale best-fit mass models. Following that, we also estimate the (tangential and radial) arc cross-section in case of both simulated and HFF cluster lenses. The source galaxy luminosity function has been taken from Cowley et al. (2018, hereafter C18). In C18, authors estimate the number of galaxies expected in the deep galaxy surveys with the JWST in different the filters for different exposure time. Although the expected galaxy population in any of the JWST filter can be used in our study, we only considered the expected galaxy population in one JWST filter, F200W, with an exposure time of $10^{4}$ seconds (please see C18 for more details).

This paper is organized as follows. In §2, we briefly revisit the basics of the singularities in the strong lensing. In §3, we present our results for simulated Irtysh clusters. The corresponding singularity maps are discussed in §3.1, followed by the discussion of the redshift distribution of singularities in §3.2. The results for the HFF clusters are presented in §4. The corresponding singularity maps, redshift distribution of singularities, and the arc cross-section are discussed in §4.1, §4.2, and §4.3, respectively. Summary and conclusions are presented in §5. The cosmological parameters used in this work to calculate the various quantities are: ${\rm H}_{0}=70\>{\rm kms}^{-1}{\rm Mpc}^{-1},\>\Omega_{\Lambda}=0.7,\>\Omega_{m}=0.3$.

In this section, we briefly review the basics of singularities in strong gravitational lensing to set the notions and notations. For a more detailed discussion about singularities, we encourage the reader to see into MB20, and Schneider, Ehlers, & Falco (1992).

For a given strong lens system, magnification factor in the lens (image) plane is given as,

where $\alpha$ and $\beta$ are the eigenvalues of the deformation tensor $\psi_{ij}$, with $\psi$ being the projected lens potential and subscript $ij$ represents its partial derivatives with respect to the image plane coordinates. $a=D_{\rm ds}/D_{\rm s}$ is the distance ratio, where $D_{\rm ds}$ and $D_{\rm s}$ represent the angular diameter distance from the lens to the source and from the observer to the source, respectively. The points in the lens plane where $\mu^{-1}=0$ are known as the singular points of the lens mapping (see MB20 for more details). The location of these singular points in the lens plane defines the critical curves and the corresponding points in the source plane are known as the caustics. The critical curves in the lens plane are smooth closed curves, and trace regions with high magnification, whereas their counterpart, caustic curves in the source are closed but not necessarily smooth. These caustics in the source plane are made of smooth segments (known as folds) meeting with each other at the cusp points.

For a given strong lens system, fold and cusp are the only stable singularity of the lens mapping as they are present for all possible source redshifts. The evolution of caustic structure in the source plane with source redshift mainly depicts the formation, destruction, or exchange of cusps between radial and tangential caustics. Other higher-order singularities (also known as point singularities) like swallowtail, purse, and pyramid are only present for specific source redshifts and are sensitive to the lens parameters. A small change in the lens system parameters can result in the complete disappearance of a point singularity. These point singularities, along with the $A_{3}$-lines constitute the so called singularity map. $A_{3}$-lines form the backbone of the singularity map and denote the points in the lens plane which correspond to the cusps in the source plane for all possible source redshifts. Different point singularities are associated with the formation, destruction, or exchange of cusp between radial and tangential caustics, hence, they also satisfy the $A_{3}$-line condition along with additional criteria. $A_{3}$-lines are given by the condition $n_{\lambda}.\nabla_{x}\lambda=0$, where $\lambda$ denotes the eigenvalue of the deformation tensor, and $n_{\lambda}$ is the corresponding eigenvector. The swallowtail singularity denotes the points in the lens plane where the eigenvector $n_{\lambda}$ is tangential to the corresponding $A_{3}$-line, and at a swallowtail singularity two extra cusps appear in the source plane. The hyperbolic and elliptic umbilics mark the points in the lens plane where both eigenvalues of the deformation tensor are equal to each other, i.e., zero shear points. We encourage the reader to see MB20 for more details about the point singularities and the corresponding characteristic image formations.

As these point singularities depend on the second and high-order derivatives of the lens potential, the effect of variation of lens parameters can be seen clearly in the corresponding singularity map. As a result, the singularity map is also very sensitive to the lens mass reconstruction method, as shown in MB21. In MB21, one can see that the cluster mass models reconstructed using non-parametric method grale give a significantly lower number of point singularities compared to the parametric mass models for the same cluster, hence, putting the lower limit on the point singularity cross-section. The singularity maps corresponding to the non-parametric mass models also show very simple $A_{3}$-line structures compared to the singularity map corresponding to the parametric mass models. Hence, the three image (tangential and radial) arc cross-section in parametric and non-parametric mass models is also expected to show significant differences.

In MB21, it has been shown that the parametric and non-parametric mass models for a given HFF cluster show a significant difference in both cusp and point singularity cross-section. However, each mass model has uncertainties associated with it due to the finite amount of observational data. As a result, these uncertainties also affect the results derived from the reconstructed mass models. In this section, we study the effect of these uncertainties on the singularity maps corresponding to reconstructed mass models for simulated galaxy clusters. We have considered two simulated clusters, Irtysh I and II from GWL20. GWL20 reconstructed free-form mass models using grale for both Irtysh I and II with three different sets of multiple images, 150, 500, 1000. The lens mass reconstruction with 150 images corresponds to the current observational scenario, whereas the 500 and 1000 multiple image cases correspond to future observations with the Hubble Space Telescope (HST) and James Webb Space Telescope (JWST). For each case, individual mass models are obtained from forty different and independent grale runs and the final best-fit mass model is an average of the mass distributions from these forty independent grale runs. While constrained by the required computational resources this number is consistent with previous works with grale and it sufficiently spans the range of uncertainties in the reconstructed mass models. We refer the reader to look into GWL20 for more details. Hereafter, for simplicity, the average mass models for 1000, 500, 150 image cases for Irtysh I/II will be referred as Irtysh IA/IIA, IB/IIB, IC/IIC.

Under the Hubble Frontier Fields (HFF) Survey²²2https://archive.stsci.edu/prepds/frontier/ program, Hubble Space Telescope observed a total of six massive merging clusters (see Lotz et al. 2017 for more details). For every HFF cluster, different groups have reconstructed mass models using parametric, non-parametric and hybrid (a combination of parametric and non-parametric) methods (e.g., Diego et al., 2007; Bradač et al., 2005; Jauzac et al., 2015; Johnson et al., 2014; Lam et al., 2014; Merten et al., 2011; Mohammed et al., 2016; Oguri, 2010; Sendra et al., 2014; Williams, Sebesta, & Liesenborgs, 2018; Zitrin et al., 2013). To study the statistical uncertainties on the singularity maps, we consider the non-parametric HFF cluster mass models. These mass models for HFF clusters are reconstructed using the free-form method grale. There are two reasons for choosing these mass models for the HFF clusters:
(1) The best-fit grale mass models for HFF clusters gives the simplest singularity maps compared to other techniques as shown in MB21 and puts the lower limit on the point singularity cross-section. Hence, it is worthwhile to check how statistical uncertainties affect the lower limit.
(2) In MB21, we only considered the central $40^{\prime\prime}\times 40^{\prime\prime}$ region for the analysis. This choice was made due to the increasing number of spurious singularities in the parametric mass models as we go away from the center of the cluster. Hence, it will not be useful to analyze the effect of uncertainties on the singularity map in the presence of artifacts. However, such a problem does not exist for these non-parametric mass models as the singularity map is very simple and a resolution of $<0.1^{\prime\prime}$ is sufficient enough as shown in MB21.

Similar to the Irtysh clusters, for every HFF cluster, grale provides forty individual mass maps and one best-fit mass map that is the average of these forty mass maps. Depending on the number of images used in the reconstruction, there are different versions of the reconstructed mass models available. Here, we are using the v4 mass models for all of the HFF clusters, and, for simplicity, we will be using the abbreviated names for the HFF cluster lenses.

In this work, we have studied the effect of statistical uncertainties associated with cluster lens mass reconstruction with grale on point singularities. We considered two simulated clusters, Irtysh I and II from GWL20 and six galaxy clusters from the HFF survey. For both simulated and real galaxy clusters, the mass models were reconstructed using grale. For each simulated galaxy cluster, we had three different mass models reconstructed with different sets of multiple images (please see GWL20 for more details). For each cluster lens, 40 reconstructions were done, and the best-fit mass models were the average of them. As a first step, we have constructed singularity maps for the original mass models (in the case of Irtysh I and II), for all of the individual mass models and for the best-fit mass models. In the next step, we have compared the singularity maps corresponding to the individual realizations with the original and average best-fit mass models. These comparisons for simulated and HFF clusters are shown in Figure 3 and 6, respectively.

We find that the singularity maps corresponding to the best-fit mass models have a relatively small number of point singularities compared to the individual runs. As a result, the best-fit mass models provide a lower bound on the number of singularities that a lens has to offer. For mass models reconstructed by grale, such a result is expected due to the fact that individual realizations contain spurious features, and the averaging to get the best-fit mass models is done in order to remove the contribution from these spurious features. Hence, the final best-fit mass model contains a low number of small scale structure and low number of point singularities. The differences in the number of point singularities between individual runs and the best-fit are maximum at source redshift $z<2$ for both simulated and real clusters. One can expect such discrepancy to arise at these redshifts as the mass models in the central regions of the clusters, at present, is not very well constrained. However, in order to be sure about this explanation, one will need to construct the mass models for simulated clusters with multiple images near the cluster center.

Such an averaging (over 40 realizations) to get the final best-fit mass model also simplifies the over $A_{3}$-line structures in the singularity maps. As $A_{3}$-lines in the singularity maps correspond to the number of cusps in the source plane at every source redshift, the final best-fit mass model also has a smaller number of cusps compared to individual realizations. Hence, the final best-fit mass models also provides a lower limit on the number of image formations near cusps. In our current work, we have also calculated the cusp cross-section for the HFF clusters using the best-fit mass models. we find that A370 is most efficient in producing three image arcs and the image formation near point singularities among the HFF clusters based on grale modeling. In the HFF clusters, one would expect to observe at least 10-20 tangential and 5-10 radial three image arcs with JWST.

The key take away from our current work is that the best-fit mass models constructed using the grale provide the lower limit on the number of cusps and point singularities for a given cluster lens. Folding in the statistical uncertainties does not decrease the numbers obtained from the best-fit mass models. As shown in MB21, one would expect to observe at least one hyperbolic umbilic and one swallowtail for every five clusters with JWST even when we include the corresponding statistical uncertainties. Apart from that, one also expects at least 10-20 tangential and 5-10 radial three image arcs with JWST in the HFF clusters.

The above results are based on the non-parametric modeling by grale, which underestimates the contribution from the galaxy scale substructures in the cluster lenses. Hence, if one considers parametric modeling, the above quoted numbers are expected to increase. As shown above, these singularity maps also provide a very nice way to compare the efficiency of different cluster lenses. This feature can also be used to compare the simulated and real clusters on qualitative (visual comparison of singularity maps) and quantitative basis (arc and point singularity cross-sections).

Apart from that, image formation corresponding to arcs or point singularities shows three or more images lie in a very small region of the lens plane. Hence, the time delay between these images is expected to be small compared to other (typical) multiple image formation scenarios. The other important feature of point singularities is their characteristic image formation, which can be very helpful in locating them. These possibilities are the subject of our ongoing work.

AKM would like to thank Council of Scientific $\&$ Industrial Research (CSIR) for financial support through research fellowship No. 524007. Authors thank the anonymous referee for useful comments. This research has made use of NASA’s Astrophysics Data System Bibliographic Service. We acknowledge the HPC@IISERM, used for some of the computations presented here.

The simulated cluster (Irtysh I and II) data from GWL20 and the high resolution HFF cluster mass models corresponding to the Williams group are available from the modelers upon reasonable request.