<figure><img src="image_1.png" /></figure>

THE EVOLUTION OF THE FAINT END OF THE UV LUMINOSITY FUNCTION DURING THE PEAK EPOCH OF STAR FORMATION ($1 < z < 3$)*

ANAHITA ALAVI¹, BRIAN SIANA¹, JOHAN RICHARD², MARC RAFELSKI³,⁴, MATHILDE JAUZAC⁵,⁶,⁷, MARCEAU LIMOUSIN⁸, WILLIAM R. FREEMAN¹, CLAUDIA SCARLATA⁹, BRANT ROBERTSON¹⁰, DANIEL P. STARK¹¹, HARRY I. TEPLITZ¹², AND VANDANA DESAI¹²

¹Department of Physics and Astronomy, University of California, Riverside, CA 92521, USA

²Centre de Recherche Astrophysique de Lyon, Université Lyon 1, 9 Avenue Charles André, F-69561 Saint Genis Laval Cedex, France

³Goddard Space Flight Center, Code 665, Greenbelt, MD 20771, USA

⁴Space Telescope Science Institute, Baltimore, MD, USA

⁵Centre for Extragalactic Astronomy, Department of Physics, Durham University, Durham DH1 3LE, UK

⁶Institute for Computational Cosmology, Durham University, South Road, Durham DH1 3LE, UK

⁷Astrophysics and Cosmology Research Unit, School of Mathematical Sciences, University of KwaZulu-Natal, Durban 4041, South Africa

⁸Aix Marseille Univ, CNRS, LAM, Laboratoire d'Astrophysique de Marseille, Marseille, France

⁹Minnesota Institute for Astrophysics, University of Minnesota, Minneapolis, MN 55455, USA

¹⁰Department of Astronomy and Astrophysics, University of California, Santa Cruz, 1156 High Street, Santa Cruz, CA 95064, USA

¹¹Department of Astronomy, Steward Observatory, University of Arizona, 933 North Cherry Avenue, Rm N204, Tucson, AZ 85721, USA

¹²Infrared Processing and Analysis Center, Caltech, Pasadena, CA 91125, USA

Received 2016 June 4; revised 2016 September 6; accepted 2016 September 13; published 2016 November 17

## ABSTRACT

We present a robust measurement of the rest-frame UV luminosity function (LF) and its evolution during the peak epoch of cosmic star formation at $1 < z < 3$. We use our deep near-ultraviolet imaging from WFC3/UVIS on the *Hubble Space Telescope* and existing Advanced Camera for Surveys (ACS)/WFC and WFC3/IR imaging of three lensing galaxy clusters, Abell 2744 and MACS J0717 from the Hubble Frontier Field survey and Abell 1689. Combining deep UV imaging and high magnification from strong gravitational lensing, we use photometric redshifts to identify 780 ultra-faint galaxies with $M_{UV} < -12.5$ AB mag at $1 < z < 3$. From these samples, we identified five new, faint, multiply imaged systems in A1689. We run a Monte Carlo simulation to estimate the completeness correction and effective volume for each cluster using the latest published lensing models. We compute the rest-frame UV LF and find the best-fit faint-end slopes of $\alpha = -1.56 \pm 0.04$, $\alpha = -1.72 \pm 0.04$, and $\alpha = -1.94 \pm 0.06$ at $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively. Our results demonstrate that the UV LF becomes steeper from $z \sim 1.3$ to $z \sim 2.6$ with no sign of a turnover down to $M_{UV} = -14$ AB mag. We further derive the UV LFs using the Lyman break “dropout” selection and confirm the robustness of our conclusions against different selection methodologies. Because the sample sizes are so large and extend to such faint luminosities, the statistical uncertainties are quite small, and systematic uncertainties (due to the assumed size distribution, for example) likely dominate. If we restrict our analysis to galaxies and volumes above $>50\%$ completeness in order to minimize these systematics, we still find that the faint-end slope is steep and getting steeper with redshift, though with slightly shallower (less negative) values ($\alpha = -1.55 \pm 0.06$, $-1.69 \pm 0.07$, and $-1.79 \pm 0.08$ for $z \sim 1.3, 1.9$, and $2.6$, respectively). Finally, we conclude that the faint star-forming galaxies with UV magnitudes of $-18.5 < M_{UV} < -12.5$ covered in this study produce the majority ($55\%–60\%$) of the unobscured UV luminosity density at $1 < z < 3$.

**Key words:** galaxies: evolution – galaxies: high-redshift – galaxies: luminosity function, mass function

## 1. INTRODUCTION

The galaxy luminosity function (LF) is a fundamental tool to study the formation and evolution of galaxies as the shape of the LF is mainly determined by the mechanisms that regulate star formation in galaxies (Rees & Ostriker 1977; White & Rees 1978; Benson et al. 2003). Comparing the LF with the underlying dark matter halo mass function reveals the importance of different modes of feedback in galaxy formation, with active galactic nucleus feedback dominating the bright end and supernova and radiation-driven winds dominating the faint end (e.g., Dekel & Birnboim 2006; Somerville et al. 2008). Furthermore, the LF is a key probe to assess the contribution of

galaxies with different luminosities to the total light budget at different redshifts.

As ultraviolet (UV) light is a tracer of recent star formation in galaxies, the UV LF can help determine the total star formation rate density at all epochs. In addition, the UV LF is one of the few galaxy observables that is directly measurable at all epochs using current telescopes. Over the past 20 years, many studies have been devoted to UV LF measurements at high redshifts with $z > 3$ (Steidel et al. 1999; Adelberger & Steidel 2000; Bunker et al. 2004; Dickinson et al. 2004; Ouchi et al. 2004; Yan & Windhorst 2004; Beckwith et al. 2006; Sawicki & Thompson 2006; Yoshida et al. 2006; Bouwens et al. 2007; Iwata et al. 2007; McLure et al. 2009; Ouchi et al. 2009; van der Burg et al. 2010; Bradley et al. 2012; Cucciati et al. 2012; McLure et al. 2013; Schenker et al. 2013; Atek et al. 2014, 2015a, 2015b; Schmidt et al. 2014b; Bouwens et al. 2015; Bowler et al. 2015; Finkelstein et al. 2015; Ishigaki et al. 2015), intermediate redshifts with $1 < z < 3$ (Dahlen

<footnote>* Some of the data presented herein were obtained at the W.M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California, and the National Aeronautics and Space Administration. The Observatory was made possible by the generous financial support of the W.M. Keck Foundation.</footnote>
et al. 2007; Reddy et al. 2008; Hathi et al. 2010; Oesch et al. 2010a; Cucciati et al. 2012; Sawicki 2012; Parsa et al. 2016) including our previous work (Alavi et al. 2014, hereafter A14), as well as low redshifts with $z < 1$ (Arnouts et al. 2005; Budavári et al. 2005; Wyder et al. 2005; Haberzettl et al. 2009; Ly et al. 2009; Cucciati et al. 2012). Taken together, these measurements suggest a rise and fall in the history of cosmic star formation from high redshifts to the present time with a peak sometime between $1 < z < 3$ (Madau & Dickinson 2014 and references therein). Therefore, the redshift range of $1 < z < 3$, known as the peak epoch of cosmic star formation, is a critical time in galaxy evolution.

Many wide and shallow surveys have probed the UV LF of rarer, luminous galaxies at $1 < z < 3$. Arnouts et al. (2005) used the WFPC2 data in the HDF-North and HDF-South and measured a faint-end slope of $\alpha = -1.5 \pm 0.2$ for the UV LF at $z = 2-3$. Later, Reddy & Steidel (2009) used a wide ground-based survey covering luminosities<fnref n="13" /> with $L > 0.05 L^*$ and measured a steep faint-end slope of $z = -1.73 \pm 0.07$ at $z = 2.3$. Following the installation of WFC3 on the Hubble Space Telescope (HST), Oesch et al. (2010a) used the wide, shallow Early Release Survey (ERS; Windhorst et al. 2011) and measured steep faint-end slopes ($-1.46 < \alpha < -1.84$) for the UV LFs at $z = 1.0 - 2.5$. However, in order to study the UV LF at fainter luminosities and accurately quantify the faint-end slope, deeper surveys were needed. In A14 (see next paragraph for more details), we used a very deep UV observation of the Abell 1689 (hereafter A1689) cluster obtained with the WFC3/UVIS channel, and we extended the $z \sim 2$ UV LF $100 \times$ fainter than previous shallower surveys ($L \sim 0.0005 L^*$). We concluded that the UV LF has a steep faint-end slope of $\alpha = -1.74 \pm 0.08$ with no evidence of a turnover down to $M_{UV} = -13$. Parsa et al. (2016) recently found galaxies as faint as $L > 0.002 L^*$ utilizing the CANDELS/GOODS-South, UltraVISTA/COSMOS, and HUDF data. However, their estimate of the faint-end slope $\alpha = -1.32 \pm 0.03$ is significantly shallower than others. A shortcoming of these two deep surveys is that they probe a single field (A1689 in A14 and HUDF dominating the faint luminosities in Parsa et al. 2016), where the field-to-field variations affect the LF measurements. In this paper, we attempt to overcome this problem by combining deep observations of three lines of sight.

Faint star-forming galaxies play a critical role in galaxy formation and evolution because they significantly contribute to IGM metal enrichment (Madau et al. 2001; Porciani & Madau 2005), are the most plausible sources of ionizing photons during the reionization epoch (Kuhlen & Faucher-Giguère 2012; Robertson et al. 2013), and maintain the ionizing background at $z > 3$ (Nestor et al. 2013). However, these faint galaxies are inaccessible at high redshifts as they lie outside of the detection limits of current surveys. One powerful way to explore these faint galaxies is to exploit the magnification of strong gravitational lensing offered by foreground massive systems and thus push the detection limits to lower luminosities. There have been many studies of high-redshift galaxies lensed by individual galaxies (e.g., Pettini et al. 2002; Siana et al. 2008b, 2009; Stark et al. 2008; Jones et al. 2010; Yuan et al. 2013; Vasei et al. 2016). However, galaxy clusters acting as gravitational lenses can magnify a

large area (e.g., Narayan et al. 1984; Kneib & Natarajan 2011), allowing a study of many highly magnified galaxies in a single pointing. In A14, combining our deep observations and magnification from strong gravitational lensing from A1689 enabled us to identify background ultra-faint galaxies.

This technique of targeting lensing galaxy clusters has been extensively used since the discovery of the first gravitationally lensed arc in the Abell 370 cluster (Soucail et al. 1987), and has culminated with recent large surveys of lensing clusters such as the CLASH (Postman et al. 2012) and Hubble Frontier Field (HFF; Lotz et al. 2016) surveys. The HFF program obtains very deep optical and near-infrared imaging over six lensing clusters using the HST/Advanced Camera for Surveys (ACS) and HST/WFC3, respectively. These deep images enable a search for the faint galaxies as opposed to the shallow CLASH data, which restrict the search to bright galaxies even in the case of high magnification. In addition, the HFF primary observations are complemented with data from Spitzer, ALMA, Chandra, XMM, VLA, VLT, and Subaru as well as our deep HST/WFC3 UV imaging in this study. Since the beginning of the HFF program, many groups have studied the faint-end of the UV LF at $z > 5$ up to $z = 9$ (Atek et al. 2014, 2015a, 2015b; Ishigaki et al. 2015; Livermore et al. 2016).

There are two primary methods of identifying high-redshift galaxies, via photometric redshifts and color-color selection of the Lyman break. Both techniques require assumptions about stellar populations, dust reddening, and star formation histories. However, each technique has its advantages. The photometric redshift method uses the full spectral energy distribution (SED) whereas the Lyman break method requires fewer filters and simpler completeness corrections. Some groups use the Lyman break technique (e.g., Hathi et al. 2010; Bouwens et al. 2015), while other groups prefer photometric redshifts (e.g., Finkelstein et al. 2015; Parsa et al. 2016). A general agreement between the UV LFs from these two methods is shown both at intermediate (Oesch et al. 2010a) and high-redshift studies (McLure et al. 2011, 2013; Schenker et al. 2013). One of the goals of this paper is to exploit the available multiwavelength imaging to provide a comparison between the UV LFs derived with these two selection techniques.

In this paper, we utilize the strong gravitational lensing magnification from three foreground galaxy clusters (two from the HFF program) in combination with our deep WFC3/UVIS imaging to construct a robust sample of faint star-forming galaxies at $1 < z < 3$. The study is similar to A14, but spanning the entire redshift range $1 < z < 3$, and measuring the LF behind three clusters instead of one. This allows us to study the evolution of the UV LF during the peak epoch of global star formation activity (i.e., $1 < z < 3$) and to compare with previous determinations. The structure of this paper is as follows. In Section 2, we summarize the available observations and the data reduction for each lensing cluster. The catalog construction and photometric redshift measurements are described in Sections 3 and 4, respectively. We briefly review the lens models and the multiple image identification in Section 5. We present our selection criteria and photometric redshift samples in Section 6. This is followed in Section 7, where we provide detailed description for the completeness simulation. We then discuss the UV LF measurements for the photometric redshift samples in Section 8 and for the dropout samples in Section 9. We compare the UV LFs obtained by different selection techniques, evolution of the UV LF and UV

<footnote n="13">¹³ To be consistent with other studies, we quote these limits in terms of $L_{z=3}^*$, i.e., $M_{1700,AB}^* = -21.07$, from Steidel et al. (1999).</footnote>
luminosity density in Section 10. Finally in Section 11, we provide a summary of our conclusions. In the appendices, we describe our color–color selection criteria, the corresponding LBG samples, and the completeness simulation for the LBG UV LF. We also provide a list of newly found multiple images of A1689.

In this paper, all distances and volumes are in comoving coordinates. All magnitudes are quoted in the AB system (Oke & Gunn 1983), and we adopt $\Omega_M = 0.3$, $\Omega_\Lambda = 0.7$, and $H_0 = 70 \text{ km s}^{-1} \text{ Mpc}^{-1}$.

## 2. DATA

In this section, we describe the data sets of three lensing fields used in this study and briefly explain the data reduction processes, as a more detailed description will be included in a future UV survey paper (B. Siana et al. 2016, in preparation). In this work we use deep HST imaging of three lensing clusters in a wide wavelength range, from UV to NIR, as described below.

### 2.1. HFF Observations and Data Reduction

The HFF survey uses the HST Director’s Discretionary time (GO/DD 13495, PI Lotz) to obtain deep WFC3/IR and ACS/ WFC images of six lensing clusters and their parallel fields (Lotz et al. 2016). The two HFF clusters analyzed here, Abell 2744 (hereafter A2744) and MACS J0717.5+3745 (hereafter MACS J0717), were observed during cycles 21 and 22, with 140 orbits of ACS/WFC and WFC3/IR imaging for each cluster/field pair. The NIR images are taken in the F105W, F125W, F140W, and F160W filters, and the optical data are obtained in the F435W, F606W, and F814W filters for each cluster.

In addition, we obtained deep near-ultraviolet images in F275W (eight orbits) and F336W (eight orbits) for three HFF clusters (including A2744 and MACS J0717) using the WFC3/ UVIS channel on board HST. These deep UV images are part of HST program ID 13389 (PI: B. Siana), which were taken between November 2013 and April 2014.

The Space Telescope Science Institute (STScI) handles the reduction and calibration of the optical and NIR images of the HFFs and releases the final mosaics in the Mikulski Archive for Space Telescopes (MAST).<fnref n="14" /> We used the version 1.0 release of the public optical and NIR mosaics with a pixel scale of 60 mas pixel⁻¹. To make these mosaics, the raw optical and NIR exposures were initially calibrated using the PYRAF/ STSDAS CALACS and CALWF3 programs, respectively. The calibrated images were then aligned and combined using the Tweakreg and AstroDrizzle (Gonzaga et al. 2012) tasks in the PYRAF/DrizzlePac package, respectively. In order to further improve the data reduction processes, the HFF team provides “self-calibrated” ACS images including more accurate dark subtraction and charge transfer efficiency (CTE) correction as well as WFC3/IR images corrected for the time-variable sky lines.

To calibrate the raw UV data, we applied two major improvements in addition to the standard WFC3/UVIS calibration approach. The first improvement is related to the CTE degradation of the UVIS CCD detectors. This degradation caused by radiation damage in the CCDs results in a loss of source flux and affects the photometry and morphology

<table><caption>**Table 1**
Observations and Image Depths</caption><thead><tr><th rowspan="2">Cluster<br>Instrument/<br>Filter</th><th colspan="2">A2744(HFF)</th><th colspan="2">MACS J0717<br>(HFF)</th><th colspan="2">A1689</th></tr><tr><th>Orbits</th><th>Depth<sup>a</sup></th><th>Orbits</th><th>Depth<sup>a</sup></th><th>Orbits</th><th>Depth<sup>a</sup></th></tr></thead><tbody><tr><td>WFC3/F225W</td><td>...</td><td>...</td><td>...</td><td>...</td><td>10</td><td>27.71</td></tr><tr><td>WFC3/F275W</td><td>8</td><td>27.80</td><td>8</td><td>27.43</td><td>30</td><td>28.14</td></tr><tr><td>WFC3/F336W</td><td>8</td><td>28.20</td><td>8</td><td>27.86</td><td>18</td><td>28.36</td></tr><tr><td>ACS/F435W</td><td>18</td><td>28.70</td><td>19</td><td>28.46</td><td>...</td><td>...</td></tr><tr><td>ACS/F475W</td><td>...</td><td>...</td><td>...</td><td>...</td><td>4</td><td>28.04</td></tr><tr><td>ACS/F606W</td><td>9</td><td>28.70</td><td>11</td><td>28.59</td><td>...</td><td>...</td></tr><tr><td>ACS/F625W</td><td>...</td><td>...</td><td>...</td><td>...</td><td>4</td><td>27.76</td></tr><tr><td>ACS/F775W</td><td>...</td><td>...</td><td>...</td><td>...</td><td>5</td><td>27.69</td></tr><tr><td>ACS/F814W</td><td>41</td><td>29.02</td><td>46</td><td>28.87</td><td>28</td><td>28.72</td></tr><tr><td>ACS/F805LP</td><td>...</td><td>...</td><td>...</td><td>...</td><td>7</td><td>27.30</td></tr><tr><td>WFC3/F105W</td><td>24.5</td><td>28.97</td><td>27</td><td>29.02</td><td>...</td><td>...</td></tr><tr><td>WFC3/F125W</td><td>12</td><td>28.64</td><td>13</td><td>28.60</td><td>...</td><td>...</td></tr><tr><td>WFC3/F140W</td><td>10</td><td>28.76</td><td>12</td><td>28.61</td><td>...</td><td>...</td></tr><tr><td>WFC3/F160W</td><td>24.5</td><td>28.77</td><td>26</td><td>28.65</td><td>...</td><td>...</td></tr></tbody></table>

**Note.**

measurements especially in low background images (e.g., UV data; Teplitz et al. 2013). To correct for these charge losses in our UV images, we used a pixel-based CTE correction tool provided on the STScI website.<fnref n="15" /> The second improvement is in the dark current subtraction from the UV images. As shown in a recent work by Teplitz et al. (2013), the standard WFCS3/ UVIS dark subtraction process is not sufficient for removing dark structures and hot pixels, mainly due to the low background level in the UV data. This regular technique leaves a background gradient and blotchy patterns in the final science image. Therefore, we used a new methodology introduced by Rafelski et al. (2015) for subtracting the dark current and masking the hot pixels properly. A detailed description of this technique is presented in Rafelski et al. (2015).

After making the modified calibrated UV images, we use the PYRAF/DrizzlePac package to drizzle these images to the same pixel scale of 60 mas and astrometrically align them with the optical and NIR data. The AstroDrizzle program subtracts the background, rejects the cosmic rays, and corrects the input images for the geometric distortion due to the nonlinear mapping of the sky onto the detector. In addition to the science output images, AstroDrizzle generates an inverse variance map (IVM) which we use later to make the weight images and to calculate the image depths. A summary of all the images and their depths is given in Table 1.

### 2.2. A1689 Observations and Data Reduction

In addition to the two HFF clusters, we observed the A1689 cluster. This cluster has been observed in three WFC3/UVIS bandpasses (F225W, F275W, and F336W) as part of program IDs 12201 and 12931 (PI: B. Siana), taken in cycle 18 in December 2010 and cycle 20 in February and March 2012, respectively. The cycle 18 data (30 orbits in F275W, 4 orbits in F336W) were used in A14 to measure the UV LF of lensed, dwarf galaxies at $z \sim 2$. In cycle 20, we added an F225W

<footnote>a 5σ limit in a 0.′′2 radius aperture.</footnote>

<footnote n="14">¹⁴ https://archive.stsci.edu/prepds/frontier/</footnote>

<footnote n="15">¹⁵ http://www.stsci.edu/hst/wfc3/tools/cte_tools</footnote>
image (10 orbits) and deeper F336W data (14 orbits, for a total of 18 orbits) to expand our redshift range from $1 < z < 3$.

The data calibration and reduction are the same as explained above for the HFF UV images. These data are corrected for the CTE degradation and dark subtraction, as well. Moreover, A1689 is observed with ACS/WFC in five optical bandpasses (F475W, F625W, F775W, F814W, and F850LP), which were calibrated and reduced as was described in A14. The A1689 images are all mapped to the same pixel scale of $40 \text{ mas pixel}^{-1}$.

### 3. OBJECT PHOTOMETRY

A detailed description for the A1689 photometry is given in A14. Here we provide the details of the photometric measurements for the HFF data. Since our HFF data cover a large range of wavelengths (from UV up to NIR), the width of the point-spread function (PSF) changes considerably. To perform multiband photometry, we match the PSF of all of the images to the F160W band, which has the largest PSF. We used the IDL routine `StarFinder` (Diolaiti et al. 2000) to stack all of the unsaturated stars in the field and extract the PSF. We fit a simple Gaussian function to each extracted PSF using the IRAF imexamine task, and then derive the PSF-matched images by convolving each band with a Gaussian kernel of appropriate width. We use `SExtractor` (Bertin & Arnouts 1996) to perform object detection and photometry. The final catalog areas are 4.81, 5.74, and 6.42 arcmin² where the WFC3 and ACS images are available for A2744, MACS J0717, and A1689, respectively.

We run `SExtractor` in dual image mode, with the F475W and F435W bands as detection images for the A1689 and HFF clusters, respectively. We use the F435W band to minimize contamination from the cluster galaxies and intracluster light, as this filter probes below the 4000 Å break, where the galaxies are considerably fainter. To improve the detection of faint objects and to avoid detecting spurious sources (i.e., overblended from very bright galaxies), for the `SExtractor` parameters, we set `DETECT MINAREA` to 4 (5) and `DETECT THRESH` to $0.9\sigma$ ($1.0\sigma$) significance for A2744 (MACS J0717). The minimum contrast parameter for deblending (`DEBLENS MINCONT`) is set to 0.02 for both cluster fields. The fluxes are measured in isophotal (ISO) apertures. The IVM images produced by the drizzling process as mentioned in Section 2.1 were converted to the rms MAPs by taking their inverse square root. `SExtractor` uses these rms MAPs to derive the flux uncertainties. We correct these rms MAPs for the correlated noise (A14; Casertano et al. 2000) from drizzling the mosaics. Finally, we correct our photometry for the Galactic extinction toward each cluster using the Schlafly & Finkbeiner (2011) IR dust maps. To account for systematic error (i.e., due to uncertainty in the Galactic extinction, the zero point values, PSF-matched photometry), we add, in quadrature, a 3% flux error (Dahlen et al. 2010; Vargas et al. 2014) in all bands for all three cluster fields.

### 4. PHOTOMETRIC REDSHIFTS

We use a template fitting function code, `EAZY` (Brammer et al. 2008), to estimate the photometric redshift of galaxies in all of our lensing fields. `EAZY` has two characteristic features that distinguish it from the other photometric redshift codes. First, it derives the optimized default template set from

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Comparison between the photometric and spectroscopic redshifts for 244 galaxies in all three lensing fields. The spectroscopic redshifts are either from our Keck/MOSFIRE and Keck/LRIS data or from the literature (for more details, see the text). The purple, orange, and cyan circles show the measurements for the A1689, MACS J0717, and A2744 cluster fields, respectively.</figcaption></figure>

semianalytical models with perfect completeness down to very faint magnitudes rather than using biased spectroscopic samples. Second, it has the ability to fit to a linear combination of basis templates rather than fitting to a single template, which is usually not a good representation of a real galaxy. We varied several `EAZY` input parameters to find the optimal values. Running `EAZY` using a variety of empirical (Coleman et al. 1980; Kinney et al. 1996) or stellar synthetic templates (Grazian et al. 2006; Blanton & Roweis 2007) allows us to find the set of models where the output photometric redshifts are in the best agreement with the spectroscopic redshifts. We use `PÉGASE` (Fioc & Rocca-Volmerange 1997) stellar synthetic templates, which provide a self-consistent treatment of nebular emission lines and include a wide variety of star formation histories (constant, exponentially declining) and a Calzetti dust attenuation curve (Calzetti et al. 2000). We do not use the template error function capability in `EAZY` because it causes poorer agreement with spectroscopic redshifts. We also do not use the magnitude priors, as these functions do not cover the faint luminosities targeted in this work. `EAZY` uses the Madau (1995) prescription for absorption from the intergalactic medium (IGM).

For the HFFs (A1689), we derive the photometric redshifts using the nine (eight) complete photometry bands of F275W, F336W, F435W, F606W, F814W, F105W, F125W, F140W, and F160W (F225W, F275W, F336W, F475W, F625W, F775W, F814W, and F850LP) with the central wavelengths covering 0.27 to 1.54 (0.24–0.91) µm. Figure 1 shows the comparison between the photometric redshifts and the spectroscopic redshifts for all three clusters. For both of the HFF clusters, we use the spectroscopic redshifts from the GLASS program, which obtained grism spectroscopy of 10 massive clusters including the HFFs (Schmidt et al. 2014a; Treu et al. 2015). We note that we only include their measurements with high quality parameter (i.e., quality >4) for a secure redshift estimate. In addition, for A2744, we also use the spectroscopic redshifts from the literature (Owers et al. 2011; Richard et al. 2014; Wang et al. 2015). For MACS J0717, we add the spectroscopic redshifts from our Keck/MOSFIRE spectral observations as well as the redshifts from the literature
(Limousin et al. 2012; Ebeling et al. 2014). Most of the spectroscopic redshifts of A1689 were described in A14, but here we also include our new measurements from our Keck/MOSFIRE spectra taken on 2015 January. A detailed study of spectroscopic data for these samples will be presented in a future paper. From the 186 galaxies with spectroscopic redshifts, 68 are within our target redshift range of $1 < z < 3$. For these galaxies with spectroscopic redshift of $1 < z_{\text{spec}} < 3$, we calculate the normalized median absolute deviation<fnref n="16" /> to be $\sigma_{\text{NMAD}} = 0.025$ (Ilbert et al. 2006) and find six outliers defined to have $\Delta z/(1+z_{\text{spec}}) > 5\sigma_{\text{NMAD}}$ (Brammer et al. 2008). The median and mean values of the fractional redshift error, $\Delta z/(1+z_{\text{spec}})$, after excluding outliers are 0.02 and 0.03, respectively.

Though the agreement between the photometric and spectroscopic redshifts is strong evidence for the reliability of our redshift estimates, it is restricted to the brighter galaxies. While our photometric redshift samples contain galaxies as faint as F606W (F625W for A1689) = 30 AB magnitudes, our spectroscopic samples cover magnitudes down to F606W (F625W for A1689) = 26.46. We note that among these objects, we have five galaxies at $1.2 < z_{\text{spec}} < 2.2$ with very faint magnitudes of $-15.4 < M_{\text{UV}} < -14$, where their spectroscopic and photometric redshifts agree well with the mean $\Delta z/(1+z_{\text{spec}}) = 0.04$. To further investigate the reliability of our photometric redshift estimates of the faint galaxies,<fnref n="17" /> where the spectroscopic redshifts are not available, we use a redshift quality parameter, $Q$.<fnref n="18" /> It is a statistical estimate of the reliability of the photometric redshift outputs of EAZY. Brammer et al. (2008) find that the photometric redshift scatter (i.e., difference between photometric redshift and spectroscopic redshift) is an increasing function of the $Q$ parameter with a sharp increase above $Q = 2-3$. We calculate the $Q$ parameter for our faint galaxies, as well as for the galaxies with spectroscopic redshifts of $1 < z < 3$. A comparison between these two sub-samples shows that the distributions of $Q$ values are similar (i.e., the faint galaxies are *not* skewed toward higher values of $Q$), such that the spectroscopic galaxies have a median $Q$ of 0.9, 1.1, and 1.5 relative to the faint galaxies with a median $Q$ of 0.6, 1.1, and 2.2 for the $z \sim 1.3, 2.2$, and 2.6 samples, respectively. We note that these values are within the safe regime for the $Q$ parameter (i.e., $Q < 3$, as explained above).

## 5. LENS MODELS

In order to estimate the intrinsic properties (i.e., luminosity) of the background lensed galaxies in our samples, we require an accurate mass model of the galaxy cluster to calculate the lensing magnification. For the HFF program, there are several groups working independently to use deep HFF optical and NIR imaging to model the mass distribution for all of the six clusters (Bradač et al. 2005; Liesenborgs et al. 2006; Diego et al. 2007; Jullo et al. 2007; Jullo & Kneib 2009; Merten

et al. 2009; Zitrin et al. 2009, 2013; Oguri 2010; Merten et al. 2011; Sendra et al. 2014). The main distinction between these models is that some groups assume light traces mass and parametrize the total mass distribution as a combination of individual cluster members and large-scale cluster halo components, while the other groups use a non-parametric mass modeling technique, avoiding any priors on the light distribution. In a recent study, Priewe et al. (2016) provide a comparison between these different lens models. All of these models are constrained by the location and the redshift of known multiply imaged systems. Besides observational constraints from strong gravitational lensing, several teams also incorporate the weak-lensing shear profile from ground-based observations. All of the HFF lens models and the methodologies adopted by each team are publicly available via the STScI website.<fnref n="19" /> In this section, we briefly review the mass models that we used for each of our lensing clusters.

### 5.1. HFF Lensing Models

For the HFF clusters, we utilize the lens models produced by the Clusters As TelescopeS (CATS) collaboration (Co-PIs J.-P. Kneib and P. Natarajan; Admin PI H. Ebeling) who use the *Lenstool* software<fnref n="20" /> (Jullo et al. 2007) to parameterize the lens mass distribution. *Lenstool* is a hybrid code that combines both strong- and weak-lensing data to constrain the lens mass model. *Lenstool* models each cluster's mass as a composition of one or more large cluster halos plus smaller subhalos associated with individual galaxies identified either spectroscopically or photometrically as cluster members. The output best model from *Lenstool* is parameterized through a Bayesian approach.

For A2744, we use the strong lensing model of Jauzac et al. (2015), which uses 61 multiply imaged systems found in the complete HFF optical and NIR data. For MACS J0717, we use the strong lensing model of Limousin et al. (2016), which uses 55 multiply imaged systems found in the complete HFF optical and NIR data.

### 5.2. A1689 Lensing Model

As in A14, the lens model that we use for A1689 is from Limousin et al. (2007). Similar to the mass reconstruction techniques for the HFF clusters, Limousin et al. (2007) optimize a parametric model implemented in the *Lenstool* using 32 multiply imaged systems behind A1689. Their optimized lens model for A1689 is a composite of two large-scale halos and the subhalos of cluster member galaxies.

### 5.3. Multiply Imaged Systems

Finding more multiply imaged systems is critical for improving a lens model, as the lens model is constrained by the location and redshift of these systems. In addition, identifying the multiple images is important as we need to remove them from the galaxy number counts. We run *Lenstool* using each previously described lens model as an input to look for the potential counter-images for each lensed galaxy in the sample.

Currently, there is no automated process for identifying multiply imaged systems. Here, we summarize the approach

<footnote n="16">¹⁶ The normalized median absolute deviation is defined as $\sigma_{\text{NMAD}} = 1.48 \times \text{median}(|\Delta z - \text{median}(\Delta z)|/(1+z))$ (Ilbert et al. 2006; Brammer et al. 2008). Unlike the usual standard deviation, $\sigma_{\text{NMAD}}$ is not sensitive to the presence of outliers.</footnote>

<footnote n="17">¹⁷ We define the faint galaxies based on the limiting magnitude used in our sample selection criteria (see Section 6). They are defined to have a signal-to-noise ratio (S/N) < 5 in either detection filter or the rest-frame 1500 Å filter.</footnote>

<footnote n="18">¹⁸ The $Q$ parameter (see Equation (8) in Brammer et al. 2008) combines the reduced-$\chi^2$ of the fitting procedure with the width of the 68% confidence interval of the redshift probability distribution function to present an estimate of the reliability of the output redshift.</footnote>

<footnote n="19">¹⁹ https://archive.stsci.edu/prepds/frontier/lensmodels/</footnote>

<footnote n="20">²⁰ https://projets.lam.fr/projects/lenstool/wiki</footnote>
that we took to find new multiply imaged systems. (1) We run Lenstool, entering the coordinates and the photometric redshift of each galaxy to predict the location of its potential counter-images. In this step, Lenstool first de-lenses the galaxy image to its original position in the source plane at the given photometric redshift, and then re-lenses it back to all of the possible multiple image positions in the image plane. (2) We search for the objects with the same color and symmetry in the morphology near the predicted positions. (3) If we find any nearby candidate from step 2, we then repeat the first step to check if the potential counter-images of the candidate match with the first object. (4) Finally, we require the same photometric redshifts (within 1σ accuracy) for all of the newly found multiple images. This final criterion exhibits the importance of covering rest-frame UV wavelengths, which enables us to identify the Lyman break to distinguish the high-redshift objects (in this case $1 < z < 3$) from the lower redshift interlopers, since both often have flat, featureless SEDs at rest-frame optical wavelengths.

Following this procedure for all of the galaxies, we find five and three new multiply imaged systems behind A1689 and MACS J0717, respectively. Our new findings in cluster MACS J0717 added new systems 21, 80, and 82 to the list reported in Limousin et al. (2016). We introduce the new A1689 multiply imaged systems in Appendix C.

## 6. SAMPLE SELECTION

We use the photometric redshift estimates to construct our galaxy samples in three redshift ranges of $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$. To ensure the reliability of our photometric redshifts and to avoid selecting spurious objects in the sample, we require $3\sigma$ detections in the detection filter and the rest-frame 1500 Å filter. The selection criteria for the lower redshift range are:

a. $1.0 < z_{phot} < 1.6$

b. S/N > 3 in the F275W and F336W bands,

selecting 70, 134, and 93 candidates in A1689, A2744, and MACS J0717, respectively. The selection criteria for the middle redshift range are

a. $1.6 < z_{phot} < 2.2$

b. S/N > 3 in the F336W and F435W (F475W) bands for the HFFs (A1689),

selecting 128, 121, and 69 candidates in A1689, A2744, and MACS J0717, respectively. And finally, the selection criteria for the higher redshift range are

a. $2.2 < z_{phot} < 3.0$

b. S/N > 3 in the F435W and F606W bands for the HFFs,

selecting 176 and 102 galaxies in the A2744 and MACS J0717 fields, respectively. We should note that we do not include data from A1689 for the highest redshift ($z \sim 2.6$) analysis because, due to the cluster redshift of $z = 0.18$, the Balmer break of faint cluster members (like globular clusters; Alamo-Martínez et al. 2013) can mimic the Lyman break at $z \sim 3$.

In total, we have 297, 318, and 278 candidates at $1 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively. As explained in Section 5.3, we must clean our samples of multiple images. Among each multiply imaged system, we keep the brightest image and remove the rest of the images from our samples. However, if the brightest image has a

magnification higher than 3.0 mag, we then select the next brightest image. This condition on magnification is considered to ensure the reliability of the magnification value predicted from the lensing models.

Furthermore, to ensure the purity of the samples, we consider different possibilities of contamination in the photometric redshift-selected samples. First, to find possible contamination from stars, we use the Pickles (1998) stellar spectra library to predict stellar colors for a variety of stars and compare with the color of our candidate galaxies. In the case of similar colors, we visually inspect the objects. We found only 1 (~0.3%), 2 (~0.6%), and 0 stars in the $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$ samples, respectively. We also visually inspect all of the galaxies to exclude objects associated with diffraction spikes and nearby bright galaxies. The contamination is only 2 (~0.7%), 3 (~0.9%), and 0 for the $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$ samples, respectively. Finally, after excluding all of the multiple images and the contamination, we have 277, 269, and 252 galaxies at $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively.

With the aim to measure the UV LF, we use the F336W, F435W, and F606W bands for the HFFs and the F336W and F475W bands for the A1689 samples to measure the absolute magnitude at rest-frame 1500 Å ($M_{UV} = M_{1500}$) at redshifts $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively. As we did in A14, we determine the intrinsic absolute magnitudes of $M_{1500}$ by applying the magnification corrections computed from the lens models discussed in Section 5:

$$M_{1500} = m + \mu_{mag} - 5 \log(d_L/10 \text{ pc}) + 2.5 \log(1+z), \quad (1)$$

where $\mu_{mag}$ is the predicted magnification in magnitude units from the lensing model of each cluster. We limit our samples to galaxies brighter than $M_{1500} < -12.5$ mag to ensure a reliable absolute magnitude measurement. All of the galaxies brighter than this limit have magnification uncertainty from lensing models below 0.5 mag with a mean value of 0.03 mag, but the galaxies fainter than $M_{1500} > -12.5$ have magnification uncertainties above 2.0 mag. This limit excludes 7 (2.5%), 10 (3.7%), and 1 (0.4%) galaxies from the $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$ samples, respectively.

Figure 2 shows the distributions of magnifications for galaxies in three photometric redshift samples. The magnification values range between $\mu_{mag} = 0.5-4.8$ (equivalent to 1.58-83 in flux density units) with median values of $\mu_{mag} = 1.27$, 1.61, and 1.24 for the $z \sim 1.3$, 1.9, and 2.6 samples, respectively. As shown in Figure 2, most of the highly magnified galaxies ($\mu_{mag} > 2.5$) in the $z \sim 1.3$ and $z \sim 1.9$ samples are from A1689. We note that because A1689 has a large Einstein radius, it provides high magnification (i.e., median $\mu_{mag} = 2.0$) over a large area in the source plane. Therefore, objects with a high magnification in A1689 are not required to be close to the critical lines, where the magnification formally diverges. For example, the galaxies with high magnification ($\mu_{mag} = 2.5-4$) in the A1689 sample are on average 15 arcsec (with median of 10 arcsec) from the critical lines whose positions are predicted with a precision of 2.87 arcsec by the lens model (Limousin et al. 2007). Therefore, these magnification estimates are not strongly affected by uncertainties in the location of the critical lines.
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Magnification distribution of galaxies expressed in magnitude units. Purple, orange, and cyan show the number of candidate galaxies for each magnification bin on A1689, A2744, and MACS J0717, respectively. These clusters provide a large range of magnifications, with higher values mostly from A1689 (see the text).</figcaption></figure>

In Figure 3, we show the histograms of absolute UV magnitudes for each lensing cluster in three redshift bins. This figure emphasizes the importance of including A1689, since it dominates the number of galaxies at the faintest magnitudes, $M_{1500} > -14.5$.

## 7. COMPLETENESS SIMULATIONS

In order to connect the observed galaxies to the underlying population of all star-forming galaxies, we need to precisely estimate the completeness of our sample. This is more critical for low luminosity bins, where the galaxies are close to the detection limits. An approach commonly used in the blank field studies to estimate the completeness (e.g., Oesch et al. 2010a; Grazian et al. 2011; Bowler et al. 2015; Finkelstein et al. 2015) is to generate artificial galaxies with properties similar to the real galaxies and then apply an identical selection technique as for the observed candidates to calculate the fraction of recovered simulated galaxies in a given magnitude and redshift bin. This technique is also applicable in gravitationally lensed studies (e.g., Atek et al. 2015a, 2015b). However, one needs to incorporate the added complexity due to the strong lensing amplification.

In this work, we adopt a Monte Carlo simulation following the methodology presented in detail in A14. Here, we briefly describe these completeness simulations, and we provide additional details where our approach deviates from what was done in A14.

We compute the completeness in a 3D grid of redshift, magnitude, and lensing magnification. For each point in this 3D space, we assign a redshifted and magnified template galaxy spectrum, which is generated by Bruzual & Charlot (2003; hearafter BC03) synthetic stellar population models assuming a 0.2 Z⊙ metallicity and an age of 100 Myr. A detailed justification for these assumptions is given in A14. The SED is dust attenuated using the Calzetti extinction curve (Calzetti et al. 2000) and a random color excess, $E(B - V)$, value taken from a Gaussian distribution centered at 0.15 as measured in A14 and other studies (Steidel et al. 1999; Reddy & Steidel 2009; Hathi et al. 2013) with a standard deviation of 0.1. In order to understand the effect of a changing reddening

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Intrinsic absolute UV magnitude (i.e., corrected for the lensing magnification) distribution of all galaxies in our three redshift slices. The colors are similar to Figure 2, and they show the number of candidate galaxies for each absolute magnitude bin. We cover a large luminosity range by combining the HFFs with A1689, which finds the faintest galaxies ($M_{1500} > -14.5$).</figcaption></figure>

distribution, we also examined the completeness for a model in which the dust reddening linearly decreases toward fainter luminosities. To derive this linear function, we measured the relation between the UV spectral slope and $M_{1500}$ magnitude for our galaxies, and we calculate the dust reddening values assuming a Calzetti reddening curve. The final completeness corrections from this examination show only negligible changes relative to our original simulations.<fnref n="21" />

We then create transmission curves (as a function of wavelength) for 300 lines of sight through the IGM at that redshift. The IGM opacity is calculated using a Monte Carlo simulation to randomly place hydrogen absorbers in each line of sight as described in A14 (see also Siana et al. 2008a). Our completeness simulation is modified relative to A14 in the following two ways.

**Updating the Size Distribution of Star-forming Galaxies:** One of the key factors in estimating the incompleteness is the assumed size distribution for galaxies. As shown in Grazian et al. (2011), the completeness correction at low luminosities depends critically on the adopted size distribution in the simulation, as using too small (large) a size distribution can cause one to over- (under-) estimate the completeness. As reported in various observational studies (e.g., Bouwens et al. 2004; Ferguson et al. 2004; Huang et al. 2013), the rest-frame UV sizes of high-redshift Lyman break galaxies follow a log-normal distribution. In a recent work, Shibuya et al. (2015) measured the size distribution of a large sample of galaxies at $0 < z < 8$, using the 3D-HST and CANDELS data. They showed that the circularized effective radius<fnref n="22" /> ($r_e$) distribution of star-forming galaxies at $0 < z < 8$ is well represented by a log-normal distribution whose median

<footnote n="21">²¹ We note that for the same experiment, the effective volumes of the LBG samples (see Appendix B) show slightly larger change at bright luminosities. This can be understood by considering that the color–color criteria select against very reddened galaxies. However, our final estimates of the best-fit LFs (for both sample selections) are robust against these different initial assumptions of dust reddening distribution.</footnote>

<footnote n="22">²² The circularized effective radius is defined as $r_e = r_{e,\text{major}} \sqrt{q}$, where $r_{e,\text{major}}$ is the half-light radius along the semimajor axis and $q$ is the axis ratio. The circularized radius has been extensively used in other high-redshift size measurements (e.g., Mosleh et al. 2012; Ono et al. 2013).</footnote>
decreases toward high redshifts (at a given luminosity) and changes with luminosity as $r_e \propto (L_{UV})^\alpha$ with $\alpha = 0.27$ for all redshifts. For our completeness simulation, we generate random galaxy sizes at each luminosity and redshift using the corresponding log-normal distribution from Table 8 in Shibuya et al. (2015). We extrapolate their measurements below $M_{UV} < -16$. Using the randomly selected $r_e$ values, we then adopt a Sérsic profile with index $n = 1.5$ as suggested by Shibuya et al. (2015). Other LF studies at both low (Oesch et al. 2010a) and high redshifts (Oesch et al. 2010b; Grazian et al. 2011; Atek et al. 2015a; Finkelstein et al. 2015) have also assumed a log-normal size distribution. Our size distribution assumption in this work is different from A14, where we assumed a normal (not a log-normal) distribution centered at 0.7 kpc with a standard deviation of 0.2 kpc (Law et al. 2012).

**Updating the The Effect of Lensing Magnification and Shear:** The next step in the simulation is to add the lensing effect by amplifying the flux and enlarging the size of the galaxies. The way that gravitational lensing distorts the image of a galaxy is a combination of convergence (i.e., $\kappa$, stretching a source isotropically) and shear (i.e., $\gamma$, stretching a source along a privileged direction). As discussed in other works (e.g., Oesch et al. 2015), it is crucial to account for the effect of lensing distortion. In A14, we did include the effect of convergence in distorting our simulated galaxies. For this work, we do a more complete and complex analysis such that the shape of the final distorted image can be described using tangential $(\mu_t = (1 - \kappa - \gamma)^{-1})$ and radial $(\mu_r = (1 - \kappa + \gamma)^{-1})$ magnification. As formulated in Bartelmann (2010), a circular source with a circularized radius of $r_e$ becomes an elliptical image with semimajor (a) and semiminor (b) axes as below:

$$a = \mu_r r_e, \qquad (2)$$

$$b = \mu_t r_e. \qquad (3)$$

We use the lensing models to construct the $\mu_t$ and $\mu_r$ maps at source plane for desired redshifts. We then use these maps to select random $(\mu_t, \mu_r)$ pairs and distort the image of our simulated galaxies. We also increase the flux with a magnification factor of $\mu = \mu_t \cdot \mu_r$. To mimic the same condition as real galaxies, we should note that we exclude the large cluster members from our source plane area reconstruction as the real galaxies behind these low-$z$ intervening galaxies cannot be observed.

The corresponding synthetic SED assigned to each simulated source is multiplied by the same filter curves used in the observations to generate artificial catalogs. We then add random photometric noise to the distorted image of each galaxy for each band. To detect the galaxies and generate the artificial catalogs, we use the same detection parameters we used in SExtractor for our real galaxies.

Finally, for each cell of the 3D grid, we have a SExtractor output catalog for 300 artificially created galaxies in random lines of sight, with random sizes and dust attenuation values sampled from the corresponding distributions explained above. We then run the EAZY code on these simulated catalogs and adopt the same selection criteria as we did for the real sources (see Section 6). Consequently, we calculate the completeness correction factor, $C(m, z, \mu_{mag})$, as a function of intrinsic apparent magnitude (i.e., before magnification, $m$), redshift ($z$), and magnification ($\mu_{mag}$) by counting the fraction of recovered artificial galaxies. Figure 4

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Completeness as a function of intrinsic apparent magnitude on the y-axis and redshift on the x-axis. The blue, orange, and purple contours show the completeness simulation values for A2744, MACS J0717, and A1689, respectively. The left, middle, and right columns represent the completeness contours for photometric redshift samples at $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively. These contours are drawn for a magnification of $\mu_{mag} = 2.0$ magnitudes.</figcaption></figure>

shows the completeness contours as a function of intrinsic apparent magnitude, $m$, on the y-axis and redshift on the x-axis for each redshift interval and for each lensing cluster with different colors. The contours are plotted for a magnification of $\mu_{mag} = 2.0$ mag. We can see the difference between HFFs and A1689 completeness values at the lower redshift range ($1.0 < z < 1.6$), where F225W photometry in A1689 helps to better constrain the redshift and avoids contamination from galaxies with input redshifts below 1.0. As seen in this figure, the recovered redshift distribution from completeness simulations is in agreement with our targeted redshift ranges for each sample.

## 7.1. The Effective Survey Volume

We incorporate the completeness corrections in the computation of the effective survey volume, $V_{eff}$, in each magnitude bin as below:

$$V_{eff}(m) = \int_{0}^{\infty} \int_{0}^{\infty} \frac{dV_{com}}{dz d\Omega} C(z, m, \mu) \Omega(\mu, z) dz d\mu \quad (4)$$

where $dV_{com}$ is the comoving volume element at redshift $z$ per unit area, $d\Omega$. In this equation, $C(z, m, \mu)$ is the completeness function that depends on redshift ($z$), intrinsic apparent magnitude ($m$), and magnification ($\mu$). $\Omega(\mu, z)$ is the area element in the source plane at $z$ which is magnified by a factor of $\mu$. We run Lenstool for each aforementioned cluster mass model to generate the de-lensed magnification maps at different redshifts. We then use these maps to estimate the $\Omega(\mu, z)$ of each cluster at each redshift. Similar to our completeness simulations (see Section 7), we subtract the area occupied by
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Effective volume estimates at each redshift slice in each field. The HFF clusters provide a large volume over faint magnitudes ($M_{UV} < -15.5$), while the A1689 cluster enables even fainter galaxies ($M_{UV} > -14.5$) beyond the HFF magnification limits to be probed.</figcaption></figure>

the large cluster members from our source plane area reconstruction.

Figure 5 represents the effective volumes versus the absolute magnitude at 1500 Å, $M_{UV}$, for each cluster at three redshift ranges. This plot clearly shows the importance of including A1689 to find the faintest galaxies ($M_{UV} > -14.5$). We should emphasize that the small volumes at faint luminosities are not necessarily due to a large incompleteness but because of the small area available at these magnitudes. For the volume calculation at each magnitude, unlike the field studies where the full area is available, here only a portion of area (i.e., effective area) with enough magnification (i.e., minimum magnification required for detection at each magnitude) is used. Therefore, at very faint luminosities, only a tiny fraction of area is available for the volume measurements.

## 8. LF OF PHOTOMETRIC REDSHIFT SAMPLES

Using the effective volumes, we construct the UV LF of our photometric redshift-selected galaxies at the peak epoch of cosmic star formation rate density. To be consistent with other studies at the same redshift ranges (e.g., Oesch et al. 2010a; Parsa et al. 2016) and at higher redshifts (e.g., Bouwens et al. 2007; Finkelstein et al. 2015), we measure the UV luminosities at rest-frame 1500 Å.

The galaxy LF is commonly fitted by a Schechter function (Schechter 1976) characterized by an exponential behavior at luminosities brighter than a characteristic magnitude, $M^*$, and a power law at the faint end with slope $\alpha$ as below:

$$ \phi(M) = 0.4 \ln(10) \phi^* 10^{-0.4(M-M^*)}(1+\alpha) e^{-10^{-0.4(M-M^*)}} \quad (5) $$

where $\phi^*$ is the normalization of this function.

In this section, we first calculate and compare binned UV LFs of each cluster field and then we find the best-fit Schechter parameters for the combined LF using a maximum likelihood (MLE) approach on the unbinned data.

### 8.1. The Binned UV LFs

The LF at each $M_{1500}$ bin is derived using the measured $V_{eff}$ values which account for the completeness corrections. This is the commonly used $V_{eff}$ method (e.g., A14, Oesch et al. 2010a)

where one calculates the number density of galaxies in each bin by dividing the number of galaxies in the corresponding absolute magnitude bin by the effective volume of that bin. But the effective volume might change significantly from one side of the magnitude bin to the other. Therefore, we estimate the effective volume for each individual galaxy and then sum up over all the galaxies within each bin, as shown below:

$$ \phi(M_i) dM_i = \sum_{j=1}^{N} \frac{1}{V_{\text{eff}}(M_j)} \quad (6) $$

As illustrated in Figure 6, we estimate the binned LFs of each lensing field separately as well as a total LF combining all of the cluster fields. For the combined LF, the $V_{eff}$ is a sum of the effective volumes over all of the cluster fields.

For each bin with a large number of galaxies ($N > 50$), we assign an uncertainty of $\frac{\phi_i}{\sqrt{N}}$ using Poisson statistics. In the case where fewer than 50 galaxies are in the bin, we compute the Poisson approximation, $\Delta_p$, from Gehrels (1986) and assign an uncertainty of $\frac{\phi_i \Delta_p}{N}$ to each bin. Each bin has a width of $\Delta M_{UV} = 1$ magnitude and our faintest magnitude bin is centered at $M_{UV} = -13$ (i.e., a magnitude cut at $M_{UV} = -12.5$, see Section 6). The values of the binned LFs and the number of galaxies at each bin are listed in Table 2.

The binned LFs are good for visualization but poor for inference because of arbitrary bin widths, bin centers, and loss of information within each bin. Therefore, instead of using binned estimators, we use an unbiased, unbinned MLE estimator as explained in the next section.

### 8.2. The Unbinned MLE Estimator

In this section, we explain our methodology to estimate the best Schechter function parameters by maximizing the likelihood function of the unbinned data. The standard MLE technique was first used by Sandage et al. (1979, STY79), and later by many other studies to derive the best-fit parameters for UV LFs at intermediate redshifts (A14), high redshifts (e.g., McLure et al. 2013; Bouwens et al. 2015), and for the H$\alpha$ LF (Mehta et al. 2015). Here, we adopt a similar approach as in A14 where we modify the standard STY79 MLE technique to account for uncertainties in the measurements of the absolute magnitude. This modified methodology is also used in Mehta et al. (2015). In the MLE technique, the best fit is found by maximizing the joint likelihood function $L$ defined as below:

$$ \mathcal{L} = \prod_{i=1}^{N} P(M_i) \quad (7) $$

where in the standard MLE, $P(M_i)$ is defined as below:

$$ P(M_i) = \frac{\phi(M_i) V_{\text{eff}}(M_i)}{\int_{-\infty}^{M_{\text{limit}}} \phi(M) V_{\text{eff}}(M) dM} \quad (8) $$

where $N$ is the total number of objects in each sample. $P(M_i)$ is the probability of finding a galaxy with absolute magnitude $M_i$ in a corresponding effective volume, $V_{eff}(M_i)$. We calculate this probability value for all of the galaxies in each of our samples. $\phi(M_i)$ is the parametric LF assuming a Schechter function. The $M_{limit}$ is defined for each sample to be the faintest absolute magnitude (i.e., corrected for the magnification). The $M_{limit}$ values are -12.88, -12.12, and -13.40 for the $z = 1.3, 1.9, 2.6$ samples, respectively.
<figure><img src="image_7.png" /><figcaption>Figure 6. Rest-frame UV luminosity function for each lensing cluster at 1.0 < z < 1.6, 1.6 < z < 2.2, and 2.2 < z < 3.0 in the left, middle, and right panels, respectively. The purple, blue, and orange stars show the binned LF of A1689, A2744, and MACS J0717, respectively. The black circles are the binned LFs after combining all of the three lensing clusters. The dashed line is the best-fit Schechter function (see Section 8.2).</figcaption></figure>

To incorporate absolute magnitude uncertainties in the LF analysis, we assume a Gaussian probability distribution $G(M|M_i, \sigma_i)$ for each object centered at the object's absolute magnitude $M_i$ and a standard deviation equal to the object's absolute magnitude uncertainty $\sigma_i$. We then modify Equation (8) as below:

$$P(M_i) = \frac{\int_{-\infty}^{+\infty} \phi(M) V_{\text{eff}}(M) G(M|M_i, \sigma_i) dM}{\int_{-\infty}^{M_{\text{limit}}} \phi(M) V_{\text{eff}}(M) dM} \quad (9)$$

with

$$G(M|M_i, \sigma_i) = \frac{1}{\sqrt{2\pi}\sigma_i} \exp\left(-\frac{(M - M_i)^2}{2\sigma_i^2}\right) \quad (10)$$

As also considered in A14, for our lensed galaxies, the total uncertainty, $\sigma_i$, of the intrinsic absolute magnitude is due to the uncertainty in photometric measurements ($\sigma_m$), photometric redshifts ($\sigma_z$), and the lens models ($\sigma_{\text{model}}$). Below, we investigate in detail these different sources of uncertainties.

a. $\sigma_m$: The photometric uncertainties are calculated using the SExtractor output of flux uncertainties.

b. $\sigma_z(\text{total})$: The photometric redshift uncertainty, $\sigma_z$, for each galaxy is computed as the $1\sigma$ confidence interval of its redshift probability distribution from EAZY. This redshift uncertainty impacts the measured intrinsic absolute magnitude in two ways. First, since the distance modulus is dependent on the redshift, we estimate the effect of redshift uncertainty on the absolute magnitude through an error propagation of Equation (1). Second, the magnification value of each galaxy is estimated by running Lenstool while incorporating its photometric redshift as an input. Therefore, a redshift uncertainty causes a magnification uncertainty, $\sigma_{\mu(z)}$. To estimate $\sigma_{\mu(z)}$ for each galaxy, we run Lenstool for 100 random redshifts generated from a Gaussian redshift distribution centered at the galaxy's photometric redshift with standard deviation equal to $\sigma_z$. The distribution of output random magnifications for each galaxy is fitted with a Gaussian function to derive $\sigma_{\mu(z)}$. Because $\sigma_z$ and $\sigma_{\mu(z)}$ are

correlated, we calculate the total redshift uncertainty as a sum over them, $\sigma_z(\text{total}) = \sigma_z + \sigma_{\mu(z)}$

c. $\sigma_{\mu(\text{model})}$: The final source of uncertainty is related to the lensing models. To estimate this uncertainty, we randomly sample the parameter space of each lens model. A detailed description of these measurements is given in A14.

We calculate the total uncertainty of the intrinsic absolute magnitude by adding all these uncertainties in quadrature.

Substituting Equation (9) in Equation (7), we calculate the likelihood function over a grid of faint-end slope ($\alpha$) and characteristic magnitude ($M^*$). The small survey areas probed in this study limits the number of bright galaxies (i.e., $M < M^*$). Therefore, to constrain $M^*$, we combine our $z \sim 1.3$ and $z \sim 1.9$ samples with the samples from a wider survey from Oesch et al. (2010a). To be consistent with our samples, we use their photometric redshift-selected galaxies at $1.0 < z < 1.5$ and $1.5 < z < 2.0$.

For our $2.2 < z < 3.0$ sample, our brightest LF bins are lower than the values from the literature. Furthermore, we do not have access to the individual galaxies from the literature. Therefore, we adopt a different approach to find the best-fit LF. We multiply the likelihood function by an $M^*$ prior to compute the posterior function. Utilizing the best Schechter parameters reported in Reddy & Steidel (2009), we define the prior as a Gaussian function centered at -20.70 with a standard deviation of 0.11. We should note that this discrepancy between the LFs at bright luminosities is not due to our completeness correction, as our $z \sim 2.6$ sample is >90% complete at these luminosities (see Figures 4 and 5). Considering that we only have two clusters at this redshift range, and consequently we probe a small area, it is not unlikely that this low number density may be due to a presence of an underdense region of galaxies. The reason that this underdensity appears to be affecting the bright end more than the faint end can be understood by the different spatial clustering of the bright galaxies relative to the faint ones (e.g., Zehavi et al. 2005).<fnref n="23" />

<footnote n="23">²³ In the future, when we complete the UV survey of HFFs, we will add four more clusters and consequently triple our sample size at $z = 2.6$. Therefore, our number density measurement for bright galaxies at this redshift will be more accurate.</footnote>
<table><caption>**Table 2**
Binned UV LFs</caption><thead><tr><th>z</th><th>M<sub>UV</sub></th><th>Number of Sources</th><th>φ (×10<sup>-2</sup> Mpc<sup>-3</sup> mag<sup>-1</sup>)</th></tr></thead><tbody><tr><td colspan="4">Photometric Redshift LFs</td></tr><tr><td rowspan="8">1.0 &lt; z &lt; 1.6</td><td>-20.0</td><td>6</td><td>0.069<sup>+0.041</sup><sub>-0.027</sub></td></tr><tr><td>-19.0</td><td>27</td><td>0.320<sup>+0.074</sup><sub>-0.061</sub></td></tr><tr><td>-18.0</td><td>41</td><td>0.490<sup>+0.089</sup><sub>-0.076</sub></td></tr><tr><td>-17.0</td><td>41</td><td>0.543<sup>+0.099</sup><sub>-0.084</sub></td></tr><tr><td>-16.0</td><td>60</td><td>1.064<sup>+0.137</sup><sub>-0.137</sub></td></tr><tr><td>-15.0</td><td>54</td><td>2.266<sup>+0.308</sup><sub>-0.308</sub></td></tr><tr><td>-14.0</td><td>28</td><td>8.275<sup>+1.877</sup><sub>-1.554</sub></td></tr><tr><td>-13.0</td><td>13</td><td>21.279<sup>+7.693</sup><sub>-5.826</sub></td></tr><tr><td rowspan="8">1.6 &lt; z &lt; 2.2</td><td>-20.0</td><td>5</td><td>0.095<sup>+0.064</sup><sub>-0.041</sub></td></tr><tr><td>-19.0</td><td>11</td><td>0.217<sup>+0.087</sup><sub>-0.064</sub></td></tr><tr><td>-18.0</td><td>31</td><td>0.615<sup>+0.131</sup><sub>-0.110</sub></td></tr><tr><td>-17.0</td><td>62</td><td>1.307<sup>+0.166</sup><sub>-0.166</sub></td></tr><tr><td>-16.0</td><td>69</td><td>2.186<sup>+0.263</sup><sub>-0.263</sub></td></tr><tr><td>-15.0</td><td>40</td><td>3.964<sup>+0.731</sup><sub>-0.624</sub></td></tr><tr><td>-14.0</td><td>35</td><td>19.223<sup>+3.828</sup><sub>-3.235</sub></td></tr><tr><td>-13.0</td><td>6</td><td>19.104<sup>+11.412</sup><sub>-7.578</sub></td></tr><tr><td rowspan="8">2.2 &lt; z &lt; 3.0</td><td>-20.0</td><td>6</td><td>0.108<sup>+0.065</sup><sub>-0.043</sub></td></tr><tr><td>-19.0</td><td>14</td><td>0.252<sup>+0.087</sup><sub>-0.066</sub></td></tr><tr><td>-18.0</td><td>27</td><td>0.488<sup>+0.113</sup><sub>-0.093</sub></td></tr><tr><td>-17.0</td><td>61</td><td>1.186<sup>+0.152</sup><sub>-0.152</sub></td></tr><tr><td>-16.0</td><td>67</td><td>2.170<sup>+0.265</sup><sub>-0.265</sub></td></tr><tr><td>-15.0</td><td>53</td><td>7.003<sup>+0.962</sup><sub>-0.962</sub></td></tr><tr><td>-14.0</td><td>21</td><td>20.915<sup>+5.637</sup><sub>-4.532</sub></td></tr><tr><td>-13.0</td><td>2</td><td>10.271<sup>+13.547</sup><sub>-6.635</sub></td></tr><tr><td colspan="4">LBG LFs</td></tr><tr><td rowspan="6">z &sim; 1.65</td><td>-18.0</td><td>3</td><td>0.527<sup>+0.513</sup><sub>-0.287</sub></td></tr><tr><td>-17.0</td><td>4</td><td>0.733<sup>+0.579</sup><sub>-0.351</sub></td></tr><tr><td>-16.0</td><td>5</td><td>1.329<sup>+0.899</sup><sub>-0.574</sub></td></tr><tr><td>-15.0</td><td>2</td><td>1.150<sup>+1.517</sup><sub>-0.743</sub></td></tr><tr><td>-14.0</td><td>3</td><td>8.566<sup>+8.332</sup><sub>-4.663</sub></td></tr><tr><td>-13.0</td><td>1</td><td>10.477<sup>+24.098</sup><sub>-8.665</sub></td></tr><tr><td rowspan="8">z &sim; 2.0</td><td>-20.0</td><td>3</td><td>0.058<sup>+0.057</sup><sub>-0.032</sub></td></tr><tr><td>-19.0</td><td>10</td><td>0.199<sup>+0.085</sup><sub>-0.062</sub></td></tr><tr><td>-18.0</td><td>34</td><td>0.768<sup>+0.156</sup><sub>-0.131</sub></td></tr><tr><td>-17.0</td><td>50</td><td>1.617<sup>+0.263</sup><sub>-0.228</sub></td></tr><tr><td>-16.0</td><td>40</td><td>3.556<sup>+0.656</sup><sub>-0.560</sub></td></tr><tr><td>-15.0</td><td>20</td><td>5.736<sup>+1.592</sup><sub>-1.271</sub></td></tr><tr><td>-14.0</td><td>13</td><td>20.106<sup>+7.269</sup><sub>-5.504</sub></td></tr><tr><td>-13.0</td><td>2</td><td>25.922<sup>+34.192</sup><sub>-16.746</sub></td></tr><tr><td rowspan="6">z &sim; 2.7</td><td>-20.0</td><td>6</td><td>0.112<sup>+0.067</sup><sub>-0.045</sub></td></tr><tr><td>-19.0</td><td>10</td><td>0.186<sup>+0.079</sup><sub>-0.058</sub></td></tr><tr><td>-18.0</td><td>21</td><td>0.409<sup>+0.110</sup><sub>-0.089</sub></td></tr><tr><td>-17.0</td><td>46</td><td>1.309<sup>+0.223</sup><sub>-0.192</sub></td></tr><tr><td>-16.0</td><td>37</td><td>4.218<sup>+0.814</sup><sub>-0.691</sub></td></tr><tr><td>-15.0</td><td>14</td><td>7.796<sup>+2.690</sup><sub>-2.060</sub></td></tr><tr><td rowspan="2">z &sim; 2.7</td><td>-14.0</td><td>5</td><td>32.106<sup>+21.716</sup><sub>-13.870</sub></td></tr><tr><td>-13.0</td><td>1</td><td>24.187<sup>+55.631</sup><sub>-20.003</sub></td></tr></tbody></table>

The best estimates of Schechter parameters are derived via marginalization of posterior functions at all redshifts. Figure 7 shows the binned LFs along with our best MLE determinations at each redshift range. The best Schechter parameters are tabulated

in Table 3. Our MLE estimates reveal steep faint-end slopes of $\alpha = -1.56 \pm 0.04$, $-1.72 \pm 0.04$, and $-1.94 \pm 0.06$ for the $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$ samples, respectively. We emphasize that our estimate of the faint-end slope at $z \sim 2.6$ is mostly independent of our choice of the $M^*$ prior, as we derive $\alpha = -1.97 \pm 0.06$ in the absence of a prior. Our steep LFs show no sign of turnover down to $M_{UV} = -12.5$ mag.

The contours in Figure 8 illustrate the correlation between the faint-end slope ($\alpha$) and the characteristic magnitude ($M^*$). The best Schechter parameters derived through marginalization are shown with filled blue, green, and red circles for the $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$ samples, respectively. We are also overplotting our best LFs (filled black circles) from the LBG samples (see Section 9). The red dashed contours show the likelihood function at $z \sim 2.6$, before incorporating the $M^*$ prior.

The systematic uncertainties—particularly in the size distribution assumption at faint luminosities—may affect the completeness corrections and thus the LF measurements at these magnitudes. This concern is also expressed in a recent paper by Bouwens et al. (2016), where they measure very small sizes (160–240 pc) for ultra-faint galaxies ($M_{UV} = -15$) at $2 < z < 8$ and then discuss the possible effects due to uncertain size assumptions on the LF measurements. We should emphasize that they present their conclusions for a large redshift range of $z = 2–6$, while we expect the lower redshift galaxies ($z \sim 2$) to be on average larger than their high-redshift counterparts (as seen at higher luminosities; Shibuya et al. 2015). Our assumed size distribution for ultra-faint galaxies is the closest to the Bouwens et al. (2016) measurements, relative to the other LF studies. We run some experiments to investigate whether the faint galaxies with large completeness corrections (i.e., where the systematic uncertainty dominates) are dictating our best-fit LFs by excluding all of our galaxies with completeness below 50%. This reduces the size of our $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$ samples by 33%, 53%, and 44%, such that our final “complete” samples have 186, 127, and 141 galaxies, respectively. To be consistent, we also remove the corresponding volumes from our total volume estimates. We then re-fit the LFs and measure faint-end slopes of $\alpha = -1.55 \pm 0.06$, $\alpha = -1.69 \pm 0.07$, and $\alpha = -1.79 \pm 0.08$ at $1 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively. These estimates are all steep and show the same trend of steeper slopes toward higher redshifts, though with slightly shallower slopes. We note that, although the $z \sim 1.3$ and $z \sim 1.9$ faint-end slopes measured from the “complete” sample are consistent with the slopes measured from the full sample, the $z \sim 2.6$ slope from the “complete” sample is significantly shallower by $\sim 1.5 \times$ the individual errors added in quadrature (although the measurements are not completely independent, so adding in quadrature will slightly overestimate the uncertainties). The probability of obtaining such a deviation in at least one of the three slope measurements is small (10%), and suggests that the systematic uncertainties are not negligible. Consequently, as also emphasized in Bouwens et al. (2016), the size measurements of very faint galaxies will need to be more accurately determined for higher quality LF measurements.

## 9. LF OF LBG SAMPLES

As discussed in Section 1, one of the goals of the present paper is to understand the effect of two widely used selection techniques. To this end, we have also performed a parallel determination of the UV LF based on the Lyman break
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Rest-frame UV luminosity functions at $z \sim 1.3$ (left), $z \sim 1.9$ (middle), and $z \sim 2.6$ (right). The blue, green, and red circles are our binned LFs combining all three lensing clusters (see Section 8.1). The blue and green squares are the LFs from Oesch et al. (2010a) at $1.0 < z < 1.5$ and $1.5 < z < 2.0$, respectively, where the individual data were used for our MLE fitting. The light blue diamonds are the LFs from Reddy & Steidel (2009) based on a BX-selected sample of star-forming galaxies. The gray stars and purple triangles are the results from Parsa et al. (2016; photometric redshift selection) and Sawicki (2012; BX selection), respectively. The solid line in all three panels shows our best Schechter fit through a MLE technique.</figcaption></figure>

“dropout” selection at equivalent redshift ranges. A complete description of our color–color selection, sample contamination, and the completeness simulation for dropouts is given in Appendix A. As explained there, our LBG samples consist of 19 F225W, 178 F275W, and 142 F336W dropouts at $z \sim 1.65$, $z \sim 2.0$, and $z \sim 2.7$, respectively. We note that our LBG samples have fewer galaxies than our photometric redshift samples, because we require a $5\sigma$ detection in the detection filter for these samples (see Appendix A), whereas the photometric redshift samples only require a $3\sigma$ detection (see Section 6). To ensure accurate detection of a break, we restrict our sample to objects where the imaging depth is sufficient to detect at least a one magnitude break (at $1\sigma$) between the dropout filter (F225W, F275W, and F336W at $z \sim 1.65$, 2.0, and 2.7, respectively) compared to the adjacent longer wavelength filter. This cut only removes two galaxies from the A2744 F336W-dropout LF and it does not change the rest of the LBG samples.

The effective volume including the completeness corrections is calculated for these samples using Equation (4). In order to estimate the binned UV LF for our LBG samples, we use the same methodology as we used for our primary photometric redshift samples. Similarly, we restrict our dropout samples to galaxies with $M_{1500} < -12.5$ for the same reasons that were mentioned before (see Section 6). This limit excludes 1 (~5.3%), 6 (~3.4%), and 0 galaxies from the F225W-, F275W-, F336W-dropout LFs, respectively. Finally, we have 18, 172, and 140 galaxies for the $z \sim 1.65$, $z \sim 2.0$, and $z \sim 2.7$ LBG LFs, respectively.

Furthermore, to constrain the bright end of our F225W and F275W-dropout LFs, we incorporate the binned measurements from the Oesch et al. (2010a) LBG samples. Here, we do not use the MLE technique because we do not have individual measurements for all of these bright-end LBG samples. We determine the best Schechter parameters only using the simple $\chi^2$ technique, considering that these two methods of fitting (MLE versus $\chi^2$) show good agreement for the photometric redshift LFs. For our F336W-dropout LF, we only fit to our binned data, keeping the characteristic magnitude $M^*$ at a fixed value of $-20.7$, similar to what we used for our $z \sim 2.6$

photometric redshift LF. The binned values and the best-fitting Schechter parameters for the LBG LFs are given in Tables 2 and 3, respectively.

## 10. DISCUSSION

### 10.1. Comparing the UV LFs of Photometric Redshift and LBG Samples

Figure 9 compares our LF results derived for the photometric redshift and UV-dropout selections. From left to right, the F225W-, F275W- and F336W-dropout LFs shown with black circles are compared with the photometric redshift LFs at $z \sim 1.3$ (blue circles), $z \sim 1.9$ (green circles), and $z \sim 2.6$ (red circles), respectively. Together with our data points for each redshift range, we also show the bright-end LFs of Oesch et al. (2010a) derived from their photometric redshift (blue, green, and red squares) and UV-dropout (black squares) samples. In addition, we include the LF results from several relevant studies (Reddy & Steidel 2009; Sawicki 2012; Parsa et al. 2016). To compare these LF measurements, we run a set of Monte Carlo simulations and estimate the $2\sigma$ confidence interval from each best-fit Schechter function. The gray and hatched regions in Figure 9 encompass the $2\sigma$ uncertainties of the LBG and photometric redshift LFs, respectively. Because our LBG samples have fewer galaxies than the photometric redshift samples, the corresponding LFs are more uncertain. Our LFs are in agreement within these confidence regions. Indeed, similar agreement between the LFs derived from these two selection techniques at higher redshifts has been shown before (McLure et al. 2013; Schenker et al. 2013). However, the lack of a robust knowledge of various systematic effects such as intrinsic size distribution and dust reddening at these faint luminosities still introduces moderate differences between these two LF measurements.

### 10.2. Evolution of the LF Schechter Parameters

In order to understand the evolution of LF parameters with redshift, we compare our best-fit Schechter parameters with other determinations of the rest-frame UV LF at higher and
<table><caption>**Table 3**
Best-fit Schechter Parameters for UV LFs</caption><thead><tr><th>z</th><th>α</th><th>M*</th><th>φ*(10<sup>-3</sup> Mpc<sup>-3</sup> mag<sup>-1</sup>)</th></tr></thead><tbody><tr><td colspan="4">Photometric redshift LF, MLE fitting</td></tr><tr><td>1.0 &lt; z &lt; 1.6<sup>a</sup></td><td>-1.56 ± 0.04</td><td>-19.74 ± 0.18</td><td>2.32 ± 0.49</td></tr><tr><td>1.6 &lt; z &lt; 2.2<sup>a</sup></td><td>-1.72 ± 0.04</td><td>-20.41 ± 0.20</td><td>1.50 ± 0.37</td></tr><tr><td>2.2 &lt; z &lt; 3.0<sup>b</sup></td><td>-1.94 ± 0.06</td><td>-20.71 ± 0.11 (prior)</td><td>0.55 ± 0.14</td></tr><tr><td colspan="4">LBG LF, χ<sup>2</sup> fitting</td></tr><tr><td>z &sim; 1.65<sup>c</sup></td><td>-1.50 ± 0.16</td><td>-19.85 ± 0.41</td><td>2.21 ± 1.32</td></tr><tr><td>z &sim; 2.0<sup>d</sup></td><td>-1.80 ± 0.06</td><td>-20.39 ± 0.31</td><td>1.46 ± 0.65</td></tr><tr><td>z &sim; 2.7<sup>e</sup></td><td>-2.01 ± 0.08</td><td>-20.70(fixed)</td><td>0.48 ± 0.15</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Maximum likelihood fit to the whole sample including individual galaxies from all three lensing clusters as well as the bright-end galaxies from Oesch et al. (2010a).

<sup>b</sup> Maximum likelihood fit to the individual galaxies from the HFF clusters assuming a Gaussian prior for M* (see Section 8.2).

<sup>c</sup> χ<sup>2</sup> fitting to the binned data from A1689 as well as the bright-end LBGs from Oesch et al. (2010a).

<sup>d</sup> χ<sup>2</sup> fitting to the binned data from all three lensing clusters as well as the bright-end LBGs from Oesch et al. (2010a).

<sup>e</sup> χ<sup>2</sup> fitting to the binned data from the HFF clusters assuming a fixed M* (see Section 9).

<figure><img src="image_9.png" /><figcaption>**Figure 8.** 68% and 95% contours of the $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$ photometric redshift LFs shown with blue, green, and red colors, respectively. The red dashed line shows the contours for the $z \sim 2.6$ LF before adding the M* prior (see Section 8.2). The marginalized probability distribution of each parameter $P(\alpha)$ and $P(M^*)$ are also plotted on the right and top sides, respectively. The best-fit values of Schechter parameters for each LF is shown with a filled circle. The black filled circles with error bars denote the best-fit values for the LBG LFs (see Section 9).</figcaption></figure>

lower redshifts in Figure 10. We summarize the evolution of LF parameters as below.

Exploiting the magnification from strong gravitational lensing and consequently extending the UV LF to very low luminosities enables a robust estimate of the faint-end slope. In the context of recent UV LF studies, there are not as many measurements at $z = 1-1.5$ to compare with our estimates, but as can be seen in Figure 10, our inferred value of the faint-end slope ($\alpha = -1.56 \pm 0.04$) at $z \sim 1.3$ is consistent with other results from Arnouts et al. (2005) and Oesch et al. (2010a)

given their large uncertainties. We should note that we are in better agreement with the Oesch et al. (2010a) estimate for their LBG LF (at $z \sim 1.5$), as their photometric redshift LF has a very steep faint-end slope. Regarding our estimate for the $z \sim 1.9$ LF, we are again in good agreement with several other results, particularly both the LBG and photometric redshift LFs of Oesch et al. (2010a) and also the $z = 2.3$ LF from Reddy & Steidel (2009). We note that we are also in agreement with our previous $z \sim 2$ LBG LF from A14 ($\alpha = -1.74 \pm 0.08$). Finally, regarding our estimate for the $z \sim 2.6$ LF, we derive a faint-end slope steeper than previous determinations and more similar to the steep faint-end slopes favored at higher redshifts. As a consequence, we conclude that there was rapid evolution in the faint-end slope toward shallower values during the 2.2 Gyr from $z = 2.6$ to $z = 1.3$, which seems to continue to $z = 0$. We also refer the reader to a recent work by Parsa et al. (2016; see gray filled stars in Figure 10) who study the UV LF between $z = 2-4$. For both $z = 1.9$ and $z = 2.8$, they derive a value of $\alpha = -1.32 \pm 0.04$, which is considerably shallower than most of the other studies, including ours. Consequently, they derive fainter $M^*$ and larger $\phi^*$ values relative to all of the other works at $z = 2-3$ in the literature. We note that they do not use the filter that samples the Lyman break at $z \sim 2$.

In addition to the observed LFs, we compare our results with the LFs from local group (LG) fossil records by Weisz et al. (2014; open green circles in Figure 10). Using the SFHs of LG galaxies, they reconstruct the UV LFs down to very faint magnitudes of $M_{UV} \sim -1.5$. Comparing to our results, they estimate shallower faint-end slopes ($\alpha > -1.4$) for their $z = 1.25$ and 2.0 LFs, but they derive steeper faint-end slopes when they restrict their calculations to the luminosities where their data are complete. For an exact comparison, however, we need to consider the different uncertainties (i.e., small sample size) and systematic errors (i.e., uncertainty in the stellar models used for the SFHs) in their results as well. This discrepancy between the faint-end slopes can be interpreted as a different evolution for the LG dwarfs relative to the field galaxies.

As seen in the middle panel of Figure 10, our characteristic UV magnitude, $M^*$, becomes brighter, $\Delta M^* = -0.7$, at $z = 1.9$ relative to $z = 1.3$. Also, our characteristic number density, $\phi^*$ (lower panel of Figure 10), decreases by a factor of
<figure><img src="image_10.png" /><figcaption>**Figure 9.** Comparing the UV LF of the LBG and photometric redshift samples. The black filled circles are the binned LBG LF for our F225W-, F275W-, and F336W- dropout samples from the left to right panels, respectively. The black dashed line in each panel represents the best-fit Schechter function for the corresponding LBG LF. The gray regions indicate the ±2σ confidence region for each fit derived via Monte Carlo simulation. We overplot the binned and best-fit LFs derived from our photometric redshift samples as shown before in Figure 7. The hatched region denotes the ±2σ confidence region of the best-fit photometric redshift LF. The black open squares are the binned LBG LFs from Oesch et al. (2010a). The rest of the colors and symbols are as in Figure 7. For each redshift range, our binned and best-fit LFs are consistent within the 2σ error bars.</figcaption></figure>

### 10.3. UV Luminosity Density

We use our best LF determinations to derive the comoving UV luminosity density, $\rho_{UV}$, as below:

$$ \rho_{UV} = \int_{L_{\text{faint}}}^{\infty} L\phi(L)dL = \int_{-\infty}^{M_{\text{faint}}} L(M)\phi(M)dM \quad (11) $$

where $\phi(L)(\phi(M))$ is the LF assuming a Schechter function. As an important consequence of the steep faint-end slope of the UV LFs at $1 < z < 3$, the faint star-forming galaxies have a significant contribution to the total unobscured UV luminosity density at these redshifts. To quantify this, we calculate the cumulative UV luminosity density down to various UV luminosity limits. Figure 11 shows these results for our three redshift ranges. We note that all of these calculations are from our photometric redshift LFs, as they have smaller statistical uncertainties. We normalized our cumulative UV luminosity densities to the corresponding value at $M_{UV} = -10$ assuming that there is no turnover in the LF down to this absolute magnitude. To estimate the $1\sigma$ uncertainty at each $M_{UV}$, we run a Gibbs sampler (i.e., Markov Chain Monte Carlo sampling) to obtain a sequence of random pairs of $(\alpha, M^*)$ using their 2D joint probability function and then calculate the distribution of UV luminosity density and the corresponding uncertainty. We also incorporate the Poisson uncertainty in quadrature. These $1\sigma$ uncertainty regions are shaded orange in Figure 11. The unobscured UV luminosity density measurements are tabulated in Table 4. To be consistent with previous studies, we also provide the UV luminosity density values integrated down to $0.04L^*_{z=3}$.

The faint dwarf galaxies with UV magnitudes of $-18.5 < M_{1500} < -12.5$ covered in this work, comprise the majority of the unobscured UV luminosity density at the redshifts of peak star formation activity (58%, 55%, and 59% of the total UV luminosity density at $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$, respectively). Therefore, these dwarf galaxies may contribute significantly to the total intrinsic UV luminosity density and thus to the star formation rate density at these

<figure><img src="image_11.png" /><figcaption>**Figure 10.** Redshift evolution of the Schechter parameters $\alpha$ (top), $M^*$ (middle), and $\phi^*$ (bottom). The present determinations are shown with black filled circles. Our measurements of $M^*$ and $\phi^*$ at $z \sim 2.6$ are shown with black open circles, as they are dependent on our choice of the $M^*$ prior. All symbols from the literature are summarized in the right-hand side of the plot. A detailed description of each parameter evolution is given in the text (see Section 10.2).</figcaption></figure>

1.5 over this time period. However, both of these measurements are dependent on data from other surveys, as our data only sample galaxies fainter than $M^*$.
<figure><img src="image_12.png" /><figcaption>**Figure 11.** Cumulative UV luminosity density at $1.0 < z < 1.6$ (left), $1.6 < z < 2.2$ (middle), and $2.2 < z < 3.0$ (right). The purple, green, and blue dashed lines show the UV-limiting magnitudes for the ERS (Oesch et al. 2010a), Hubble Ultra Deep Field (Parsa et al. 2016), and our samples. The orange region represents the $1\sigma$ uncertainty measured at each $M_{UV}$. We have measured the LF of galaxies responsible for 58%, 55%, and 59% of total UV luminosity density at $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$, respectively.</figcaption></figure>

<table><caption>**Table 4**
UV Luminosity Density<sup>a</sup></caption><thead><tr><th>z</th><th>M &lt; -10</th><th>M &lt; -17.475(0.04L<sup>*</sup><sub>z-3</sub>)<sup>b</sup></th><th>-18.5 &lt; M &lt; -12.5</th></tr></thead><tbody><tr><td>1.0 &lt; z &lt; 1.6</td><td>1.57<sup>+0.08</sup><sub>-0.08</sub> &times; 10<sup>26</sup></td><td>0.90<sup>+0.06</sup><sub>-0.06</sub> &times; 10<sup>26</sup></td><td>0.90<sup>+0.04</sup><sub>-0.05</sub> &times; 10<sup>26</sup></td></tr><tr><td>1.6 &lt; z &lt; 2.2</td><td>2.84<sup>+0.15</sup><sub>-0.15</sub> &times; 10<sup>26</sup></td><td>1.50<sup>+0.09</sup><sub>-0.09</sub> &times; 10<sup>26</sup></td><td>1.57<sup>+0.08</sup><sub>-0.10</sub> &times; 10<sup>26</sup></td></tr><tr><td>2.2 &lt; z &lt; 3.0</td><td>3.13<sup>+0.22</sup><sub>-0.24</sub> &times; 10<sup>26</sup></td><td>1.03<sup>+0.14</sup><sub>-0.19</sub> &times; 10<sup>26</sup></td><td>1.84<sup>+0.13</sup><sub>-0.15</sub> &times; 10<sup>26</sup></td></tr></tbody></table>

**Notes.**

<sup>a</sup> These values are in units of erg s<sup>-1</sup> Hz<sup>-1</sup> Mpc<sup>-3</sup>.

<sup>b</sup> For $L_{z=3}^*$ we use the measurement from Reddy & Steidel (2009).

However, to quantify this, we need to incorporate the effect of dust reddening and its dependence on galaxy luminosity.

In order to understand the evolution of the UV luminosity density, we compare our $\rho_{UV}$ measurements with other studies at various redshifts. As the value of $\rho_{UV}$ depends on the limiting luminosity, i.e., $L_{faint}$ in Equation (11), we use the published Schechter parameters from the literature and calculate the UV luminosity densities and corresponding uncertainties by integrating down to the same $M_{faint} = -10$. We should note that there is no straightforward way to estimate the $\rho_{UV}$ uncertainties as necessary information for these measurements such as covariance between Schechter parameters are not usually provided in the literature. But to assign uncertainty to each $\rho_{UV}$ measurement in a consistent way, we use the same methodology as Madau & Dickinson (2014). We assume that the fractional error, i.e., $\Delta\rho_{UV}/\rho_{UV}$, provided by each author is fixed and thus derive the corresponding uncertainty for our $\rho_{UV}$ value with $M_{faint} = -10$. Figure 12 illustrates these measurements. As seen in many previous studies (e.g., Cucciati et al. 2012), the unobscured (i.e., uncorrected for dust extinction) UV luminosity density rises from $z = 0$ to $z = 2.0$ where it reaches its peak and starts to decline after $z = 3$ (e.g., Oesch et al. 2010a; Finkelstein

et al. 2015). As shown in Figure 12, our $\rho_{UV}$ points (black filled circles) follow the similar trend as seen by previous determinations. However, our measurements show a more rapid evolution from $z = 1.3$ to $z = 1.9$ followed by a slower evolution up to $z = 2.6$.

We emphasize that the unobscured $\rho_{UV}$ evolution rate and the exact location of the peak depends on the wavelength (Trenti et al. 2012) where $\rho_{UV}$ is being measured, and the limiting luminosity, i.e., $L_{faint}$ in Equation (11). Therefore, to find the best-fitting function describing the evolution of unobscured UV luminosity density between $z = 0-2.6$, we only include the results from the literature at the same wavelength (1500 Å) and integrated down to the same magnitude of $M_{UV} = -10$ through our own compilation. Fitting a power law, we find $\rho_{UV} = 25.58 \times (1+z)^{1.74}$, incorporating the data points from Schiminovich et al. (2005; red filled circle), Dahlen et al. (2007; pink open diamond), Oesch et al. (2010a; for the photometric redshift sample, orange filled square), and Cucciati et al. (2012; blue filled triangle).

In addition, to study the evolution of $\rho_{UV}$ for the whole redshift range from $z = 0$ to $z = 8$, we fit a function used by Madau & Dickinson (2014) as shown below. For the higher redshifts, we incorporate the data points from McLure et al. (2009, 2013; green open diamond), Bouwens et al. (2015;
<figure><img src="image_13.png" /><figcaption>Figure 12. Redshift evolution of the unobscured UV luminosity density measured at rest-frame wavelength of 1500 Å. To estimate the $\rho_{UV}$ values, the LFs are integrated down to $M_{UV} = -10$ at all redshifts. The uncertainty for each data point is derived by retaining the fractional error of published $\rho_{UV}$ values from each author. The black filled circles are derived from our photometric redshift LF estimates. Similar to Figure 10, we show our $\rho_{UV}$ measurement at $z = 2.6$ with a black open circle as it depends on the choice of $M^*$ prior. The rest of the symbols are similar to Figure 10, except the red filled circles which are from Schiminovich et al. (2005) using the LF estimates from Arnouts et al. (2005) and Wyder et al. (2005) (shown with red filled circle and blue asterisk in Figure 10, respectively). The dashed line indicates the best-fitting power law to the data points at $z < 2.6$. The solid line shows the best-fitting function (see Equation (12)) for the redshift range of $0 < z < 8$.</figcaption></figure>

green filled star) and Parsa et al. (2016; gray filled star), as well as the data points that we used for the power law.

$$ \rho_{UV}(z) = a \frac{(1+z)^b}{1 + [(1+z)/d]^c} \text{ erg s}^{-1} \text{ Hz}^{-1} \text{ Mpc}^{-3} \quad (12) $$

where we derive $a = 0.34 \pm 0.04$, $b = 2.14 \pm 0.27$, $c = 3.41 \pm 0.23$ and $d = 3.86 \pm 0.63$. We emphasize that these best-fit values describe the $\rho_{UV}$ evolution assuming a limiting magnitude of $M_{UV} = -10$, dramatically fainter than typical limits used in previous studies (~-17.5; Madau & Dickinson 2014). Because we do not account for an increase in the uncertainty of $\rho_{UV}$ at low luminosities, we add in quadrature 12% uncertainty to all of the data points to keep the reduced chi-squared equal to one.

## 10.4. No Turnover in the UV LF

Our observations have now reached the very faint luminosities where some simulations predict a turnover in the UV LF. Though our steep LFs extend down to $M_{UV} = -12.5$ and we showed that the faint bins with large completeness corrections are not affecting the faint-end slope fit (see Section 8.2), they may affect our interpretation of whether or not there is a turnover. As seen in Figure 7, we can rule out the possibility of a turnover in the LF at magnitudes brighter than $M_{UV} < -14$, because one would need an unphysically large systematic effect to cause a turnover at this magnitude. This conclusion is in conflict with the results of the recent cosmological hydrodynamical simulation by Kuhlen et al. (2013), who predict a turnover at $M_{UV} = -16$ in the $z \sim 2.5$ UV LF. Implementing an H$_2$-regulated star formation model, they predict that the star formation is suppressed in dwarf galaxies ($M_{UV} > -16$)

because their gas surface density is below what is required to build a substantial molecular fraction. A similar tension between the observed UV LF and the turnover predicted by recent theoretical results has also been seen at higher redshifts (e.g., Jaacks et al. 2013; Livermore et al. 2016).

The presence of a turnover in the UV LF might also be used to constrain warm dark matter (WDM) models. Menci et al. (2016) provide a limit on the WDM particle mass by comparing the WDM halo mass function and the number density of ultra-faint galaxies derived from the UV LF in A14. The constraints can now be significantly improved given the much larger sample in this survey.

## 11. CONCLUSION

We have obtained deep near-UV imaging of three lensing clusters, two from the HFF surveys (A2744 and MACS J0717) and A1689, to study the evolution of the UV LF during the peak epoch of cosmic star formation at $1 < z < 3$. Combining deep data with strong gravitational lensing magnification, we obtain a large sample (780) of ultra-faint galaxies with $M_{UV} < -12.5$ at $1 < z < 3$ using the photometric redshift selection. We perform an extensive set of simulations to compute the completeness correction required for the LF measurements. We summarize our conclusions below:

1. We derive the best Schechter fit to each UV LF using an MLE technique considering various sources of uncertainty, including the lensing models. Thanks to the lensing magnification, we can extend the UV LF measurements down to very faint luminosities of $M_{UV} = -12.5$ at $1 < z < 3$. Consequently, we find a robust estimate of the UV LF faint-end slope to be $\alpha = -1.56 \pm 0.04$, $\alpha = -1.72 \pm 0.04$, and $\alpha = -1.94 \pm 0.06$ at $1.0 < z < 1.6$, $1.6 < z < 2.2$, and $2.2 < z < 3.0$, respectively. Our $\alpha$ measurements at $z \sim 1.3$ and $z \sim 1.9$ are consistent with previous studies of Reddy & Steidel (2009) and Oesch et al. (2010a). But for $z \sim 2.3$, we have a steeper faint-end slope than previous studies. Our determinations of the UV LFs show a rapid evolution in the faint-end slope toward steeper values from $z = 1.3$ to $z = 2.6$. In addition, when we run a test to minimize the systematic effects by excluding galaxies and volumes with $<50\%$ completeness, we still derive steep faint-end slopes of $\alpha = -1.55 \pm 0.06$, $-1.69 \pm 0.07$, and $-1.79 \pm 0.08$ at $z \sim 1.3$, $1.9$, and $2.6$, respectively. However, for a better determination of the LF parameters, we need a better understanding of the size and color distribution of these faint galaxies.

2. To understand the effect of different selection techniques on the UV LF, we use a two-color “dropout” selection of Lyman break galaxies at redshifts similar to our photometric redshift samples. After correcting for incompleteness and then finding the best-fit Schechter parameters, our LBG and photometric redshift LFs are in $2\sigma$ agreement.

3. We integrate our UV LFs down to a magnitude limit of $M_{UV} = -10$ and find the UV luminosity density to be $\rho_{UV} = 1.57_{-0.08}^{+0.08} \times 10^{26}$, $2.84_{-0.15}^{+0.15} \times 10^{26}$, and $3.13_{-0.24}^{+0.22} \times 10^{26}$ erg s$^{-1}$ Hz$^{-1}$ Mpc$^{-3}$ at $z \sim 1.3$, $z \sim 1.9$, and $z \sim 2.6$, respectively. We show that the faint star-forming galaxies with $-18.5 < M_{UV} < -12.5$
<figure><img src="image_14.png" /><figcaption>**Figure 13.** Color-color selection of Lyman break galaxies for F225W dropouts (left), F275W dropouts (middle), and F336W dropouts (right). A1689 and HFF LBG selections are shown separately in the upper and lower panels, respectively. Gray circles are all of the objects in the corresponding catalogs. The orange asterisks are stars from Pickles (1998). The green dashed line shows the color track of low-redshift ($0 < z < 1$) elliptical galaxies from Coleman et al. (1980). The violet lines are star-forming tracks with different dust reddening values of $E(B-V) = [0, 0.1, 0.2, 0.3]$. The blue region shows the selection criteria. The blue symbols are the LBG candidates with $5\sigma$ detection in two bands redward of the Lyman limit (see text). The blue arrows represent $1\sigma$ lower limits.</figcaption></figure>

contribute the majority of the total unobscured UV luminosity density during the peak epoch of cosmic star formation.

4. Though some models of warm dark matter and some prescriptions for H₂-regulated star formation predict a turnover in the UV LF, we see no evidence for a turnover down to $M_{UV} = -14$ at $1 < z < 3$.

This study highlights the power of gravitational lensing to produce a robust constraint on the faint-end of the LF. However, as mentioned in Section 10.1, this analysis still

suffers from uncertainties that are systematic, rather than statistical. To overcome these uncertainties, in the future, we require precise measurements of size distribution and dust reddening at low luminosities.

We thank the referee for a careful reading and useful comments that improved this paper. The authors are grateful to the STScI and HFF teams for obtaining and reducing the *HST* images. A.A. would like to thank Pascal Oesch for providing his individual measurements for photometric redshift LFs, Takatoshi Shibuya for sending size measurements, Jose Diego
<figure><img src="image_15.png" /><figcaption>**Figure 14.** Photometric redshift distribution of LBGs. The gray histogram shows the photometric redshift distribution of all objects in the fields. The blue histogram highlights the photometric redshift distribution of those galaxies selected as LBGs. The solid and dashed lines show the simulated completeness distribution over redshift for a galaxy with $\mu_{mag} = 1.0$ and $m_{UV} = 27$ and 28, respectively. The purple, blue, and orange colors present the completeness distribution for A1689, A2744, and MACS J0717, respectively. The right-hand axis shows the completeness values.</figcaption></figure>

for providing us with a list of A1689 multiple images, as well as John Blakeslee and Karla Kalamo for providing us with a list of globular clusters of A1689. A.A. also thanks Daniel Weisz for his valuable comments as well as Naveen Reddy, Nader Shakibay Senobari, Mario De Leo, Ali Ahmad Khostovan, and Kaveh Vasei for useful conversations. M.J. acknowledges support from the Science and Technology Facilities Council (grant number ST/L00075X/1 & ST/F001166/1). M.L. acknowledges the Centre National de la Recherche Scientifique (CNRS) for its support. This work is based on observations with the NASA/ESA Hubble Space Telescope, obtained at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555.

## APPENDIX A
### LYMAN BREAK SELECTION

In this section, we outline our selection criteria to find Lyman break galaxies (Steidel et al. 1999). We adopted a standard color-color diagram to sample the UV continuum break in the SED of high-redshift galaxies. As shown in Figure 13, the selection region is defined based on the location of tracks of star-forming galaxies in the color-color plot. The star-forming tracks are predictions from Bruzual & Charlot (2003) synthetic stellar population models assuming a constant star formation history, 0.2 $Z_\odot$, and an age of 100 Myr with different color excess of $E(B-V) = [0, 0.1, 0.2, 0.3]$. In the next subsection, we summarize the selection criteria we use to identify the $z \sim 1-3$ galaxies. As A1689 is
<table><caption>**Table 5**
New Multiply Imaged Systems in A1689</caption><thead><tr><th>System</th><th>R.A.(J2000)</th><th>Decl.(J2000)</th><th>z<sub>phot</sub></th><th>Δz<sub>phot</sub></th><th>F336W–F475W</th><th>F475W–F625W</th><th>F625W–F775W</th><th>F775W–F814W</th></tr></thead><tbody><tr><td>a.1</td><td>197.875427</td><td>-1.353059</td><td>2.33</td><td>0.13</td><td>0.72 ± 0.12</td><td>-0.19 ± 0.12</td><td>-0.05 ± 0.15</td><td>-0.01 ± 0.12</td></tr><tr><td>a.2</td><td>197.886612</td><td>-1.344707</td><td>1.90</td><td>0.07</td><td>0.30 ± 0.09</td><td>-0.09 ± 0.10</td><td>-0.08 ± 0.12</td><td>0.03 ± 0.10</td></tr><tr><td>a.3</td><td>197.858505</td><td>-1.337834</td><td>2.33</td><td>0.10</td><td>0.67 ± 0.11</td><td>-0.11 ± 0.11</td><td>-0.12 ± 0.13</td><td>-0.02 ± 0.11</td></tr><tr><td>a.4</td><td>197.880814</td><td>-1.335316</td><td>2.20</td><td>0.17</td><td>0.65 ± 0.15</td><td>-0.08 ± 0.14</td><td>-0.01 ± 0.17</td><td>-0.13 ± 0.14</td></tr><tr><td>b.1</td><td>197.886612</td><td>-1.352190</td><td>1.87</td><td>0.27</td><td>0.48 ± 0.29</td><td>-0.21 ± 0.33</td><td>-0.06 ± 0.41</td><td>-0.08 ± 0.33</td></tr><tr><td>b.2</td><td>197.873810</td><td>-1.333846</td><td>2.02</td><td>0.10</td><td>-0.03 ± 0.18</td><td>-0.36 ± 0.28</td><td>0.06 ± 0.34</td><td>-0.13 ± 0.27</td></tr><tr><td>c.1</td><td>197.875732</td><td>-1.350983</td><td>2.47</td><td>0.10</td><td>1.56 ± 0.22</td><td>0.07 ± 0.11</td><td>0.02 ± 0.12</td><td>-0.07 ± 0.10</td></tr><tr><td>c.2</td><td>197.858337</td><td>-1.333123</td><td>2.40</td><td>0.06</td><td>1.30 ± 0.17</td><td>0.23 ± 0.10</td><td>0.17 ± 0.10</td><td>-0.21 ± 0.08</td></tr><tr><td>d.1</td><td>197.879654</td><td>-1.342635</td><td>2.61</td><td>0.10</td><td>1.92 ± 0.30</td><td>0.01 ± 0.12</td><td>-0.16 ± 0.14</td><td>0.01 ± 0.12</td></tr><tr><td>d.2</td><td>197.878067</td><td>-1.342786</td><td>2.57</td><td>0.19</td><td>2.08 ± 0.54</td><td>0.02 ± 0.17</td><td>-0.02 ± 0.19</td><td>0.15 ± 0.15</td></tr><tr><td>d.3</td><td>197.855209</td><td>-1.339240</td><td>2.30</td><td>0.09</td><td>0.62 ± 0.12</td><td>0.11 ± 0.11</td><td>-0.04 ± 0.12</td><td>-0.19 ± 0.10</td></tr><tr><td>e.1</td><td>197.877899</td><td>-1.354296</td><td>3.19</td><td>0.19</td><td>>3.32</td><td>0.67 ± 0.09</td><td>0.10 ± 0.08</td><td>-0.07 ± 0.07</td></tr><tr><td>e.2</td><td>197.885437</td><td>-1.349361</td><td>3.23</td><td>0.18</td><td>>3.65</td><td>0.64 ± 0.08</td><td>0.05 ± 0.07</td><td>0.05 ± 0.06</td></tr><tr><td>e.3</td><td>197.879929</td><td>-1.335399</td><td>3.11</td><td>0.26</td><td>>2.39</td><td>0.76 ± 0.17</td><td>0.13 ± 0.13</td><td>0.04 ± 0.10</td></tr><tr><td>e.4</td><td>197.857147</td><td>-1.340662</td><td>3.32</td><td>0.23</td><td>>2.86</td><td>0.89 ± 0.13</td><td>0.02 ± 0.10</td><td>0.08 ± 0.08</td></tr></tbody></table>

observed with different sets of filters than the HFFs, it is not possible to use the same color criteria for all of these fields. Therefore, we construct analogous selection criteria as below.

For F225W dropout sources, considering that A1689 is the only field where F225W images are available, the selection criteria are as below:

$$ 
\begin{aligned}
& \text{F225W} - \text{F275W} > 0.75 \\
& \text{F275W} - \text{F336W} < 1.4 \\
& \text{F225W} - \text{F275W} > 1.67(\text{F275W} - \text{F336W}) - 0.42 \\
& \text{S/N(F275W)} > 5, \text{S/N(F336W)} > 5.
\end{aligned} 
\qquad (13) $$

These color criteria, which are identical to Oesch et al. (2010a), find 31 galaxy candidates in A1689.

For F275W dropout sources, the selection criteria for A1689 with F625W-band imaging are

$$ 
\begin{aligned}
& \text{F275W} - \text{F336W} > 1 \\
& \text{F336W} - \text{F625W} < 1.3 \\
& \text{F275W} - \text{F336W} > 2.67(\text{F336W} - \text{F625W}) - 1.67 \\
& \text{S/N(F336W)} > 5, \text{S/N(F625W)} > 5.
\end{aligned} 
\qquad (14) $$

These color criteria, which are identical to what we used before in A14, find 99 galaxy candidates in A1689. For F275W dropout sources, the selection criteria for HFFs where F606W-band imaging are available instead of F625W-band data, we use selection criteria identical to A1689. In total, these color criteria find 230 galaxy candidate over three clusters. Ninety-nine of these candidates are from A1689 in comparison with 84 candidates in A14, because we added 14 orbits to the 4 orbits of the F336W image that we used in A14.

For F336W dropout sources, the selection criteria for HFFs are

$$ 
\begin{aligned}
& \text{F336W} - \text{F435W} > 1 \\
& \text{F435W} - \text{814W} < 1.2 \\
& \text{F336W} - \text{F435W} > 2.4(\text{F435W} - \text{F814W}) - 0.68 \\
& \text{S/N(F435W)} > 5, \text{S/N(F814W)} > 5.
\end{aligned} 
\qquad (15) $$

These criteria find 189 galaxy candidates over HFFs. Similar to our $z \sim 2.6$ photometric redshift sample (see Section 6), we do not include A1689 in our F336W-dropout sample as there is contamination from cluster members. As discussed in Section 6, we remove all of the multiple images corresponding to a single source except the brightest image. We then use the same identification scheme as our photometric redshift samples to remove contamination from stars, stellar spikes, and spurious detections from low-redshift bright galaxies. As seen in Figure 13, the stellar sequence (orange asterisks) enters the selection region of the F225W and F336W dropouts, resulting in a contamination of 3.2% and 1.6% of stars, respectively. We also excluded a low contamination of 1.3% and 4.2% from the stellar spikes and spurious objects in the F275W and F336W dropouts, respectively.

In addition, our photometric redshift measurements show that the fraction of low-redshift interlopers in the LBG samples is low. We find that only 9.7%, 7.8%, and 5.8% of our F225W-, F275W-, and F336W-dropout samples are low-redshift interloper with $z < 1.0$, $z < 1.3$, and $z < 1.5$, respectively. As illustrated in Figure 14 and explained in the next section, we derive these redshift cuts using the expected redshift distribution from our completeness simulations for dropout samples.

## APPENDIX B
COMPLETENESS SIMULATION FOR LBG LF

Following the analysis for our photometric redshift sample, we run the same Monte Carlo simulation to assess the completeness values, $C(m, z, \mu_{mag})$, for the Lyman break samples. As described in detail in Section 7, after generating random artificial galaxies with similar properties to observed sources, we require the same selection criteria (Equations (13)–(15)) as we used for the observed LBGs. Figure 14 illustrates the completeness distribution for two intrinsic apparent magnitudes (i.e., before magnification) of $m_{UV} = 27, 28$ and magnification $\mu_{mag} = 1.0$ mag. To compare with the observed galaxies, we overplot the photometric redshift distribution of the whole catalog together with the subsample selected as LBGs. As seen in this figure, the redshift distribution of LBGs (blue histograms) relative to the redshift distribution of all
<figure><img src="image_16.png" /><figcaption>**Figure 15.** Newly identified multiply imaged systems in A1689. The arrow in each stamp denotes the position of each multiple image. The color image is a combination of the F475W (blue), F625W (green), and F814W (red) filters. We note that the reddening in the d.2 image is due to a nearby cluster member. The size of each cutout is 2".</figcaption></figure>

galaxies in the field (gray histograms) is consistent with the completeness calculations.

## APPENDIX C
NEW MULTIPLY IMAGED SYSTEMS

As described in Section 5.3, we identify five new multiply lensed system candidates in A1689. Table 5 summarizes these

systems, where we provide their photometric redshift estimates and color measurements as well. Because one of the primary indicators of multiple images is their uniform colors (i.e., magnification is independent of wavelength), we show their RGB composite image (see Figure 15) combining F814W, F625W, and F475W data as red, green, and blue filters. In Figure 16, we show the positions of all five new
<figure><img src="image_17.png" /><figcaption>**Figure 16.** The color image is a composite of the F475W (blue), F625W (green), and F814W (red) filters. The contours show the critical lines for sources at z = 2.5. The circles denote the positions of the newly found multiple images. A compass provides the orientation, and the lengths of the arrows show the 15" scale. Some of the labels have been offset slightly.</figcaption></figure>

systems on a color-composite image of A1689. We also overplot the critical lines at z = 2.5. The RGB colors are similar to Figure 15.

## REFERENCES

* Adelberger, K. L., & Steidel, C. C. 2000, ApJ, 544, 218

* Alamo-Martínez, K. A., Blakeslee, J. P., Jee, M. J., et al. 2013, ApJ, 775, 20

* Alavi, A., Siana, B., Richard, J., et al. 2014, ApJ, 780, 143

* Arnouts, S., Schiminovich, D., Ilbert, O., et al. 2005, ApJL, 619, L43

* Atek, H., Richard, J., Jauzac, M., et al. 2015a, ApJ, 814, 69

* Atek, H., Richard, J., Kneib, J.-P., et al. 2014, ApJ, 786, 60

* Atek, H., Richard, J., Kneib, J.-P., et al. 2015b, ApJ, 800, 18

* Bartelmann, M. 2010, CQGra, 27, 233001

* Beckwith, S. V. W., Stiavelli, M., Koekemoer, A. M., et al. 2006, AJ, 132, 1729

* Benson, A. J., Bower, R. G., Frenk, C. S., et al. 2003, ApJ, 599, 38

* Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393

* Blanton, M. R., & Roweis, S. 2007, AJ, 133, 734

* Bouwens, R. J., Illingworth, G. D., Blakeslee, J. P., Broadhurst, T. J., & Franx, M. 2004, ApJL, 611, L1

* Bouwens, R. J., Illingworth, G. D., Franx, M., & Ford, H. 2007, ApJ, 670, 928

* Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2015, ApJ, 803, 34

* Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2016, arXiv:1608.00966

* Bowler, R. A. A., Dunlop, J. S., McLure, R. J., et al. 2015, MNRAS, 452, 1817

* Bradač, M., Schneider, P., Lombardi, M., & Erben, T. 2005, A&A, 437, 39

* Bradley, L. D., Trenti, M., Oesch, P. A., et al. 2012, ApJ, 760, 108

* Brammer, G. B., van Dokkum, P. G., & Coppi, P. 2008, ApJ, 686, 1503

* Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000

* Budavári, T., Szalay, A. S., Charlot, S., et al. 2005, ApJL, 619, L31

* Bunker, A. J., Stanway, E. R., Ellis, R. S., & McMahon, R. G. 2004, MNRAS, 355, 374

* Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682

* Casertano, S., de Mello, D., Dickinson, M., et al. 2000, AJ, 120, 2747

* Coleman, G. D., Wu, C.-C., & Weedman, D. W. 1980, ApJS, 43, 393

* Cucciati, O., Tresse, L., Ilbert, O., et al. 2012, A&A, 539, A31

* Dahlen, T., Mobasher, B., Dickinson, M., et al. 2007, ApJ, 654, 172

* Dahlen, T., Mobasher, B., Dickinson, M., et al. 2010, ApJ, 724, 425

* Dekel, A., & Birnboim, Y. 2006, MNRAS, 368, 2

* Dickinson, M., Stern, D., Giavalisco, M., et al. 2004, ApJL, 600, L99

* Diego, J. M., Tegmark, M., Protopapas, P., & Sandvik, H. B. 2007, MNRAS, 375, 958

* Diolaiti, E., Bendinelli, O., Bonaccini, D., et al. 2000, in ASP Conf. Ser. 216, Astronomical Data Analysis Software and Systems IX, ed. N. Manset, C. Veillet, & D. Crabtree (San Francisco, CA: ASP), 623

* Ebeling, H., Ma, C.-J., & Barrett, E. 2014, ApJS, 211, 21

* Ferguson, H. C., Dickinson, M., Giavalisco, M., et al. 2004, ApJL, 600, L107

* Finkelstein, S. L., Ryan, R. E., Jr., Papovich, C., et al. 2015, ApJ, 810, 71

* Fioc, M., & Rocca-Volmerange, B. 1997, A&A, 326, 950

* Gehrels, N. 1986, ApJ, 303, 336

* Gonzaga, S., Hack, W., Fruchter, A., & Mack, J. 2012, The DrizzlePac Handbook (Baltimore, MD: STScI)

* Grazian, A., Castellano, M., Koekemoer, A. M., et al. 2011, A&A, 532, A33

* Grazian, A., Fontana, A., de Santis, C., et al. 2006, A&A, 449, 951

* Haberzettl, L., Williger, G. M., Lauroesch, J. T., et al. 2009, ApJ, 702, 506

* Hathi, N. P., Cohen, S. H., Ryan, R. E., Jr., et al. 2013, ApJ, 765, 88

* Hathi, N. P., Ryan, R. E., Jr., Cohen, S. H., et al. 2010, ApJ, 720, 1708

* Huang, K.-H., Ferguson, H. C., Ravindranath, S., & Su, J. 2013, ApJ, 765, 68

* Ilbert, O., Arnouts, S., McCracken, H. J., et al. 2006, A&A, 457, 841

* Ishigaki, M., Kawamata, R., Ouchi, M., et al. 2015, ApJ, 799, 12

* Iwata, I., Ohta, K., Tamura, N., et al. 2007, MNRAS, 376, 1557

* Jaacks, J., Thompson, R., & Nagamine, K. 2013, ApJ, 766, 94

* Jauzac, M., Richard, J., Jullo, E., et al. 2015, MNRAS, 452, 1437

* Jones, T., Ellis, R., Jullo, E., & Richard, J. 2010, ApJL, 725, L176

* Jullo, E., & Kneib, J.-P. 2009, MNRAS, 395, 1319

* Jullo, E., Kneib, J.-P., Limousin, M., et al. 2007, NJPh, 9, 447

* Kinney, A. L., Calzetti, D., Bohlin, R. C., et al. 1996, ApJ, 467, 38

* Kneib, J.-P., & Natarajan, P. 2011, A&ARv, 19, 47

* Kuhlen, M., & Faucher-Giguère, C.-A. 2012, MNRAS, 423, 862

* Kuhlen, M., Madau, P., & Krumholz, M. R. 2013, ApJ, 776, 34

* Law, D. R., Steidel, C. C., Shapley, A. E., et al. 2012, ApJ, 759, 29

* Liesenborgs, J., De Rijcke, S., & Dejonghe, H. 2006, MNRAS, 367, 1209

* Limousin, M., Ebeling, H., Richard, J., et al. 2012, A&A, 544, A71

* Limousin, M., Richard, J., Jullo, E., et al. 2007, ApJ, 668, 643

* Limousin, M., Richard, J., Jullo, E., et al. 2016, A&A, 588, A99

* Livermore, R. C., Finkelstein, S. L., & Lotz, J. M. 2016, arXiv:1604.06799

* Lotz, J. M., Koekemoer, A., Coe, D., et al. 2016, arXiv:1605.06567

* Ly, C., Malkan, M. A., Treu, T., et al. 2009, ApJ, 697, 1410
Madau, P. 1995, ApJ, 441, 18

Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415

Madau, P., Ferrara, A., & Rees, M. J. 2001, ApJ, 555, 92

McLure, R. J., Cirasuolo, M., Dunlop, J. S., Foucaud, S., & Almaini, O. 2009, MNRAS, 395, 2196

McLure, R. J., Dunlop, J. S., Bowler, R. A. A., et al. 2013, MNRAS, 432, 2696

McLure, R. J., Dunlop, J. S., de Ravel, L., et al. 2011, MNRAS, 418, 2074

Mehta, V., Scarlata, C., Colbert, J. W., et al. 2015, ApJ, 811, 141

Menci, N., Sanchez, N. G., Castellano, M., & Grazian, A. 2016, ApJ, 818, 90

Merten, J., Cacciato, M., Meneghetti, M., Mignone, C., & Bartelmann, M. 2009, A&A, 500, 681

Merten, J., Coe, D., Dupke, R., et al. 2011, MNRAS, 417, 333

Mosleh, M., Williams, R. J., Franx, M., et al. 2012, ApJL, 756, L12

Narayan, R., Blandford, R., & Nityananda, R. 1984, Natur, 310, 112

Nestor, D. B., Shapley, A. E., Kornei, K. A., Steidel, C. C., & Siana, B. 2013, ApJ, 765, 47

Oesch, P. A., Bouwens, R. J., Carollo, C. M., et al. 2010a, ApJL, 725, L150

Oesch, P. A., Bouwens, R. J., Illingworth, G. D., et al. 2010b, ApJL, 709, L16

Oesch, P. A., Bouwens, R. J., Illingworth, G. D., et al. 2015, ApJ, 808, 104

Oguri, M. 2010, PASJ, 62, 1017

Oke, J. B., & Gunn, J. E. 1983, ApJ, 266, 713

Ono, Y., Ouchi, M., Curtis-Lake, E., et al. 2013, ApJ, 777, 155

Ouchi, M., Mobasher, B., Shimasaku, K., et al. 2009, ApJ, 706, 1136

Ouchi, M., Shimasaku, K., Okamura, S., et al. 2004, ApJ, 611, 660

Owers, M. S., Randall, S. W., Nulsen, P. E. J., et al. 2011, ApJ, 728, 27

Parsa, S., Dunlop, J. S., McLure, R. J., & Mortlock, A. 2016, MNRAS, 456, 3194

Pettini, M., Rix, S. A., Steidel, C. C., et al. 2002, Ap&SS, 281, 461

Pickles, A. J. 1998, yCat, 611, 863

Porciani, C., & Madau, P. 2005, ApJL, 625, L43

Postman, M., Coe, D., Benítez, N., et al. 2012, ApJS, 199, 25

Priewe, J., Williams, L. L. R., Liesenborgs, J., Coe, D., & Rodney, S. A. 2016, arXiv:1605.07621

Rafelski, M., Teplitz, H. I., Gardner, J. P., et al. 2015, AJ, 150, 31

Reddy, N. A., & Steidel, C. C. 2009, ApJ, 692, 778

Reddy, N. A., Steidel, C. C., Pettini, M., et al. 2008, ApJS, 175, 48

Rees, M. J., & Ostriker, J. P. 1977, MNRAS, 179, 541

Richard, J., Jauzac, M., Limousin, M., et al. 2014, MNRAS, 444, 268

Robertson, B. E., Furlanetto, S. R., Schneider, E., et al. 2013, ApJ, 768, 71

Sandage, A., Tammann, G. A., & Yahil, A. 1979, ApJ, 232, 352

Sawicki, M. 2012, MNRAS, 421, 2187

Sawicki, M., & Thompson, D. 2006, ApJ, 642, 653

Schechter, P. 1976, ApJ, 203, 297

Schenker, M. A., Robertson, B. E., Ellis, R. S., et al. 2013, ApJ, 768, 196

Schiminovich, D., Ilbert, O., Arnouts, S., et al. 2005, ApJL, 619, L47

Schlafly, E. F., & Finkbeiner, D. P. 2011, ApJ, 737, 103

Schmidt, K. B., Treu, T., Brammer, G. B., et al. 2014a, ApJL, 782, L36

Schmidt, K. B., Treu, T., Trenti, M., et al. 2014b, ApJ, 786, 57

Sendra, I., Diego, J. M., Broadhurst, T., & Lazkoz, R. 2014, MNRAS, 437, 2642

Shibuya, T., Ouchi, M., & Harikane, Y. 2015, ApJS, 219, 15

Siana, B., Polletta, M. d. C., Smith, H. E., et al. 2008a, ApJ, 675, 49

Siana, B., Smail, I., Swinbank, A. M., et al. 2009, ApJ, 698, 1273

Siana, B., Teplitz, H. I., Chary, R.-R., Colbert, J., & Frayer, D. T. 2008b, ApJ, 689, 59

Somerville, R. S., Hopkins, P. F., Cox, T. J., Robertson, B. E., & Hernquist, L. 2008, MNRAS, 391, 481

Soucail, G., Fort, B., Mellier, Y., & Picat, J. P. 1987, A&A, 172, L14

Stark, D. P., Swinbank, A. M., Ellis, R. S., et al. 2008, Natur, 455, 775

Steidel, C. C., Adelberger, K. L., Giavalisco, M., Dickinson, M., & Pettini, M. 1999, ApJ, 519, 1

Teplitz, H. I., Rafelski, M., Kurczynski, P., et al. 2013, AJ, 146, 159

Trenti, M., Perna, R., Levesque, E. M., Shull, J. M., & Stocke, J. T. 2012, ApJL, 749, L38

Treu, T., Schmidt, K. B., Brammer, G. B., et al. 2015, ApJ, 812, 114

van der Burg, R. F. J., Hildebrandt, H., & Erben, T. 2010, A&A, 523, A74

Vargas, C. J., Bish, H., Acquaviva, V., et al. 2014, ApJ, 783, 26

Vasei, K., Siana, B., Shapley, A. E., et al. 2016, ApJ, 831, 38

Wang, X., Hoag, A., Huang, K.-H., et al. 2015, ApJ, 811, 29

Weisz, D. R., Johnson, B. D., & Conroy, C. 2014, ApJL, 794, L3

White, S. D. M., & Rees, M. J. 1978, MNRAS, 183, 341

Windhorst, R. A., Cohen, S. H., Hathi, N. P., et al. 2011, ApJS, 193, 27

Wyder, T. K., Treyer, M. A., Milliard, B., et al. 2005, ApJL, 619, L15

Yan, H., & Windhorst, R. A. 2004, ApJL, 612, L93

Yoshida, M., Shimasaku, K., Kashikawa, N., et al. 2006, ApJ, 653, 988

Yuan, T.-T., Kewley, L. J., & Richard, J. 2013, ApJ, 763, 9

Zehavi, I., Zheng, Z., Weinberg, D. H., et al. 2005, ApJ, 630, 1

Zitrin, A., Broadhurst, T., Umetsu, K., et al. 2009, MNRAS, 396, 1985

Zitrin, A., Meneghetti, M., Umetsu, K., et al. 2013, ApJL, 762, L30