# Flexion in Abell 2744

Flexion in Abell 2744 Abstract We present the first flexion-focused gravitational lensing analysis of the Hubble Frontier Field observations of Abell 2744 (z = 0.308). We apply a modified Analytic Image Model technique to measure source galaxy flexion and shear values at a final number density of 82 arcmin−2. By using flexion data alone, we are able to identify the primary mass structure aligned along the heart of the cluster in addition to two major substructure peaks, including an NE component that corresponds to previous lensing work and a new peak detection offset 1.43 arcmin from the cluster core towards the east. We generate two types of non-parametric reconstructions: flexion aperture mass maps, which identify central core, E, and NE substructure peaks with mass signal-to-noise contours peaking at 3.5σ, 2.7σ, and 2.3σ, respectively; and convergence maps derived directly from the smoothed flexion field. For the primary peak, we find a mass of (1.62 ±  0.12) × 1014 h−1 M⊙ within a 33 arcsec (105 h−1 kpc) aperture, a mass of (2.92 ±  0.26) × 1013 h−1 M⊙ within a 16 arcsec (50 h−1 kpc) aperture for the north-eastern substructure, and (8.81 ± 0.52) × 1013 h−1 M⊙ within a 25 arcsec (80 h−1 kpc) aperture for the novel eastern substructure. gravitational lensing: weak, galaxies: clusters: general, galaxies: clusters: individual: A2744, dark matter, cosmology: observations 1 INTRODUCTION Over the last two decades, gravitational lensing has served as one of the most successful tools in the development of our standard model of cosmology and extragalactic astronomy (e.g. Blandford & Narayan 1992; Kaiser & Squires 1993; Kaiser, Squires & Broadhurst 1995; Klypin et al. 1999; Van Waerbeke et al. 2000; Wittman et al. 2000; Bartelmann & Schneider 2001; Kravtsov et al. 2004a; Lewis & Challinor 2006; Mandelbaum et al. 2006; Munshi et al. 2008; Tinker et al. 2008; Massey, Kitching & Richard 2010; Kneib & Natarajan 2011; Hoekstra et al. 2013). Galaxy clusters are ideal candidates for observing gravitational lensing effects in background sources, as their high-mass density and broad angular extent ensure a large population of lensed images. While most weak lensing (WL) reconstructions focus on questions of total mass, radial profiles, and the measurement of cluster ellipticities, it is well-known from both numerical simulations (e.g. Klypin et al. 1999; Diemand et al. 2008) and semi-analytical work (Press & Schechter 1974; Sheth & Tormen 1999) that significant substructure is expected in dark matter (DM) haloes at approximately the 2 per cent level. Consequently, the advent of larger and deeper data sets now allows precise and accurate mass reconstructions in both the inner core and outer regions of galaxy clusters through strong lensing (SL) and combined lensing approaches (e.g. Bradac et al. 2005; Leonard et al. 2007; Jauzac et al. 2015a). With improved data sets, new methods, and the knowledge and results gleaned through a multitude of shear testing programmes and challenges (e.g. Heymans et al. 2006; Massey et al. 2007a; Bridle et al. 2010; Kitching et al. 2012; Mandelbaum et al. 2015), the last few years have seen an explosion of new strong (+weak) gravitational lensing analyses, particularly through the data obtained in the recent Hubble Frontier Fields (HFF) initiative (e.g. Medezinski et al. 2013, 2016; Zitrin et al. 2013, 2015; Jauzac et al. 2014, 2015a,b, 2016a,b; Johnson et al. 2014; Lam et al. 2014; Richard et al. 2014; Sendra et al. 2014; Diego et al. 2015a,b, 2016a,b; Grillo et al. 2015, 2016; Sharon & Johnson 2015; Wang et al. 2015; Caminha et al. 2016a,b, 2017; Mo et al. 2016; Mohammed et al. 2016; Sebesta et al. 2016; Treu et al. 2016; Lotz et al. 2017). Although highly effective where possible, SL is limited to regions where the lens projected mass approaches or exceeds a critical density, which extends only a few tens of arcseconds from the cluster centre. Conversely, WL techniques can generally be applied outside of this region, but must smooth over a large-enough region to encompass the number of galaxies required for reasonable signal, making detailed shape information of detectable substructure difficult to obtain. In general, it is very difficult to measure precise shape information for galaxy cluster substructure or even be able to detect substructure in the (1012–1013) M⊙ regime at all, particularly so when sufficiently removed from the centre and/or in complex or merging systems. Simulations of ΛCDM cluster formation suggest that substructure in even well-behaved clusters can be significantly extended away from the central core. Numerous N-body results indicate that the radial distribution of DM subhaloes in galaxy clusters is more extended than that of luminous galaxies (Okamoto & Nagashima 2001; Springel et al. 2001; Gao et al. 2004; Kravtsov, Gnedin & Klypin 2004b), while Nagai & Kravtsov (2005) in agreement with previous studies find that the radial distribution of subhaloes is significantly less concentrated than the overall cluster DM distribution. The subhaloes near the cluster centre lose more than 70 per cent of their initial mass due to tidal stripping, while haloes at much farther distances (approaching the virial radius) lose only ∼ 30 per cent on average, resulting in a subhalo population that is biased away from the central cluster regions. On scales relevant to lensing, Cain, Bradač & Levinson (2016) indicates that at least half of these DM subhaloes would exist outside of the Einstein radius of the main cluster halo, significantly impacting the probability of a multiply-imaged system occurring. We propose that one of the most-promising untapped ways to uncover the detailed structure of dynamic clusters is to look at the higher-order lensing statistic known as ‘flexion’ (Goldberg & Natarajan 2002; Goldberg & Bacon 2005; Bacon et al. 2006). Flexion induces arcing effects in lensed galaxies – similar to but smaller in scale than the giant arcs of SL – forming ‘arclets’ near local mass overdensities. The advantage to using this higher-order information is that it provides a new avenue to detect local substructure independently of shear or SL measurements. Flexion signal has the potential to be much more sensitive to small-scale perturbations of the convergence field than shear while remaining theoretically viable outside the SL region, with utility directly proportional to galaxy-image-substructure proximity. Flexion as an indicator of local structure will provide significantly more information if there is a high density of sources, and thus the HFFs are in a unique position to calibrate flexion measurement and reconstruction techniques for other, wider surveys and applications. In this work, we present a gravitational lensing flexion analysis of the widely-studied galaxy cluster Abell 2744. Our paper is structured accordingly. Section 2 provides a background of the lensing formalism used in this work along with information concerning the data set, observational pipeline, and data reduction process. In Section 3, we summarize our analysis methodology, including details and implementation of flexion measurement and our mass estimators. Section 4 presents the results of our analysis on Abell 2744 and its associated parallel field. Our derived convergence maps, best-fitting parametrized lens models, and flexion aperture mass signal maps are presented with comparisons of results from other lensing groups, and we conclude with an extended summary and discussion. We adopt the standard conventions of ΛCDM cosmology with Ωm = 0.3, ΩΛ = 0.7, Hubble constant H0 = 100 h km s−1 Mpc−1 unless otherwise specified, and with all magnitudes listed in the AB system. Likewise, we adopt complex notation for denoting directional vectors and define the gradient operator as ∂ = ∂1 + i∂2. 2 BACKGROUND 2.1 Lensing and flexion formalism Gravitational flexion was originally developed as a way to quantify gradients in the dimensionless surface mass density field, κ, through highly distorted, gravitationally lensed ‘arclet’ images (Goldberg & Natarajan 2002; Goldberg & Bacon 2005; Bacon et al. 2006; Schneider & Er 2008). As the flexion field is the third derivative of the underlying gravitational potential, it is far more sensitive than shear to substructure in cluster haloes while not requiring the multiple images necessary for a full SL analysis. Accurately, sampling the flexion field at source positions can give critical information about the local convergence field and ultimately any nearby (potentially non-luminous) structural peaks. To reconstruct the lens-plane convergence, κ, with measurable galaxy-image properties, we can establish a foreground lensing potential ψ through the Poisson equation:   $$\mathrm{\partial} ^*\mathrm{\partial} \psi = 2\kappa$$ (1)with derivatives in angular units in the lens (galaxy cluster) plane. The physical mass of the lens system may be determined by relating the line-of-sight projected surface mass density Σ to the dimensionless convergence field:   $$\Sigma (\boldsymbol {\theta }) = \kappa (\boldsymbol {\theta }) \cdot \Sigma _{{\rm crit}}\,,$$ (2)and the critical mass density defined in the usual way. The deflection potential may now be related to a number of observables, including the complex shear field,   $$\gamma = \gamma _1 + i\gamma _2 = \frac{1}{2}\mathrm{\partial} \mathrm{\partial} \psi ,$$ (3)which induces tangential elongation in a source image. Spin-1 flexion, $${\mathcal {F}}$$, has a natural interpretation as the gradient of the convergence, inducing distortion that ‘skews’ source images (coma), while spin-3 flexion, $${\mathcal {G}}$$, has m = 3 rotational symmetry and acts as the gradient of the shear, producing ‘arcing’ effects (trefoil) in an image (see fig. 1 of Bacon et al. 2006):   \begin{eqnarray} {\mathcal {F}}= |{\mathcal {F}}|{\rm e}^{i\phi } = \frac{1}{2}\mathrm{\partial} \mathrm{\partial} ^*\mathrm{\partial} \psi = \mathrm{\partial} \kappa , \end{eqnarray} (4)  \begin{eqnarray} {\mathcal {G}}= |{\mathcal {G}}|{\rm e}^{i3\phi } = \frac{1}{2}\mathrm{\partial} \mathrm{\partial} \mathrm{\partial} \psi = \mathrm{\partial} \gamma . \end{eqnarray} (5)In the very WL regime (κ ≪ 1 ,  |γ| ≪ 1), a lensing coordinate transformation may be sufficiently approximated through linear combinations of these distortion terms. The regions where flexion signal is expected to provide the most benefit however typically have fairly large convergence values that are a significant fraction of unity, and therefore we must consider the effects that magnification and the mass-sheet degeneracy impose on our lensing measurements. As first recognized in Falco, Gorenstein & Shapiro (1985) and expanded by Gorenstein, Shapiro & Falco (1988) and Schneider & Seitz (1995), every lensing observable defined through differential shape distortions will remain unchanged under application of a fixed mass sheet: $$\kappa _{\lambda }^{\prime } = \lambda \kappa + (1-\lambda ),$$ where λ can take any real non-trivial value. The same is true with shear and Flexion:   \begin{eqnarray} \gamma ^{\prime } = \lambda \gamma \,, \quad {\mathcal {F}}{}^{\prime } = \lambda {\mathcal {F}}{}\,, \quad {\mathcal {G}}{}^{\prime } = \lambda {\mathcal {G}}{} \,. \end{eqnarray} (6)To avoid the difficulties associated with degenerate measurements, it is possible to reformulate the lensing fields into new quantities that are invariant to the mass-sheet degeneracy:   $$g = \frac{\gamma }{1-\kappa }\,, \quad \Psi _1 = \frac{{\mathcal {F}}{}}{1-\kappa }\,, \quad \Psi _3 = \frac{{\mathcal {G}}{}}{1-\kappa }\,.$$ (7)where g is the reduced shear, Ψ1 is the reduced spin-1 flexion, and Ψ3 is the reduced spin-3 flexion (Schneider & Er 2008). In this work, references to shear and flexion will refer to their reduced versions unless otherwise noted. 2.2 Observations and data With ∼140 Hubble Space Telescope(HST) orbits for each of six massive cluster lenses and their respective parallel fields, the goal of the Hubble Frontier Field programme is to probe the early universe at redshifts up to and including the z = (6–8) regime, with a limiting magnitude of magAB ∼ 28.5–29 in seven passbands for a 5σ point source detection in a 0.4 arcsec aperture. The first cluster to be fully imaged was Abell 2744 at a redshift of z = 0.308 with a field size of 3.5 arcmin × 3.5 arcmin. At this redshift, assuming the previously stated cosmology gives a scale of 3.175 h−1 kpc arcsec−1 (or 4.536 kpc arcsec−1 if assuming h = 0.7). The angular scale of the data (0.05 arcsec scale, 0.03 arcsec after drizzling), the abundance of colour information, the considerable mass (e.g. Jauzac et al. (2015b) measure 2.76 × 1014 M⊙ within 250 kpc of the core of Abell 2744), and line-of-sight orientations of the cluster lenses create favourable conditions for strong flexion signal, and it is to this first cluster that we apply our flexion pipeline. We take advantage of the self-calibrated v1.0 high-level science images1 released to the Mikulski Archive for Space Telescope (MAST) archives from the HST Frontier Fields Science Products Team. These mosaics are designed as a final product to be used directly for scientific analysis, having been distortion-corrected, combined, and registered to a common pixel scale for the highest possible quality. We use the ACS F435w, F606w, and F814w filters for colour magnitudes, and make shape measurements in the F814w filter, with total respective filter integration times corresponding to 24, 14, and 46 orbits, respectively. 3 METHODS 3.1 Data reduction Flexion analyses pose unique challenges as compared to WL shear. If we consider the simplified case of a singular isothermal sphere (SIS) lens profile, the magnitudes of the shear and flexion field strengths fall off as 1/θ and 1/θ2, respectively. The much sharper distance scale on which flexion image-distortions begin to become prominent results in a scenario where only a small handful of proximal source images are able to strongly contribute to the detection of each local substructure. The majority of galaxy source images that are not readily adjacent to apparent substructure will generally exhibit only a weak or noisy flexion signal. These weaker and more radially distant sources are able to contribute either statistically or by constraining or excluding the possibility of local substructure through their absence of significant flexion distortion, in a similar manner as used for predicting multiply-lensed image positions in SL analysis. This dual-regime-spanning behaviour of flexion signal is likely to lead to a much more individualistic treatment of flexion sources when compared to shear, as any one galaxy image may be a significant part of a substructure detection. While shear reconstructions may require averaging over hundreds of galaxies to produce a significant signal, a flexion reconstruction might only depend on a handful of significant sources near each local structure, and thus would require a much higher signal-to-noise (S/N) ratio per galaxy. To optimize the quality of our galaxy sources and the accuracy of their associated uncertainty maps in preparation for flexion measurement, further handling and manipulation of the Space Telescope Science Institute (STScI)-calibrated HST images is necessary. The process of identifying source galaxies down to the background level, deblending a dense cluster field, extracting individual galaxies while masking other contaminating sources, and targeting the correct sources for measurement is non-trivial. We create a post-processing pipeline to optimize the specifics of our data creation, combining it with the Source Extractor (Bertin & Arnouts 1996) utility which can be used for identifying background signatures that match specified criteria, namely involving pixel thresholding levels, aperture size, deblending, and noise parameters. Alongside the identification of galaxy objects in the field, there are selection cuts that can be made both before and after measurements are made to decrease possible contaminants in the flexion signal, and which depend on the inherent and apparent properties of each source. For example, flexion signal depends on the gradient of the lens convergence and is therefore very sensitive to the apparent size of the lensed source, causing larger objects to be preferentially selected. Additionally, the intrinsic shape of a galaxy can play a significant role in a measurement of its ‘intrinsic flexion’ – the flexion that would be measured if the object was not gravitationally lensed at all. Irregular galaxies in particular can mimic a strong external flexion signal. As the flexion signal can depend on far fewer objects than shear, data reduction and source selection is an extremely important component of this analysis. Our data reduction strategy involves the following general steps: We crop all images to exclude regions that have poor stacking coverage, decreasing background estimation bias in Source Extractor and increasing performance. Using redshift and magnitude catalogues (Owers et al. 2011; Merten et al. 2011), we position-match and identify any established objects belonging to either the cluster itself or the foreground. We identify the foreground-associated pixels down to background level and replace them with a randomized noise map set to the mean of the local background. We implement a ‘hot-cold’ strategy (Rix et al. 2004; Leonard et al. 2007) running Source Extractor in dual-image mode across all appropriate bands – first clearing the field of large bright objects not previously removed, then targeting the expected source population with a minimum of 15 pixels over a 2σ detection threshold. As spectroscopic and photometric redshifts are available for only a fraction of our detected objects, we also use a series of criteria to exclude probable source contaminants or low S/N objects. Namely, we reject bright or large galaxies with mag F814w < 24 or full width at half-maximum (FWHM) > 0.9 arcsec as probable member or foreground contaminants, establish a low signal bound by excluding galaxies with mag F814w > 28, FWHM < (2 × FWHM PSF), or S/N ratio <20. Additionally, we remove objects which are flagged as incomplete (FLAG ≥8) in the Source Extractor catalogues, indicating closely associated or fragmented objects, and mark objects which have FLAG totals <8 to use as a discriminant in later visual inspection (see Section 3.2.7). We generate square postage stamps of selected sources. These are centred at the galaxy centroids and have a windowed radial extent set to 4.5 times the calculated half-light radii, chosen to be large enough to ensure that the flexion-susceptible galaxy wings and background zero-constraints are included. Alongside the image stamp, we generate an associated 1σ noise stamp that includes both background sky noise and Poisson noise. 3.2 Flexion measurement – Analytic Image Modeling While research into the potential power of flexion measurement continues to remain popular (Er, Tereno & Mao 2012; Viola, Melchior & Bartelmann 2012; Cain, Bradač & Levinson 2016; Cardone et al. 2016; Rexroth, Natarajan & Kneib 2016), there have only been a handful of flexion analyses applied to real data to date (Leonard et al. 2007; Okura, Umetsu & Futamase 2008; Cain, Schechter & Bautz 2011; Leonard, King & Goldberg 2011), with several using the same data set and the majority in the widely studied rich cluster Abell 1689. Flexion measurements of real, individual sources have been achieved through a few techniques, including by measuring a combination of various third (and higher) order moments of the light distribution (as originally suggested by Goldberg & Natarajan 2002; Okura et al. 2008), decomposing the projected galaxy shape on to a polar orthonormal ‘shapelets’ basis set and truncating the series at a particular threshold (Goldberg & Leonard 2007; Massey et al. 2007b), and exploring the local potential field through parametrized ray-tracing (Cain et al. 2011), known as Analytic Image Modeling (AIM). We utilize a modified version of the AIM technique to determine the localized flexion field values for individual sources. AIM is distinct in that instead of measuring derived quantities (such as weighted surface brightness moments), it fits the lensed galaxy objects using a parametric model. By comparing the observed image to the uncertainty-weighted model image, we are able to optimize the parameters over reasonable bounds and thus constrain the flexion fields (and other shape information.) 3.2.1 Model parametrization Galaxy profiles are largely fit to parametric models with a radial luminosity distribution, and for this work we implement an elliptical Sérsic intensity profile. The Sérsic profile is particularly useful in that it encompasses a range of different models, including ones already well established through galaxy-luminosity profiling, including exponential, Gaussian, and de Vaucoulers profiles. Our model takes the form   $$I(\theta ) = I_e \exp {\left\lbrace -b_n\left[\left(\frac{\theta }{\theta _s}\right)^{1/n}-1\right]\right\rbrace },$$ (8)where θ is the circularized radius, θs is the radius of the isophote containing half the total flux of the galaxy (the half-light radius), Ie is the brightness at this effective radius, and n is the Sérsic index, which controls the steepness of the profile. We also include a complex ellipticity such that   $$|e| \,=\, \frac{q^2 - 1}{q^2 + 1}\,,$$ (9)and where q represents the axis ratio. Finally, bn is empirically determined to set θe to the half-light isophotal radius, and we adopt an approximated functional form   $$b_n = 1.992n - 0.3271, \quad 0.5 < n < 8.0$$ (10)valid within the given Sérsic index range (Capaccioli 1989). 3.2.2 Flexion, shear, and intrinsic ellipticity coupling As shown in Viola et al. (2012), measurements of the spin-1 and spin-3 distortions in a lensed image cannot be straightforwardly related to the $${\mathcal {F}}{}$$ and $${\mathcal {G}}{}$$ flexions, in general. In the intermediate regime between very weak and SL, coupling terms between the shear and flexion produce non-negligible spin-1 and spin-3 terms that contribute to the distortions originating from flexion alone. A significant inherent ellipticity in a source will also produce a similar effect through flexion–ellipticity coupling. In particular, moment-based approaches can introduce large biases in flexion measurements when failing to take these terms in into account, especially when accompanied by a centroid shift that depends on the ellipticity and galaxy morphology itself. Forward model-fitting techniques such as AIM avoid some of these issues by including an ellipticity/shear optimization term that accounts for some of the parameter coupling, but inherent source ellipticity can still introduce bias during parameter optimization – in general there is no appropriate application of a spin-1, spin-2, and spin-3 distortion to a circular source that is able to reproduce the same shape as when these distortions are applied to a source with inherent ellipticity. Attempting to fit both inherent ellipticity and shear simultaneously to an observed galaxy will likely result in non-convergence, as first identified in Cain et al. (2011) using simulated data. We approach this issue by fixing the lensing shear γ to zero and allowing the intrinsic ellipticity to absorb the two degenerate parameters into one. Although it is possible to set the lensing shear to that of a predetermined model, we instead aim to obtain an unbiased estimate of the shear field alongside the targeted flexion fields. As noted, it is possible that bias may be introduced by neglecting the flexion-coupling terms from inherent source ellipticity. We investigate this issue by rerunning our analysis with the source lensing shear values fixed to model data. We utilize an independent shear map drawn from public HFF SL cluster reconstructions (Lotz et al. 2017) and interpolate shear priors for our flexion sources before fitting with AIM. The resulting convergence and aperture mass reconstructions are very similar to those produced with null prior flexion, with peak locations and enclosed mass deviations within measurement uncertainties, |δθpeak| < 2 arcsec, |δMenclosed| < 1013 M⊙. The signs of these deviations are not consistent across all peak measurements and do not clearly show any evidence of systematic trends towards a particular direction. From this, we confirm the trend found in Cain et al. (2011) that shear priors on average do not drastically influence measured flexion values when using the AIM measurement technique. Removing this parameter degeneracy, there are seven model parameters for the unlensed galaxy light profile, and four effective parameters for the lensing fields. While we can reasonably expect this model to work well for simple source galaxies, it is also important to explore typical parameter spaces for existing data and to consider any known limitations. 3.2.3 Typical model suitability One of the biggest problems for flexion measurement, in general, is the inherent inhomogeneity of galaxy structure, particularly when combined with the emphasis on a select few sources (again as a result of strong local mass sensitivity). A complete analysis of ‘inherent flexion’ in unlensed populations has yet to be fully investigated at this point, though previous analyses have done so in cluster lensed populations (Leonard et al. 2007; Okura et al. 2008; Cain et al. 2011). Most galaxies are not simply circular or elliptical, and many can be expected to contain structural irregularities such as arms, bars, or star-forming regions. Elliptical models tend to have a large range of Sérsic indices, whose values depend strongly on luminosity (Blanton & Moustakas 2009). Furthermore, ellipticals can also be separated into ‘boxy’ and ‘discy’ subtypes. For the most part, however, a Sérsic model can be used to model the overall structure of most galaxy types. For elliptical galaxies, the Sérsic index n generally reflects what is apparently a single component galaxy. Spiral galaxies can also be described by a Sérsic index, but in this case n reflects a balance between the disc and the bulge, two clearly distinct components. Morphology typically focuses on the separation of the disc from the bulge, usually treating the disc as an exponential profile and fitting the bulge to a Sérsic profile, with n dependent on the type of the central component. Along with lenticular galaxies, these form the basis of our unlensed galaxy reconstruction efforts. The range of possible parameters create a broad distribution of shapes; drawing a good statistical likelihood of inherent shape is difficult. In this work, we aim to partially investigate inherent flexion by performing a simultaneous flexion analyses on the offset parallel Abell 2744 field. 3.2.4 Distortions and PSF corrections Our images consist of many offset, rotated, and stacked individual exposures over noncontiguous time periods (with a total average integration time of over 104 000 s, or ∼29 h in the F814w band) increasing the object S/N by a factor roughly equivalent to the number of orbits, at the cost of possible geometric/elliptical distortions and a more complicated effective Point Spread Function (PSF) that is non-trivial to model. The 90-min HST orbital period introduces thermal breathing fluctuations, while the telescope focus deviates over periods of weeks, and our data spans a range of these cycles. Furthermore, the off-axis position of the ACS introduces a spatially varying geometric distortion across the WFC chips. Despite the complexities inherent in the exact shape of the ACS PSF, Schrabback et al. (2010) finds that ∼ 97 per cent of the total PSF ellipticity variation in random pointings can be attributed to a single modelled parameter – the telescope focus. While there are techniques that have been successfully applied to model the PSF of stacked ACS images (e.g. Bacon et al. 2003; Rhodes et al. 2007; Harvey et al. 2015; Jauzac et al. 2016b), these have mostly found use in shear measurement, with accuracy requirements at the ∼1 per cent level. As the induced flexion in a source from the PSF is   $${\mathcal {F}}_{\mathrm{induced}} \sim {\mathcal {F}}_{\mathrm{PSF}} \frac{a_{\mathrm{PSF}}^4}{a_{\mathrm{source}}^4 + a_{\mathrm{PSF}}^4}$$ (11)as derived in Leonard et al. (2007), for typical WFC values ($${\mathcal {F}}_{{\rm PSF}} \sim 10^{-3} \, {\rm arcsec^{-1}}, a_{{\rm PSF}} \sim 0.1125{\rm \,arcse}{\rm c}$$) there is not a significant flexion contribution from the PSF, provided the source is sufficiently large. As one of our pre-measurement cuts requires Source Extractor sources to have an FWHM greater than twice that of the ACS PSF (2 × 0.1125 arcsec), this minimizes the amount of any possible PSF-induced flexion. To test the flexion anisotropy across the ACS WFC chips, we used the tinytim2 software to simulate an ACS PSF under typical instrument parameters. By varying the input pixel coordinates across the two chips, we created a grid of spatially varying PSFs representing a single exposure, and were then able to measure the flexion signal using our AIM implementation. The resultant $${\mathcal {F}}$$-flexion vectors vary smoothly in magnitude and direction across the chips with a maximum magnitude of $$|{\mathcal {F}}| = 0.005$$, much smaller than the lensing flexion signal expected in Abell 2744, particularly in sources larger than the PSF. In addition, stacking the offset and rotated individual exposures into our integrated data image reprojects and averages the directional biases of the flexion anisotropies, further reducing possible measurement biases. In practice, as an approximation we use a simple Gaussian convolution with FWHM $$= 2\sqrt{2\mathrm{ln}\,2}\sigma = 0.1125{\rm \,arcsec}$$ to account for the general smoothing of a PSF. We expect the PSF to have relatively small effects on the overall shape and shape parameters of our selected galaxies, but as a self-consistency check we perform the best fit again without any applied PSF. If the fit shape or model parameters are markedly different we exclude that object from our analyses. 3.2.5 Implementation To determine the best-fitting model parametrization, we use a Levenberg–Marquardt least-squares minimization scheme, the core of which is implemented through a modified version of the python translation3 of the idl code mpfit. Our pipeline combines this minimization technique with galaxy models simulated through the GalFlex4 module, which ray-trace second-order weak gravitational lensing effects (flexion), along with shear and convergence, on simulated source images in real space. While a full exploration of the 11-dimensional parameter space is possible through combined downhill and sampling techniques (e.g. simulated annealing), for this work we take advantage of initial galaxy parameter estimates to aid global minimization. 3.2.6 Goodness of fit and parameter uncertainties Already noted in Section 3.1, the source-selection process is a critical part of flexion analysis. In a lensing reconstruction of Abell 1689 using ACS data, Leonard et al. (2007) combines SL, shear, and flexion analysis to confirm substructure in the form of a second mass peak apart from the central cluster. However, the central mass peak is not recovered in their flexion analysis despite having a large mean number density of flexion sources ($$\bar{n}_g \sim 75$$ arcmin−2). A later flexion-only analysis of A1689 by Cain et al. (2011) is able to infer mass structure consistent with previous measurements with a mean flexion source density of $$\bar{n}_g \sim 26$$ arcmin−2. Okura et al. (2008) finds that they can recover significant structure in the core of A1689 with only nine galaxies using Subaru data, through a meticulous source-selection procedure resulting in a final source density of $$\bar{n}_g \sim 7.75$$ arcmin−2. The lack of a significantly stronger signal in the Leonard et al. (2007) analysis despite the larger source density indicates the necessity of selecting which galaxy sources effectively ‘make the cut.’ Our final flexion source catalogues is constructed by excluding the galaxies that have either poor or unphysical best-fitting parameters. A commonly used metric for determining a ‘good’ fit is the reduced chi-squared statistic $$\chi ^2_{\mathrm{red}}$$, defined as the weighted sum of squared deviations (the χ2) divided by the number of independent degrees of freedom of the system. Although we use χ2 as the optimized statistic in AIM, we do not use the $$\chi ^2_{\mathrm{red}}$$ as the main selector for cutting our post-fit sources as in Cain et al. (2011) – highly correlated pixel data combined with a highly non-linear model render the actual meaning of the $$\chi ^2_{\mathrm{red}}$$ uncertain. We point the reader to Andrae, Schulze-Hartung & Melchior (2010) for further details. Nevertheless, this statistic is useful to broadly separate groups, as large values do indicate a likely poor fit. Fig. 1 shows an example of accepted, uncertain, and rejected galaxy source fits and their residuals, separated by $$\chi ^2_{\mathrm{red}}$$ for convenience. Figure 1. View largeDownload slide Best-fitting Sérsic galaxy models in the Abell 2744 HST F814w band. Each row shows a typical reduced χ2 value used in part to group objects by a goodness-of-fit statistic. As systematic uncertainties in inherent model suitability, best-fitting model error analyses, and model-parameter-error-size correlation are not yet well-understood, this statistic is not the defining accepted/rejected source criterion, but rather a weighted component of an extended source-selection pipeline. Figure 1. View largeDownload slide Best-fitting Sérsic galaxy models in the Abell 2744 HST F814w band. Each row shows a typical reduced χ2 value used in part to group objects by a goodness-of-fit statistic. As systematic uncertainties in inherent model suitability, best-fitting model error analyses, and model-parameter-error-size correlation are not yet well-understood, this statistic is not the defining accepted/rejected source criterion, but rather a weighted component of an extended source-selection pipeline. Instead of the reduced chi-squared statistic, we employ an impartial stamp inspection Graphical Universe Interface (GUI) as described in Section 3.2.7, in combination with the following parameter cuts: We reject fits where the reduced flexions {|Ψ1|, |Ψ3|} ≥ 1.0 arcsec−1. These large values indicate either a poor or unphysical fit, or that the source is in a region where WL assumptions are not valid. We also impose the constraint that $$\lbrace \sigma _{\Psi _1}, \sigma _{\Psi _3}\rbrace > 0.001$$ arcsec−1 to prevent overfitting the associated parameters {Ψ1, Ψ3}, as in Cain et al. (2011). We reject object fits with axis ratio q = 1.0 in the same vein, as a best fit which is exactly symmetric is likely to either be a star or overfitted. Finally, for practical convenience, we also eliminate objects whose best-fit χred2 is large (≳5), in most cases representing either an improperly deblended or foreground galaxy. 3.2.7 Interactive-limited stamp inspection While the Source Extractor utility turns the non-trivial procedure of source identification, deblending, and shape estimation into something more routine for most purposes, the various strategy and parameter combinations – particularly for crowded fields with diverse redshift populations – can still prove a need for artful navigation. Parameters must be selected so as to identify both the small and faint sources, as well as the larger, brighter, and closer ones; as described in Section 3.1, a common technique in WL studies is to perform a ‘hot-cold’ routine. However, deblending issues still can persist throughout the pipeline. These issues can entail multiple galaxies masquerading as one, or vice versa, especially with extended galaxies containing luminous star-forming regions. Therefore, an additional layer of inspection can benefit a signal dependent on an accurate measure of a single galaxy's shape. While performing this kind of source inspection in a shear analysis would be time-prohibitive, this is not necessarily true for flexion analysis. As shown in Viola et al. (2012), accurate flexion measurements are particularly dependent on high S/N sources, and so by selecting objects according to suitable criteria the number of usable galaxies for flexion measurement is narrowed down to a much smaller and more feasible subset. Our stamp/fit inspection routine is implemented through a GUI as part of our flexion pipeline. Built with the python bindings for Qt (PyQt4), we compile the base code and all relevant modules into an executable to create a stand-alone application. Seen in Fig. 2, the GUI can be used to quickly and efficiently scroll through input stamp data sets for acceptance/rejection. We show a slightly extended region around the selected stamp in order to inspect its immediate neighbours for blending issues or galaxy fragmentation. We purposefully restrict the field of view in an attempt to minimize any observational bias. We also inspect the best-fitting parameters to check for problematic behaviours, including having values on the established parameter bounds or exactly null, or having particular combinations of values that clearly indicate a poor fit (e.g. a best-fit half-light-radius larger than the stamp itself). Figure 2. View largeDownload slide Inspecting our image and best-fitting stamp and fit parameters is critical for an analyses that can depend on just a few handful of strong signal images. We use a custom GUI to streamline the process of accepting/rejecting/flagging these sources for use in the mass reconstruction. Figure 2. View largeDownload slide Inspecting our image and best-fitting stamp and fit parameters is critical for an analyses that can depend on just a few handful of strong signal images. We use a custom GUI to streamline the process of accepting/rejecting/flagging these sources for use in the mass reconstruction. 3.3 Mass estimators A major component of flexion analysis is the non-trivial process of incorporating flexion measurements into an overall mass reconstruction. Ideally, both shear and flexion measurements could be used to constrain the mass profile more than either could alone. As we've noted, the shear probes large scales more effectively, while flexion is particularly suited to investigate smaller-scale structure. However, the flexion field alone can give valuable information about both small-scale structure in a lensing cluster and the overall mass profile. To do so, we employ two independent techniques to get a measure of a mass structure based solely on positional flexion measurements: a direct calculation of the convergence field, stemming directly from the fact that flexion is the gradient of the convergence; and the construction of a mass-signal map through the aperture mass technique (Schneider 1996), which relates the weighted integral of flexion measurements within an aperture to the respective weighted integral of the convergence within that aperture (Cain et al. 2011; Leonard et al. 2011). Each technique offers a different interpretation of the flexion measurements, and the utilization of both can provide a broader image of any substructure present within the galaxy cluster. 3.3.1 Direct flexion-convergence convolution Using the various directional derivative relations between lensing distortions and a gravitational potential, we are able to relate (up to the mass-sheet degeneracy) the convergence field κ with the flexion fields $${\mathcal {F}}, {\mathcal {G}}$$ through the convolution:   $$\kappa (\boldsymbol {\theta }) = \int _{\mathbb {R}^2} d^2 \theta ^{\prime } \, \mathcal {D}_M (\boldsymbol {\theta } - \boldsymbol {\theta }^{\prime }) \, {\mathcal {F}}(\boldsymbol {\theta }^{\prime }),$$ (12)where $$\mathcal {M}$$ is the $${\mathcal {F}}$$ or $${\mathcal {G}}$$ flexion field in the lens plane, and the respective convolution kernels are given by   \begin{eqnarray} \mathcal {D}_{{\mathcal {F}}} = \frac{\theta _1 + i\theta _2}{2\pi |\theta |^2}, \end{eqnarray} (13)  \begin{eqnarray} \mathcal {D}_{{\mathcal {G}}} = \frac{(\theta _1 + i\theta _2)^3}{2\pi |\theta |^4}. \end{eqnarray} (14)as outlined in Leonard et al. (2011). To avoid the effects of infinite sampling noise (Bartelmann & Schneider 2001), the irregularly spaced flexion data is first smoothed with an uncertainty-weighted Gaussian kernel before convolving the resulting finely spaced flexion grid with the appropriate flexion convolution kernel. The resulting convergence map is a soft function of the selected smoothing scale, which generally reflects the scale of structure that is emphasized. To obtain estimates of the uncertainty in our convergence reconstructions, we investigate the corresponding uncertainties in both the values of the measured flexions (or shear) and the shot noise as a result of the scattered spatial positions of selected flexion (shear) sources. We simulate the prior probability distributions for these observables and ultimately apply Monte Carlo resampling methods to determine the variance in the derived convergence map as well as observational derivatives such as our estimates of enclosed peak masses. In order to sample the probability distribution, we make an assumption that our flexion measurement uncertainties {$$\sigma _{{\mathcal {F}}{}_i}, \sigma _{{\mathcal {G}}{}_i}$$} are approximately normally distributed and follow the relation:   \begin{eqnarray} \sigma _{{\mathcal {F}}{}_i ({\mathcal {G}}{}_i)} = \sqrt{ \left(\frac{\sigma _{\mathcal {F}\cdot a}}{a}\right)^2 + (\sigma _{\mathcal {F}_i, {\rm AIM}})^2}, \end{eqnarray} (15)where the term $$\sigma _{\mathcal {F}\cdot a}$$ (or $$\sigma _{\mathcal {G}\cdot a}$$) indicates the 1σ variation in the measured values of intrinsic flexion inherent in the unlensed galaxy population, added in quadrature to $$\sigma _{\mathcal {F}_i, {\rm AIM}}$$ (or $$\sigma _{\mathcal {G}_i, {\rm AIM}}$$), the individual flexion measurement uncertainty derived from the best-fitting model's covariance matrix. For clarity, we note that the population dimensionless flexion dispersions $$\sigma _{\mathcal {F}\cdot a}$$ and $$\sigma _{\mathcal {G}\cdot a}$$ are calculated from the combined sets of real and imaginary components {$$\mathcal {F}_1 \cdot a\,, \, \mathcal {F}_2 \cdot a$$} , {$$\mathcal {G}_1 \cdot a\,, \mathcal {G}_2 \cdot a$$} of spin-1 and spin-3 flexions, respectively, which are assumed to be normally distributed with intrinsic mean value of zero. For each iteration in a large set of successive convergence reconstructions, the resampled flexion values are drawn from a normal distribution with mean and standard deviation equal to the measured flexion values and corresponding 1σ uncertainties (including both model-based and population uncertainty). By stacking the value-resampled reconstructions, we are able to construct a two-dimensional variance map along the stacked axis. Our peak mass uncertainties are obtained from the variance of the enclosed aperture mass as measured in each of these frames. While this convolution relation offers a direct linear inversion map with a straightforward interpretation, it suffers from the same difficulties as the shear convolution, namely that: The integral in equation (12) extends over $$\mathbb {R}^2$$, while usable flexion data are constrained to a few square arcminutes. Smoothing over the lensing fields is necessary to avoid infinite noise from discrete sampling, which is a particular concern as the strength of the flexion signal falls off much more quickly than shear further from mass structure. Source-position shot noise can significantly affect local field values depending on smoothing scale. A low source density and large variance of nearest-neighbour separation distance can strongly influence the relative importance of a single source in any particular smoothing aperture. Using multiple scales can mitigate potential local field extremes by allowing narrow-band peak outliers to be identified and removed. The cluster regions where flexion signal is expected to be strong are usually those where the convergence cannot be assumed to be small, and thus the observables are better approximated by the reduced flexion Ψ1 and Ψ3. To offer an alternative and complementary method to the direct convergence inversion, we also utilize our implementation of the aperture mass technique. 3.3.2 Aperture mass statistic The aperture mass statistic, $$\mathcal {F}_{{\rm aper}}$$, is able to relate measured flexion to the underlying lens convergence through more general weighting kernels than the direct flexion-convergence convolution. Apertures of characteristic radius R are laid down in a grid-like pattern over the field. Within each aperture the convergence is convolved with a filter function w(r). This convolution is then related to the measured flexion convolved with an appropriate filter function Q(r). The role of the aperture mass technique is to evaluate mass structure detections. By randomly rotating the flexion vectors in the image and running the reconstruction multiple times, it is possible to create a S/N map of detected mass signal, analogous to the flux S/N ratio of a galaxy, for example. The aperture statistic is defined by relating the spin-1 flexion and convergence as follows:   \begin{eqnarray} \mathcal {F}_{{\rm aper}}{}(\theta _0;R) = \int _{|\theta |\le R} {\rm d}^2\, \theta \, \kappa (\theta +\theta _0)w(|\theta |), \end{eqnarray} (16)  \begin{eqnarray} \mathcal {F}_{{\rm aper}}{}(\theta _0;R) = \int _{|\theta |\le R} {\rm d}^2\, \theta \, \mathcal {F}_E(\theta ;\theta _0)Q(|\theta |), \end{eqnarray} (17)where $$\mathcal {F}_E$$ is the component of the first flexion oriented towards the centre of the aperture. The weight and filter functions are free to be selected, while obeying the relations   $$Q(r) = -\frac{1}{r} \int _0^r w(r^{\prime }) r^{\prime } {\rm d}r^{\prime } \quad w(r) = -\frac{1}{r} Q(r) - \frac{{\rm d}Q}{{\rm d}r},$$ (18)where the variable r represents |θ − θ0|. The weight function is constrained to go to zero smoothly at the aperture boundary and have a mean value of zero,   $$\int _0^R w(r^{\prime }) r^{\prime } {\rm d}r^{\prime } = 0.$$ (19)In this work, we use a family of polynomial weighting functions as described and used in the literature (Schneider et al. 1998; Cain et al. 2011; Leonard et al. 2011). These filter functions are chosen to test the persistence of mass estimates across a variety of scales. The convergence weighting function is defined as   $$w(r) = A_l \frac{(2 + l)^2}{\pi } \left(1 - \frac{r^2}{R^2}\right)^l \left(\frac{1}{2+l} - \frac{r^2}{R^2}\right)$$ (20)along with the flexion kernel   $$Q(r) = -A_l \frac{2 + l}{2\pi } r \left(1 - \frac{r^2}{R^2}\right)^{1+l}$$ (21)and normalization   $$A_l = \frac{4}{\sqrt{\pi }}\frac{\Gamma (\frac{7}{2} + l)}{\Gamma (3 + l)}$$ (22)as a function of the polynomial index, ℓ, and the characteristic aperture radius, R. A larger polynomial index ℓ makes the kernel more sensitive to smaller scale structure within the aperture, at the expense of having larger noise fluctuations. Lower indices smooth over a larger area, reducing noise at the cost of resolution. The aperture radius R has a similar effect: as the radius increases the flexion signal associated with the centre of the aperture falls off quickly, smoothing over smaller scale structures. To implement this method, we turn the flexion convolution into a discrete sum of kernel-weighted flexion sources at each point θ0 over a regular coordinate grid. We include flexion uncertainty in the kernel Q(r) by introducing a normalized weighting term qn applied to each flexion source:   $$Q(r) \rightarrow Q_i(r) \, \,q_i ; \quad q_i = \frac{1}{\sigma _{{\mathcal {F}}{}_i}^2},$$ (23)where the flexion uncertainty, $$\sigma _{{\mathcal {F}}{}_i}$$, refers to the total combined uncertainty from both the AIM measurement and the inherent population dispersion (equation 15). In practical application, the methods used to implement both the flexion-convergence inversion and the aperture mass statistic techniques are quite similar, though the noise and smoothing properties as well as the physical meanings of the reconstructions are different. 3.3.3 Flexion profile fitting While the two techniques as described above are able to give more than sufficient information for high-quality mass reconstructions, they are both subject to similar types of restrictions and assumptions. As a third perspective we directly fit gravitational lens models to our flexion data to further predict and constrain the position, shape, and size of cluster substructure. Parametric lens modelling has the advantage of using observed data more directly than by using a smoothed data set along with simultaneously being capable of breaking the mass-sheet degeneracy by constraining the slope of the shear/flexion data. We fit a SIS profile to each of the identified peak locations and allow the centre coordinates to moderately vary. The best-fitting values of the parameter set {θ1, 0, θ2, 0, θE} are optimized through a comparison between our flexion data and the model flexion,   $${\mathcal {F}}{} = -\frac{\theta _E}{2(\theta -\theta _0)^2} \, e^{i\phi }\qquad {\mathcal {G}}{}(\theta ) = \frac{2\theta _E}{3(\theta -\theta _0)^2} \, e^{3i\phi },$$ (24)where the preceding constant θE is the Einstein radius of the lens, which is related to the lensing depth and velocity dispersion, $$\sigma _{\nu _{{\rm 1D}}}$$, as   $$\theta _E = 4\pi \left(\frac{\sigma _{\nu _{{\rm 1D}}}}{c}\right)^2 \left(\frac{D_{{\rm ls}}}{D_{{\rm s}}}\right).$$ (25) We utilize the non-reduced formulations of the flexion in this model fitting – this avoids the direction-reversing behaviour that the SIS model tends to predict for the measured flexion within appreciable distances of the origin. The source galaxies that would be affected by this behaviour would exhibit radial shear with highly distorted profiles according to this model, and would be selected out of our source sample before flexion measurement. A result of this approximation is that the best-fitting model estimates for the Einstein radii and velocity dispersions may be somewhat larger than if a more accurate flexion profile and mass-degeneracy treatment were utilized. While fitting to more complicated profiles is possible, the SIS lens is commonly referenced in both past and current lensing work as a simple and analytically tractable model with well-understood model parameters. As a complement to our non-parametric reconstruction maps, the SIS flexion model should prove sufficient for the purposes of this work. 4 RESULTS AND SUMMARY 4.1 Source selection results Our extraction and pre-measurement selection routine returned 1344 out of the 4249 detected galaxies in the F814w-observed 3.5 arcmin × 3.5 arcmin field of Abell 2744. The after-measurement selection cuts give us our final catalogue of 969 sources, corresponding to 82.1 sources arcmin−2 or about 23 per cent of the initially detected galaxy population. In the cluster's associated parallel field, our pre-selection routine returned 1001 out of 4142 detected galaxies. Post-measurement cuts reduce our final catalogues to 923 usable galaxies, corresponding to 78.2 sources arcmin−2. As the observational details of the parallel field are identical to those of the cluster (e.g. a field of the exact same size, with the same total integration time in each band, in the same number orbits), a source density comparable to the cluster field is in line with the expectation that our cluster selection routine is able to exclude contaminating cluster members. We measure the spin-1 and spin-3 dimensionless flexion dispersions in our cluster source population to be $$\sigma _{\mathcal {F}\cdot a} = 0.09$$ and $$\sigma _{\mathcal {G}\cdot a} = 0.12$$, with similar quantities calculated for the selected parallel field population. The mean dimensionless flexions $$\langle {\mathcal {F}}{}\cdot a\rangle , \langle {\mathcal {G}}{}\cdot a\rangle$$ were found to be consistent with zero in both the cluster and parallel populations as well. Previous flexion analyses (Leonard et al. 2007; Okura et al. 2008; Cain et al. 2011) have found population dimensionless flexion dispersions of $$\sigma _{\mathcal {F}\cdot a}\simeq (0.03 {\rm -}0.05)$$ and $$\sigma _{\mathcal {G}\cdot a} \simeq (0.04 {\rm -}0.06)$$, roughly a factor of 2–3 times smaller than our calculated values. It is unclear if this population difference is simply due to chance or if possible systematic differences exist between the different flexion procedures. As all other previous flexion analyses have focused on the galaxy cluster Abell 1689, it is possible that differences in observational details (redshifts, integration time, etc.) or properties of the selected source populations may account for some or all of the variance in measured flexion. In addition to our flexion statistics, by deriving shear estimators from measured ellipticity $$\widetilde{\gamma } \approx e$$, we also note the AIM-associated shear population statistics: $$\sigma _{\gamma _1} = 0.30, \sigma _{\gamma _2} = 0.32,$$ 〈γ1〉 = 0.0037, 〈γ2〉 = −0.0361. These shear dispersions are typical of WL studies for the most part, although there is a non-negligible mean shear component γ2. As the high-resolution HST data set used in this work extends only a few arcminutes from the cluster centre, this systematic mean shear bias is much less than the expected lensing shear signal (γ ∼ 0.1) at fields bounds. However, most wide field shear analyses commonly attempt to measure cluster lensing signal at extended distances much farther than those explored here, and a similar mean shear residual would need to be accounted for in order to observe an actual lensing signal in the range of a few percentages. We approximate our effective background source redshifts by averaging the lensing depths determined in two previous lensing studies using the same HFF images. The source sample used in Medezinski et al. (2016) was determined to have 〈zeff〉 = 1.27 ± 0.08 with a mean lensing depth of Dls/Ds = 0.68 ± 0.11 while the sample used in Jauzac et al. (2016b) had an effective source redshift of 〈zeff〉 = 1.586 with a mean lensing depth of Dls/Ds = 0.729. Taking a simple average yields an effective critical density of Σcrit = 3.6 × 1010 h−1M⊙ arcsec−2. 4.2 Mass reconstruction results Following the non-parametric methods as described in Section 3.3, we recreate the mass distribution of Abell 2744 over a 5 arcmin × 5 arcmin square plane on which our 3.5 arcmin × 3.5 arcmin data field is diagonally circumscribed. We use an equirectangular projection of the sky plane as a planar coordinate system, oriented such that declination, δ, increases along the positive y-axis to the north; and right ascension, α, decreases along the positive x-axis towards the west. The origin of our relative sky coordinate system is set to the equatorial coordinates: {α = 3.5890213°, δ = −30.397196°} (deg, J2000.0), as the centre position in our stacked HST pixel data image. The flexion measurements from the lensed images can be visualized in a number of ways. In Fig. 3, we plot a cluster-scale reconstruction of the convergence map of Abell 2744 as well as the parallel field noise reconstruction, using $${\mathcal {F}}$$-flexion data exclusively. At a smoothing scale of 15 arcesc, only two of three substructures are detected – a large central peak that follows the elliptical luminosity structure of the central cluster galaxies, declining smoothly towards the western edges, and a second peak centred just south of the bright spider-like foreground galaxy to the east, offset 1.43 arcmin (273 h−1 kpc) from the cluster core. Figure 3. View largeDownload slide Independent $${\mathcal {F}}$$-flexion data are used to create cluster-scale mass field reconstructions of Abell 2744 and noise in its associated parallel field. The top left figure shows the cluster convergence (up to a mass-sheet degeneracy scaling) using a 15 arcsec smoothing kernal, while the top right shows a parallel field noise reconstruction under the same parameters. The bottom left and right figures are the flexion aperture mass S/N reconstructions for the cluster and parallel fields, respectively, at aperture radius 60 arcsec and index ℓ = 5. The axes are in units of arcseconds. Figure 3. View largeDownload slide Independent $${\mathcal {F}}$$-flexion data are used to create cluster-scale mass field reconstructions of Abell 2744 and noise in its associated parallel field. The top left figure shows the cluster convergence (up to a mass-sheet degeneracy scaling) using a 15 arcsec smoothing kernal, while the top right shows a parallel field noise reconstruction under the same parameters. The bottom left and right figures are the flexion aperture mass S/N reconstructions for the cluster and parallel fields, respectively, at aperture radius 60 arcsec and index ℓ = 5. The axes are in units of arcseconds. While the 15 arcsec smoothing scale presented in Fig. 3 only hints at an offshoot of increased convergence extending from the central region towards the north-east, and a second towards the south-east, finer smoothing scales in addition to our mass aperture reconstructions reveal apparent mass peaks within both tendrils, with a more substantial north-eastern peak and a lesser south-eastern detection. Using a 60 arcsec smoothing kernel, on the other hand, reveals a third structure, towards the north-east. We can get a better sense of the detailed shape of the substructures – and in particular the central substructure – by using a complementary map, the aperture mass reconstruction, in Fig. 4. Finally, in Fig. 5 presents the results of our shear-only convergence reconstructions, with noise maps and two different scaling kernels. We describe more detailed interpretations of these maps below. Figure 4. View largeDownload slide Magnified section of flexion aperture mass reconstructions centred around the detected core substructure. The aperture characteristic truncation scale $$\mathcal {R} = \lbrace 30{\rm },45,60,75\, {\rm arcsec}\rbrace$$ varies along the vertical axis while the flexion weighting function's polynomial index ℓ = {3, 5, 7} spans the horizontal. The plots are arranged such that the effective resolution of the aperture maps is the broadest towards the lower left corner of the figure, and the finest towards the top right corner, with approximate resolution scaling $$\Delta \sim \mathcal {R}_{{\rm aper}} / \ell$$. At higher effective resolutions, the detected core substructure appears to begin to deblend into multiple components – the most significant peak focuses towards the region associated with the southern BCG pair, while mass signature begins to show both slightly north-west of the BCG pair and also around the slightly easterly spiral galaxy. Our peak-moment analysis at $$\mathcal {R} = 60{\rm \,arcsec}, \ell = 5$$ shows a 13 arcsec offset between the maximum value peak and the centroid, indicating either a very skewed or dimodal structure. Figure 4. View largeDownload slide Magnified section of flexion aperture mass reconstructions centred around the detected core substructure. The aperture characteristic truncation scale $$\mathcal {R} = \lbrace 30{\rm },45,60,75\, {\rm arcsec}\rbrace$$ varies along the vertical axis while the flexion weighting function's polynomial index ℓ = {3, 5, 7} spans the horizontal. The plots are arranged such that the effective resolution of the aperture maps is the broadest towards the lower left corner of the figure, and the finest towards the top right corner, with approximate resolution scaling $$\Delta \sim \mathcal {R}_{{\rm aper}} / \ell$$. At higher effective resolutions, the detected core substructure appears to begin to deblend into multiple components – the most significant peak focuses towards the region associated with the southern BCG pair, while mass signature begins to show both slightly north-west of the BCG pair and also around the slightly easterly spiral galaxy. Our peak-moment analysis at $$\mathcal {R} = 60{\rm \,arcsec}, \ell = 5$$ shows a 13 arcsec offset between the maximum value peak and the centroid, indicating either a very skewed or dimodal structure. Figure 5. View largeDownload slide Independent AIM-derived shear data are used in a cluster-scale mass field reconstruction of Abell 2744 and a noise reconstruction in its associated parallel field (both up to a mass-sheet degeneracy scaling). The top left and right figures show the convergence field resulting from a 15 arcsec smoothing kernal applied to the cluster and parallel fields, respectively, while the bottom left and right figures display the results while using a 30 arcsec kernal instead. The axes are in units of arcseconds. Figure 5. View largeDownload slide Independent AIM-derived shear data are used in a cluster-scale mass field reconstruction of Abell 2744 and a noise reconstruction in its associated parallel field (both up to a mass-sheet degeneracy scaling). The top left and right figures show the convergence field resulting from a 15 arcsec smoothing kernal applied to the cluster and parallel fields, respectively, while the bottom left and right figures display the results while using a 30 arcsec kernal instead. The axes are in units of arcseconds. 4.2.1 Convergence-flexion inversion reconstructions There is a degree to which the smoothing scale for a mass reconstruction is somewhat arbitrary, though we typically set it so that the noise does not dominate. Further, when turning a reconstruction into concrete mass estimates, there is bound to be some ambiguity. For instance, the mass sheet degeneracy allows us to adjust the cluster up and down, and is dependent on the boundary conditions(e.g. Bradač, Lombardi & Schneider 2004). However, the effect on our measured peak mass within an aperture is limited – for example, increasing or reducing the assumed boundary convergence by κ ± 0.1 results in a contained mass difference of ∓0.04 × 1014 h−1 M⊙ for our primary structure, or 2 per cent of our calculated mass, with similar results (∼3 per cent) for the secondary and tertiary peaks. We focus on a smoothing scale of 15 arcsec (except as discussed below), and focus on comparisons to the literature within particular apertures. Our reconstructed peak positions and ellipticities, compared with those found previous researchers, can be found in Table 1, while the estimated masses within apertures (selected for comparison with previous researchers) can be found in Table 2. Table 1. Detected substructure geometries. The listed structural properties are derived from the moments of each peak's associated pixels. Included pixels are determined through nearest-neighbour proximity and have a minimum threshold of 2.5σ above the peak-subtracted background level. The detection geometries from three different reconstruction methods are compared: flexion aperture mass inversion, $$\mathcal {F}_{{\rm aper}}{}$$; direct flexion-convergence inversion, $$\kappa _{{\mathcal {F}}{}}$$; and optimization of an SIS model to measured flexion data. Coordinates are specified in arcseconds and relative to the centre of our HST pixel image, at α = 3.589 0213, δ = −30.397 196 (deg, J2000.0). Method/  Xpeak  Ypeak  a  b  e  ϕ  Peak  structure  (arcsec)  (arcsec)  (arcsec)  (arcsec)    (deg)  value  $$\quad {\mathcal {F}}{}_{{\rm aper}}$$            S/N  (1) Core  11.47 ± 13.3  −9.98 ± 17.1  31.36 ± 3.6  16.13 ± 1.3  0.58 ± 0.1  67.05 ± 8.5  3.50 ± 0.84  (2) NE  −55.15 ± 6.7  72.26 ± 14.0  18.03 ± 1.8  9.78 ± 0.6  0.55 ± 0.1  105.04 ± 6.3  2.68 ± 0.85  (3) E  −68.26 ± 9.9  −0.59 ± 10.0  17.16 ± 1.9  10.47 ± 0.8  0.46 ± 0.1  136.00 ± 8.0  2.34 ± 0.73  $$\quad \kappa _{{\mathcal {F}}{}}$$            $$\kappa \qquad$$  (1) Core  3.82 ± 4.0  −1.76 ± 4.4  30.1 ± 0.3  19.19 ± 0.1  0.42 ± 0.01  54.02 ± 0.6  2.50 ± 0.29  (2) NE  –  –  –  –  –  –  –  (3) E  −61.7 ± 2.3  −5.28 ± 2.4  26.0 ± 0.2  23.86 ± 0.2  0.08 ± 0.02  71.02 ± 2.9  2.65 ± 0.30  Lensmodel (SIS)            σν (km s−1)  (1) Core  3.14 ± 4.3  −11.97 ± 6.5  –  –  –  –  1904 ± 346  (2a) NE  −47.16 ± 2.2  41.41 ± 3.4  –  –  –  –  747 ± 380  (2b) NE  −35.99 ± 2.7  68.83 ± 4.2  –  –  –  –  594 ± 690  (3a) E  −52.79 ± 2.5  −21.78 ± 2.2  –  –  –  –  1161 ± 293  (3b) E  −65.44 ± 3.2  0.21 ± 6.2  –  –  –  –  881 ± 552   κγ            $$\kappa \qquad$$  (1) Core  3.82 ± 1.99  −8.81 ± 4.57  50.22 ± 2.58  37.51 ± 1.29  0.28 ± 0.06  0.50 ± 81.6  0.32 ± 0.03  (1) NE  −40.95 ± 3.33  55.81 ± 3.91  6.71 ± 0.24  3.06 ± 0.11  0.65 ± 0.03  128.74 ± 2.20  0.036 ± 0.024  Method/  Xpeak  Ypeak  a  b  e  ϕ  Peak  structure  (arcsec)  (arcsec)  (arcsec)  (arcsec)    (deg)  value  $$\quad {\mathcal {F}}{}_{{\rm aper}}$$            S/N  (1) Core  11.47 ± 13.3  −9.98 ± 17.1  31.36 ± 3.6  16.13 ± 1.3  0.58 ± 0.1  67.05 ± 8.5  3.50 ± 0.84  (2) NE  −55.15 ± 6.7  72.26 ± 14.0  18.03 ± 1.8  9.78 ± 0.6  0.55 ± 0.1  105.04 ± 6.3  2.68 ± 0.85  (3) E  −68.26 ± 9.9  −0.59 ± 10.0  17.16 ± 1.9  10.47 ± 0.8  0.46 ± 0.1  136.00 ± 8.0  2.34 ± 0.73  $$\quad \kappa _{{\mathcal {F}}{}}$$            $$\kappa \qquad$$  (1) Core  3.82 ± 4.0  −1.76 ± 4.4  30.1 ± 0.3  19.19 ± 0.1  0.42 ± 0.01  54.02 ± 0.6  2.50 ± 0.29  (2) NE  –  –  –  –  –  –  –  (3) E  −61.7 ± 2.3  −5.28 ± 2.4  26.0 ± 0.2  23.86 ± 0.2  0.08 ± 0.02  71.02 ± 2.9  2.65 ± 0.30  Lensmodel (SIS)            σν (km s−1)  (1) Core  3.14 ± 4.3  −11.97 ± 6.5  –  –  –  –  1904 ± 346  (2a) NE  −47.16 ± 2.2  41.41 ± 3.4  –  –  –  –  747 ± 380  (2b) NE  −35.99 ± 2.7  68.83 ± 4.2  –  –  –  –  594 ± 690  (3a) E  −52.79 ± 2.5  −21.78 ± 2.2  –  –  –  –  1161 ± 293  (3b) E  −65.44 ± 3.2  0.21 ± 6.2  –  –  –  –  881 ± 552   κγ            $$\kappa \qquad$$  (1) Core  3.82 ± 1.99  −8.81 ± 4.57  50.22 ± 2.58  37.51 ± 1.29  0.28 ± 0.06  0.50 ± 81.6  0.32 ± 0.03  (1) NE  −40.95 ± 3.33  55.81 ± 3.91  6.71 ± 0.24  3.06 ± 0.11  0.65 ± 0.03  128.74 ± 2.20  0.036 ± 0.024  View Large Table 2. Sky position and enclosed mass for substructure identified in this work, Merten et al. (2011), Medezinski et al. (2016), and Jauzac et al. (2016b). The coordinates Xpeak and Ypeak are specified in arcseconds and are relative to the centre of our HST pixel image, at α = 3.589 0213, δ = −30.397 196 (deg, J2000.0). These refer to the positions of either the peak maximum value (as in this work), or the peak centroid (as in Merten et al. 2011), as reported in each respective analysis. In order to directly compare our mass and size estimates with those in the relevant literature we assume here the dimensionless Hubble parameter h = 0.7. Uncertainties in coordinates and enclosed masses for the different techniques used in this work, $$\mathcal {F}_{{\rm aper}}{}, \kappa _{\mathcal {F}}{}, \kappa _{\gamma }, {\rm SIS}_{(1)}, {\rm SIS}_{(2)}$$, are derived from the variance of peak-associated pixel values throughout independent sets of bootstrap realizations as discussed in the text. Errors are derived from the standard deviation of parameter values as calculated for the selected pixels applied to each bootstrap realization.   Source  Xpeak (arcsec)  Ypeak (arcsec)  Renclosing (kpc)  M( < Renclosing) (1013 M⊙)  Core  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  11.47 ± 13.3  −9.98 ± 17.1   n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  3.82 ± 4.0  −1.76 ± 4.4  150  16.2 ± 1.2    κγ (this work)  3.81 ± 2.0  −8.81 ± 4.6  150 (250)  6.27 ± 0.35 (13.3 ± 0.8)    SIS (this work)  3.14 ± 4.3  −11.97 ± 6.5  150  19.9 ± 7.2    Jauzac et al. (2016b)  8.58  −10.72  150  13.55 ± 0.09    Medezinski et al. (2016)  12.91 ± 10.8  −15.84 ± 9.6  200c  77 ± 34    Merten et al. (2011)  $$-2.24^{+11}_{-8}$$  $$-16.71^{+6}_{-13}$$  250  22.4 ± 5.5  E  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −68.26 ± 9.9  −0.59 ± 10.0  $$\varnothing$$  $$\varnothing$$    $$\kappa _{\mathcal {F}}{}$$ (this work)  −61.7 ± 2.3  −5.28 ± 2.4  110  8.81 ± 0.52    SIS(1) (this work)  −52.79 ± 2.5  −21.78 ± 2.2  150  7.4 ± 3.7    SIS(2) (this work)  −65.44 ± 3.2  0.21 ± 6.2  150  4.3 ± 5.3    Jauzac et al. (2016b)  ..  ..  ..  ..    Medezinski et al. (2016)  ..  ..  ..  ..    Merten et al. (2011)  ..  ..  ..  ..  NE  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −55.15 ± 6.7  72.26 ± 14.0  n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  −  −  68  2.92 ± 0.26b    SIS(1) (this work)  −47.16 ± 2.2  41.41 ± 3.4  150  3.1 ± 3.1    SIS(2) (this work)  −35.99 ± 2.7  68.83 ± 4.2  150  1.9 ± 4.5    Jauzac et al. (2016b) a   −46.91  81.16  150  5.00 ± 0.40    Medezinski et al. (2016)  −52.02 ± 19.8  78.25 ± 31.8  200c  28 ± 16    Merten et al. (2011)  ×  ×  ×  ×    Source  Xpeak (arcsec)  Ypeak (arcsec)  Renclosing (kpc)  M( < Renclosing) (1013 M⊙)  Core  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  11.47 ± 13.3  −9.98 ± 17.1   n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  3.82 ± 4.0  −1.76 ± 4.4  150  16.2 ± 1.2    κγ (this work)  3.81 ± 2.0  −8.81 ± 4.6  150 (250)  6.27 ± 0.35 (13.3 ± 0.8)    SIS (this work)  3.14 ± 4.3  −11.97 ± 6.5  150  19.9 ± 7.2    Jauzac et al. (2016b)  8.58  −10.72  150  13.55 ± 0.09    Medezinski et al. (2016)  12.91 ± 10.8  −15.84 ± 9.6  200c  77 ± 34    Merten et al. (2011)  $$-2.24^{+11}_{-8}$$  $$-16.71^{+6}_{-13}$$  250  22.4 ± 5.5  E  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −68.26 ± 9.9  −0.59 ± 10.0  $$\varnothing$$  $$\varnothing$$    $$\kappa _{\mathcal {F}}{}$$ (this work)  −61.7 ± 2.3  −5.28 ± 2.4  110  8.81 ± 0.52    SIS(1) (this work)  −52.79 ± 2.5  −21.78 ± 2.2  150  7.4 ± 3.7    SIS(2) (this work)  −65.44 ± 3.2  0.21 ± 6.2  150  4.3 ± 5.3    Jauzac et al. (2016b)  ..  ..  ..  ..    Medezinski et al. (2016)  ..  ..  ..  ..    Merten et al. (2011)  ..  ..  ..  ..  NE  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −55.15 ± 6.7  72.26 ± 14.0  n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  −  −  68  2.92 ± 0.26b    SIS(1) (this work)  −47.16 ± 2.2  41.41 ± 3.4  150  3.1 ± 3.1    SIS(2) (this work)  −35.99 ± 2.7  68.83 ± 4.2  150  1.9 ± 4.5    Jauzac et al. (2016b) a   −46.91  81.16  150  5.00 ± 0.40    Medezinski et al. (2016)  −52.02 ± 19.8  78.25 ± 31.8  200c  28 ± 16    Merten et al. (2011)  ×  ×  ×  ×  Notes.aThe NE structure is consistent with the Jauzac et al. (2016b) S1 substructure. bTo directly compare our derived masses and distances in physical units (M⊙ and kpc), we temporarily adopt the shared common cosmological model with Hubble constant h = 0.7. View Large As the NE peak is not clearly detected at this scale, the NE aperture is set to the centroid of the $$\mathcal {R} = 60{\rm \,arcsec}$$$$\mathcal {F}_{{\rm aper}}{}$$ reconstruction. The aperture sizes were selected such that the they approximately encircled the identified peaks down to a value around the calculated 1σ background level for both the $$\mathcal {F}_{{\rm aper}}{}$$ and κ reconstructions, with consideration for a representative size that could be easily standardized or replicated and that worked reasonably well for each reconstruction type. 4.2.2 Aperture mass reconstructions While the convergence–flexion relation allows us to calculate physical properties of the cluster or peaks inside a given aperture (provided the convergence is normalized to the mass-sheet degeneracy), the aperture mass method gives a more qualitative result of where mass peaks lie. The lower left and lower right graphs in Fig. 3 show the cluster/parallel aperture mass/noise reconstructions using aperture R = 60 arcsec and polynomial index ℓ = 5. Mass S/N contours in the cluster indicate a strong central peak aligned and centred on the Brightest Cluster Galaxy (BCG) and extending along the cluster's luminosity semimajor axis. A corresponding lesser peak is identified to the east along with another to the north-east. The signal tapers to zero near the edge of the data, with the exception of a few structures along the edge. It is important to note, however, that it is not straightforward to compare the significance values of these peak detections between works directly. Our significance values are derived through different means (our aperture mass method) than those of our comparison papers, and therefore the statistical meaning of these sets of values is not necessarily the same. While direct summation techniques (e.g. Fourier estimates of κ) are essential in providing mass estimates of individual structures, we stress the importance of the aperture mass statistic in this work because of its usefulness in identifying individual structures. Some care must be taken with the choice of smoothing radius and order of the polynomial index, however, in that they yield an effective smoothing scale of $$\sim \mathcal {R}_{{\rm aper}} / \ell$$. In Fig. 3, this gives an approximate smoothing scale of ∼12 arcsec. However, in Fig. 4, we see that the location and number of central peak(s) is(are) highly sensitive to these parameters. It is significant to note that a north-eastern local maxima is not clearly identified in our convergence map at the depicted scale size, but is well-pronounced within all explored $$\mathcal {F}_{{\rm aper}}{}$$ parametrizations. This may suggest that the aperture mass statistic might be better suited for substructure detection than a flexion-convergence inversion, by using the data more directly or through more suitable weighting profiles. Further study would be necessary in order to draw valid conclusions in this area, however. 4.2.3 Shear reconstructions We also compare the scales of structural predictions between flexion and shear by producing convergence reconstructions through direct shear inversion. Although our treatment of the PSF (using a simple Gaussian) would likely introduce ellipticity bias in an analysis on much wider scales, we expect a much stronger shear signal in the inner 3.5 arcmin × 3.5 arcmin cluster field. Fig. 5 shows our convergence reconstructions using just the averaged galaxy ellipticities measured through AIM fitting, again up to a mass-sheet degeneracy. Following a similar format, the top two graphs show the cluster and parallel field κ maps using a 15 arcsec shear smoothing kernel, while the bottom graphs utilize a 30 arcsec shear kernel. The clear broad central structure of the cluster peak at each scale is apparent, and both estimate the cluster centroid to be very close to what the flexion maps predict. However, both shear smoothing scales do not show any strong indications of local substructure. In the 15 arcsec map, the convergence extends slightly towards and hints at the eastern and north-eastern structures highlighted by our flexion analysis, but does not reproduce the flexion map's well-defined structural peaks, while the 30 arcsec reconstruction shows no indication of possible local structure at all. Employing our peak-finding algorithm to the 15 arcsec shear-convergence reconstruction with a σdetect = 2.5 detection threshold and a 1σbg background pixel-association threshold, we found the resulting region to have semimajor axis a = 50.22 ± 2.58 and semiminor axis b = 37.51 ± 1.29. As noted, the geometry values and uncertainties here refer to the reconstructed peak as detected in the κ map as reconstructed from 15 arcsec-smoothed shear data. Through the SL+WL analysis in Jauzac et al. (2016b), a WL model centred on the cluster core and SL region was found to have semimajor axis a = 60 arcsec, semiminor axis b = 42 arcsec, and orientation angle θ = 60°. This result is quite similar to our own reconstruction that finds an axis ratio of approximately 2:1, with an orientation angle of approximately 67°. To test the PSF assumptions and reconstruction method, we produced noise reconstructions of the shear in the parallel field in the same manner as the flexion procedure. Neither parallel field shear smoothing scale shows any appreciable structural significance, particularly when compared with the strong central profile of the cluster field, and the overall noise signal is centred at 0.0 similar to the flexion parallel field reconstructions. 4.2.4 Parallel reconstruction The parallel field convergence noise reconstruction under the same parameters is shown in the top right graph of Fig. 3. Compared to the cluster mass reconstruction, the derived magnitudes are larger but the noise peaks are much broader and do not correlate with any luminous matter. In addition, the magnitude of the highest value is approximately the same as that of the lowest, unlike our derived cluster convergence map. As we expect the parallel field to have 〈κ〉 ∼ 0, the observed |κ| parity around 0.0 is an indicator of the magnitude of systematic noise in an unconstrained flexion–convergence reconstruction. The parallel field map does not appear to have any significant peaks, beyond apparent edge effects (specifically at the north and northeast edges), nor do the minor peaks detected within the field (S/N ∼ 1) correlate with any luminous structure present. 4.3 Comparing flexion reconstructions to previous results One of the reasons that the galaxy cluster Abell 2744 was chosen to be observed as part of the HFF programme was that there already existed a wealth of observational information and a substantial amount of work detailed in the literature. One of the advantages of this strong existing framework is that it allows new methods and science to rest on, build from, and compare to studies and techniques that have already withstood academic rigour. While flexion signal would ideally be incorporated into larger and more encompassing combined approaches, it is also important to understand at a practical level how well flexion alone can reproduce the results of previous work on detecting complex and extended mass distributions. We compare our structure analysis with three works in particular: Merten et al. (2011), whose SL+WL reconstruction first detected four substructures arrayed around the core of the cluster; a later study by Medezinski et al. (2016) that independently detects four substructure components as well in a purely WL analysis, although not all structures detected by Merten were detected by Medezinski, and vice versa; and Jauzac et al. (2016b) who combine the techniques and results of the previous studies of Jauzac et al. (2015b) and Eckert et al. (2015) with additional information to form a broad view of the currently-known mass structure of Abell 2744, the summary of which can be found in Tables 1 and 2. Comparing with other estimates of Abell 2744, in a SL analysis Jauzac et al. (2015b) found a contained mass of M(<200 kpc) = (2.162 ± 0.005) × 1014 M⊙ within a 200 kpc aperture, which is similar to our estimate of the central potential. Although our relatively simple and straightforward mass calculation does not lie within the error bars of the precise SL estimate, Jauzac et al. stress that the precision of the mass models from gravitational-lensing studies depends strongly on the mass modelling technique, and mass estimates from different groups using different SL algorithms will find different results. As there are only a handful of flexion analyses of real galaxy data up to this point, it is difficult to determine if flexion reconstructions consistently over- or underestimate the mass of detected structure. Okura et al. (2008) obtain a convergence map of A1689 but do not estimate any enclosed masses. Cain et al. (2011) use aperture mass signal exclusively for reconstruction. Leonard et al. (2011) models aperture peaks with NFW and SIS models and obtain M200 estimates for these, but do not directly compare these masses with masses for those peaks attained in other works. As higher-order distortion measurements become more standard in gravitational lensing analyses we expect a clearer view of any biases to emerge. Reconstructions differ not only in terms of mass but also in terms of number of substructures. Using a joint SL+WL+X-ray analysis, Jauzac et al. (2016b) detect eight separate substructure peaks over the range of both the cluster core centre as well as the extended field that is outside of the observed HFF range, while our analysis only identifies three. Likewise, Merten et al. (2011) detect four distinct substructures within a field size of 600 arcsec × 600 arcsec, though the statistical significance of each is not made clear. While considering the question of matching the presence or absence of our detected substructure with other work, we identify a third major peak in our analysis that we previously ignored for its proximity to our field edge. However, the peak is over 25 arcsec from the edge, and the structure's had an elongated radial extent of over 50 arcsec, which brings a large portion of it well within an acceptable field of view. A particular motivation for the inclusion is the spatial correlation between our detected peak and an SL peak as identified in both Medezinski et al. (2016) and Jauzac et al. (2016b). 5 DISCUSSION We close with a few thoughts related to the differences between our work and previous analyses of A2744 (and other clusters). There are a few reasons why flexion analyses may not be able to pick up previously known mass peaks at any particular position, with the most prominent being that there are simply no (or very few) galaxy sources that both are close to the mass substructure and are suitable sources (e.g. large apparent size and good and realistic model fit) from which to measure a strong flexion signal. Because the flexion field is most prominent close to mass peaks and is buried by inherent shape noise farther away, if there are no chance source detections within a few radial arcseconds of the peak then it is likely that smaller substructure (on the 1012–1013 M⊙ scale) will remain undetected, without a substantial improvement in flexion measurement capability. Besides the structure correlated with potential DM peaks, our mass reconstructions also appear to capture a large portion of the gas distribution and X-ray surface brightness as examined and investigated independently in M11, M16, J16 using Chandra and XMM–Newton maps. The three analyses show generally consistent X-ray luminosity structure, with the global X-ray peak located just SW of the core BCG that expands broadly to the north-west. This structural orientation gives evidence towards the conclusion that the cluster had undergone a merger along the SE-NW axis that is supported by all three works along with other types of supporting evidence. In our own analysis, our derived mass distribution contours run parallel with J16's X-ray contours on a slanted line running from the core BCG NE towards the group of cluster galaxies located at approximately α = 3.596 5217°, δ = −30.388 257° (deg, J2000). While the gravitational shear field can be used effectively to determine the overall mass structure of galaxy clusters, its extended nonlocal effects as well as its inherent ellipticity-degeneracy limit its use to broad mass distributions, and thus it is not a viable candidate for higher resolution substructure detections. Additionally, while SL can and has led to precise characterizations of inner cluster cores, multiply-imaged source systems are not guaranteed to be located near substructure and quickly become sparse outside the dense core. As an intermediary signal spanning the strong and WL regimes, gravitational flexion signal has the ability to effectively probe significant cluster substructure on scales and at angular extents that cannot be practically detected through other means. As inherent flexion noise and systematic bias become more well-understood, flexion signal has the potential to be a key component in both exploring the behaviour of galaxy cluster formation and evolution as well as understanding the nature of DM structural dynamics. In this work, we have used an AIM implementation to measure flexion signal in the Abell 2744 galaxy cluster and inherent flexion in Abell 2744's associated parallel field. We show the efficacy of using flexion alone as an indicator of structure, exploring a much deeper view into the inner core of the cluster than shear would allow, and investigate the role of different mass estimators in both the cluster and parallel field. We identify and obtain mass estimates for both the central core of the cluster and a detected substructure offset 1.43 arcsec to the east of the core. Finally, we demonstrate that we are able to make simultaneous measurements of the shear field while measuring flexion through the AIM technique, and reconstruct the broader cluster mass structure while finding no such signal in the parallel field. Acknowledgements J.P.B. would like to thank both Austen Groener and Markus Rexroth for valuable discussions and dialog throughout this project. He would also like to recognize Mike Jarvis for his early coordination and support of the conceptual development and integration of the GalFlex flexion module as an affiliate into the GalSim software package. Finally, we thank the anonymous referee and editors for their constructive critiques, beneficial suggestions, and patience. The image data used in this paper are based on observations obtained with the NASA/ESA Hubble Space Telescope, retrieved from the MAST at the STScI. STScI is operated by the Association of Universities for Research in Astronomy, Inc. under NASA contract NAS 5-26555. This work also utilizes a subset of the gravitational lensing models produced by PIs Bradac, Natarajan & Kneib (CATS), Merten & Zitrin, Sharon, Williams, and the GLAFIC and Diego groups. These lens models have been partially funded by the HST Frontier Fields programmes conducted by STScI and obtained through MAST, as described above. Footnotes 1 https://archive.stsci.edu/pub/hlsp/frontier/abell2744/images/hst/ 2 http://www.stsci.edu/hst/observatory/focus/TinyTim 3 https://code.google.com/archive/p/astrolibpy/ 4 http://physics.drexel.edu/jbird/galflex/index.html REFERENCES Andrae R., Schulze-Hartung T., Melchior P., 2010, preprint (arXiv:1012.3754) Bacon D. J., Massey R. J., Refregier A. R., Ellis R. S., 2003, MNRAS , 344, 673 https://doi.org/10.1046/j.1365-8711.2003.06877.x CrossRef Search ADS   Bacon D. J., Goldberg D. M., Rowe B. T. P., Taylor A. N., 2006, MNRAS , 365, 414 https://doi.org/10.1111/j.1365-2966.2005.09624.x CrossRef Search ADS   Bartelmann M., Schneider P., 2001, Phys. Rep. , 340, 291 https://doi.org/10.1016/S0370-1573(00)00082-X CrossRef Search ADS   Bertin E., Arnouts S., 1996, A&AS , 117, 393 CrossRef Search ADS   Blandford R. D., Narayan R., 1992, ARA&A , 30, 311 https://doi.org/10.1146/annurev.aa.30.090192.001523 CrossRef Search ADS   Blanton M. R., Moustakas J., 2009, ARA&A , 47, 159 https://doi.org/10.1146/annurev-astro-082708-101734 CrossRef Search ADS   Bradac M., Schneider P., Lombardi M., Erben T., 2005, A&A , 437, 39 Bradač M., Lombardi M., Schneider P., 2004, A&A , 424, 13 https://doi.org/10.1051/0004-6361:20035744 CrossRef Search ADS   Bridle S. et al.  , 2010, MNRAS , 405, 2044 Cain B., Schechter P. L., Bautz M. W., 2011, ApJ , 736, 43 https://doi.org/10.1088/0004-637X/736/1/43 CrossRef Search ADS   Cain B., Bradač M., Levinson R., 2016, MNRAS  Caminha G. B. et al.  , 2016a, A&A , 587, A80 https://doi.org/10.1051/0004-6361/201527670 CrossRef Search ADS   Caminha G. B. et al.  , 2016b, A&A , 595, A100 https://doi.org/10.1051/0004-6361/201527995 CrossRef Search ADS   Caminha G. B. et al.  , 2017, A&A , 600, A90 https://doi.org/10.1051/0004-6361/201629297 CrossRef Search ADS   Capaccioli M., 1989, in Corwin H. G. Jr, Bottinelli L., eds, World of Galaxies (Le Monde des Galaxies) . Springer-Verlag, Berlin, p. 208 Google Scholar CrossRef Search ADS   Cardone V. F., Vicinanza M., Er X., Maoli R., Scaramella R., 2016, MNRAS , 462, 4028 https://doi.org/10.1093/mnras/stw1803 CrossRef Search ADS   Diego J. M., Broadhurst T., Molnar S. M., Lam D., Lim J., 2015a, MNRAS , 447, 3130 https://doi.org/10.1093/mnras/stu2660 CrossRef Search ADS   Diego J. M., Broadhurst T., Zitrin A., Lam D., Lim J., Ford H. C., Zheng W., 2015b, MNRAS , 451, 3920 https://doi.org/10.1093/mnras/stv1168 CrossRef Search ADS   Diego J. M. et al.  , 2016a, MNRAS , 473, 4279 https://doi.org/10.1093/mnras/stx2609 CrossRef Search ADS   Diego J. M., Broadhurst T. W J., Silk J., Lim J., Zheng W., Lam D., Ford H., 2016b, MNRAS , 459, 3447 https://doi.org/10.1093/mnras/stw865 CrossRef Search ADS   Diemand J., Kuhlen M., Madau P., Zemp M., Moore B., Potter D., Stadel J., 2008, Nature , 454, 735 https://doi.org/10.1038/nature07153 CrossRef Search ADS PubMed  Eckert D. et al.  , 2015, Nature , 528, 105 https://doi.org/10.1038/nature16058 CrossRef Search ADS PubMed  Er X., Tereno I., Mao S., 2012, MNRAS , 421, 1443 https://doi.org/10.1111/j.1365-2966.2012.20408.x CrossRef Search ADS   Falco E. E., Gorenstein M. V., Shapiro I. I., 1985, ApJ , 289, L1 https://doi.org/10.1086/184422 CrossRef Search ADS   Gao L., De Lucia G., White S. D. M., Jenkins A., 2004, MNRAS , 352, L1 https://doi.org/10.1111/j.1365-2966.2004.08098.x CrossRef Search ADS   Goldberg D. M., Bacon D. J., 2005, ApJ , 619, 741 https://doi.org/10.1086/426782 CrossRef Search ADS   Goldberg D. M., Leonard A., 2007, ApJ , 660, 1003 https://doi.org/10.1086/513137 CrossRef Search ADS   Goldberg D. M., Natarajan P., 2002, ApJ , 564, 65 https://doi.org/10.1086/324202 CrossRef Search ADS   Gorenstein M. V., Shapiro I. I., Falco E. E., 1988, ApJ , 327, 693 https://doi.org/10.1086/166226 CrossRef Search ADS   Grillo C. et al.  , 2015, ApJ , 800, 38 https://doi.org/10.1088/0004-637X/800/1/38 CrossRef Search ADS   Grillo C. et al.  , 2016, ApJ , 822, 78 https://doi.org/10.3847/0004-637X/822/2/78 CrossRef Search ADS   Harvey D., Massey R., Kitching T., Taylor A., Tittley E., 2015, Science , 347, 1462 https://doi.org/10.1126/science.1261381 CrossRef Search ADS PubMed  Heymans C. et al.  , 2006, MNRAS , 368, 1323 https://doi.org/10.1111/j.1365-2966.2006.10198.x CrossRef Search ADS   Hoekstra H., Bartelmann M., Dahle H., Israel H., Limousin M., Meneghetti M., 2013, Space Sci. Rev. , 177, 75 https://doi.org/10.1007/s11214-013-9978-5 CrossRef Search ADS   Jauzac M. et al.  , 2014, MNRAS , 443, 1549 https://doi.org/10.1093/mnras/stu1355 CrossRef Search ADS   Jauzac M. et al.  , 2015a, MNRAS , 446, 4132 https://doi.org/10.1093/mnras/stu2425 CrossRef Search ADS   Jauzac M. et al.  , 2015b, MNRAS , 452, 1437 https://doi.org/10.1093/mnras/stv1402 CrossRef Search ADS   Jauzac M. et al.  , 2016a, MNRAS , 457, 2029 https://doi.org/10.1093/mnras/stw069 CrossRef Search ADS   Jauzac M. et al.  , 2016b, MNRAS , 463, 3876 https://doi.org/10.1093/mnras/stw2251 CrossRef Search ADS   Johnson T. L., Sharon K., Bayliss M. B., Gladders M. D., Coe D., Ebeling H., 2014, ApJ , 797, 48 https://doi.org/10.1088/0004-637X/797/1/48 CrossRef Search ADS   Kaiser N., Squires G., 1993, ApJ , 404, 441 https://doi.org/10.1086/172297 CrossRef Search ADS   Kaiser N., Squires G., Broadhurst T., 1995, ApJ , 449, 460 https://doi.org/10.1086/176071 CrossRef Search ADS   Kitching T. D. et al.  , 2012, MNRAS , 423, 3163 https://doi.org/10.1111/j.1365-2966.2012.21095.x CrossRef Search ADS   Klypin A., Kravtsov A. V., Valenzuela O., Prada F., 1999, ApJ , 522, 82 https://doi.org/10.1086/307643 CrossRef Search ADS   Kneib J.-P., Natarajan P., 2011, A&AR , 19, 47 CrossRef Search ADS   Kravtsov A. V., Berlind A. A., Wechsler R. H., Klypin A. A., Gottlöber S., Allgood B., Primack J. R., 2004a, ApJ , 609, 35 https://doi.org/10.1086/420959 CrossRef Search ADS   Kravtsov A. V., Gnedin O. Y., Klypin A. A., 2004b, ApJ , 609, 482 https://doi.org/10.1086/421322 CrossRef Search ADS   Lam D., Broadhurst T., Diego J. M., Lim J., Coe D., Ford H. C., Zheng W., 2014, ApJ , 797, 98 https://doi.org/10.1088/0004-637X/797/2/98 CrossRef Search ADS   Leonard A., Goldberg D. M., Haaga J. L., Massey R., 2007, ApJ , 666, 51 https://doi.org/10.1086/520109 CrossRef Search ADS   Leonard A., King L. J., Goldberg D. M., 2011, MNRAS , 413, 789 https://doi.org/10.1111/j.1365-2966.2010.18171.x CrossRef Search ADS   Lewis A., Challinor A., 2006, Phys. Rep. , 429, 1 https://doi.org/10.1016/j.physrep.2006.03.002 CrossRef Search ADS   Lotz J. M. et al.  , 2017, ApJ , 837, 97 https://doi.org/10.3847/1538-4357/837/1/97 CrossRef Search ADS   Mandelbaum R., Seljak U., Kauffmann G., Hirata C. M., Brinkmann J., 2006, MNRAS , 368, 715 https://doi.org/10.1111/j.1365-2966.2006.10156.x CrossRef Search ADS   Mandelbaum R. et al.  , 2015, MNRAS , 450, 2963 https://doi.org/10.1093/mnras/stv781 CrossRef Search ADS   Massey R. et al.  , 2007a, MNRAS , 376, 13 https://doi.org/10.1111/j.1365-2966.2006.11315.x CrossRef Search ADS   Massey R., Rowe B., Refregier A., Bacon D. J., Bergé J., 2007b, MNRAS , 380, 229 https://doi.org/10.1111/j.1365-2966.2007.12072.x CrossRef Search ADS   Massey R., Kitching T., Richard J., 2010, Rep. Prog. Phys. , 73, 086901 https://doi.org/10.1088/0034-4885/73/8/086901 CrossRef Search ADS   Medezinski E. et al.  , 2013, ApJ , 777, 43 https://doi.org/10.1088/0004-637X/777/1/43 CrossRef Search ADS   Medezinski E., Umetsu K., Okabe N., Nonino M., Molnar S., Massey R., Dupke R., Merten J., 2016, ApJ , 817, 24 https://doi.org/10.3847/0004-637X/817/1/24 CrossRef Search ADS   Merten J. et al.  , 2011, MNRAS , 417, 333 https://doi.org/10.1111/j.1365-2966.2011.19266.x CrossRef Search ADS   Mo W. et al.  , 2016, ApJ , 818, L25 https://doi.org/10.3847/2041-8205/818/2/L25 CrossRef Search ADS   Mohammed I., Saha P., Williams L. L. R., Liesenborgs J., Sebesta K., 2016, MNRAS , 459, 1698 https://doi.org/10.1093/mnras/stw727 CrossRef Search ADS   Munshi D., Valageas P., van Waerbeke L., Heavens A., 2008, Phys. Rep. , 462, 67 https://doi.org/10.1016/j.physrep.2008.02.003 CrossRef Search ADS   Nagai D., Kravtsov A. V., 2005, ApJ , 618, 557 https://doi.org/10.1086/426016 CrossRef Search ADS   Okamoto T., Nagashima M., 2001, ApJ , 547, 109 https://doi.org/10.1086/318375 CrossRef Search ADS   Okura Y., Umetsu K., Futamase T., 2008, ApJ , 680, 1 https://doi.org/10.1086/587676 CrossRef Search ADS   Owers M. S., Randall S. W., Nulsen P. E. J., Couch W. J., David L. P., Kempner J. C., 2011, ApJ , 728, 27 https://doi.org/10.1088/0004-637X/728/1/27 CrossRef Search ADS   Press W. H., Schechter P., 1974, ApJ , 187, 425 https://doi.org/10.1086/152650 CrossRef Search ADS   Rexroth M., Natarajan P., Kneib J.-P., 2016, MNRAS , 460, 2505 https://doi.org/10.1093/mnras/stw1017 CrossRef Search ADS   Rhodes J. D. et al.  , 2007, ApJS , 172, 203 https://doi.org/10.1086/516592 CrossRef Search ADS   Richard J. et al.  , 2014, MNRAS , 444, 268 https://doi.org/10.1093/mnras/stu1395 CrossRef Search ADS   Rix H.-W. et al.  , 2004, ApJS , 152, 163 https://doi.org/10.1086/420885 CrossRef Search ADS   Schneider P., 1996, MNRAS , 283, 837 https://doi.org/10.1093/mnras/283.3.837 CrossRef Search ADS   Schneider P., Er X., 2008, A&A , 485, 363 https://doi.org/10.1051/0004-6361:20078631 CrossRef Search ADS   Schneider P., Seitz C., 1995, A&A , 294, 411 Schneider P., van Waerbeke L., Jain B., Kruse G., 1998, MNRAS , 296, 873 https://doi.org/10.1046/j.1365-8711.1998.01422.x CrossRef Search ADS   Schrabback T. et al.  , 2010, A&A , 516, A63 https://doi.org/10.1051/0004-6361/200913577 CrossRef Search ADS   Sebesta K., Williams L. L. R., Mohammed I., Saha P., Liesenborgs J., 2016, MNRAS , 461, 2126 https://doi.org/10.1093/mnras/stw1433 CrossRef Search ADS   Sendra I., Diego J. M., Broadhurst T., Lazkoz R., 2014, MNRAS , 437, 2642 https://doi.org/10.1093/mnras/stt2076 CrossRef Search ADS   Sharon K., Johnson T. L., 2015, ApJ , 800, L26 https://doi.org/10.1088/2041-8205/800/2/L26 CrossRef Search ADS   Sheth R. K., Tormen G., 1999, MNRAS , 308, 119 https://doi.org/10.1046/j.1365-8711.1999.02692.x CrossRef Search ADS   Springel V., White S. D. M., Tormen G., Kauffmann G., 2001, MNRAS , 328, 726 https://doi.org/10.1046/j.1365-8711.2001.04912.x CrossRef Search ADS   Tinker J., Kravtsov A. V., Klypin A., Abazajian K., Warren M., Yepes G., Gottlöber S., Holz D. E., 2008, ApJ , 688, 709 https://doi.org/10.1086/591439 CrossRef Search ADS   Treu T. et al.  , 2016, ApJ , 817, 60 https://doi.org/10.3847/0004-637X/817/1/60 CrossRef Search ADS   Van Waerbeke L. et al.  , 2000, A&A , 358, 30 Viola M., Melchior P., Bartelmann M., 2012, MNRAS , 419, 2215 https://doi.org/10.1111/j.1365-2966.2011.19872.x CrossRef Search ADS   Wang X. et al.  , 2015, ApJ , 811, 29 https://doi.org/10.1088/0004-637X/811/1/29 CrossRef Search ADS   Wittman D. M., Tyson J. A., Kirkman D., Dell'Antonio I., Bernstein G., 2000, Nature , 405, 143 https://doi.org/10.1038/35012001 CrossRef Search ADS PubMed  Zitrin A., Menanteau F., Hughes J. P., Coe D., Barrientos L. F., Infante L., Mandelbaum R., 2013, ApJ , 770, L15 https://doi.org/10.1088/2041-8205/770/1/L15 CrossRef Search ADS   Zitrin A. et al.  , 2015, ApJ , 801, 44 https://doi.org/10.1088/0004-637X/801/1/44 CrossRef Search ADS   © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Monthly Notices of the Royal Astronomical Society Oxford University Press

# Flexion in Abell 2744

, Volume 476 (1) – May 1, 2018
15 pages

/lp/ou_press/flexion-in-abell-2744-A8CNZwVlV0
Publisher
The Royal Astronomical Society
ISSN
0035-8711
eISSN
1365-2966
D.O.I.
10.1093/mnras/sty300
Publisher site
See Article on Publisher Site

### Abstract

Abstract We present the first flexion-focused gravitational lensing analysis of the Hubble Frontier Field observations of Abell 2744 (z = 0.308). We apply a modified Analytic Image Model technique to measure source galaxy flexion and shear values at a final number density of 82 arcmin−2. By using flexion data alone, we are able to identify the primary mass structure aligned along the heart of the cluster in addition to two major substructure peaks, including an NE component that corresponds to previous lensing work and a new peak detection offset 1.43 arcmin from the cluster core towards the east. We generate two types of non-parametric reconstructions: flexion aperture mass maps, which identify central core, E, and NE substructure peaks with mass signal-to-noise contours peaking at 3.5σ, 2.7σ, and 2.3σ, respectively; and convergence maps derived directly from the smoothed flexion field. For the primary peak, we find a mass of (1.62 ±  0.12) × 1014 h−1 M⊙ within a 33 arcsec (105 h−1 kpc) aperture, a mass of (2.92 ±  0.26) × 1013 h−1 M⊙ within a 16 arcsec (50 h−1 kpc) aperture for the north-eastern substructure, and (8.81 ± 0.52) × 1013 h−1 M⊙ within a 25 arcsec (80 h−1 kpc) aperture for the novel eastern substructure. gravitational lensing: weak, galaxies: clusters: general, galaxies: clusters: individual: A2744, dark matter, cosmology: observations 1 INTRODUCTION Over the last two decades, gravitational lensing has served as one of the most successful tools in the development of our standard model of cosmology and extragalactic astronomy (e.g. Blandford & Narayan 1992; Kaiser & Squires 1993; Kaiser, Squires & Broadhurst 1995; Klypin et al. 1999; Van Waerbeke et al. 2000; Wittman et al. 2000; Bartelmann & Schneider 2001; Kravtsov et al. 2004a; Lewis & Challinor 2006; Mandelbaum et al. 2006; Munshi et al. 2008; Tinker et al. 2008; Massey, Kitching & Richard 2010; Kneib & Natarajan 2011; Hoekstra et al. 2013). Galaxy clusters are ideal candidates for observing gravitational lensing effects in background sources, as their high-mass density and broad angular extent ensure a large population of lensed images. While most weak lensing (WL) reconstructions focus on questions of total mass, radial profiles, and the measurement of cluster ellipticities, it is well-known from both numerical simulations (e.g. Klypin et al. 1999; Diemand et al. 2008) and semi-analytical work (Press & Schechter 1974; Sheth & Tormen 1999) that significant substructure is expected in dark matter (DM) haloes at approximately the 2 per cent level. Consequently, the advent of larger and deeper data sets now allows precise and accurate mass reconstructions in both the inner core and outer regions of galaxy clusters through strong lensing (SL) and combined lensing approaches (e.g. Bradac et al. 2005; Leonard et al. 2007; Jauzac et al. 2015a). With improved data sets, new methods, and the knowledge and results gleaned through a multitude of shear testing programmes and challenges (e.g. Heymans et al. 2006; Massey et al. 2007a; Bridle et al. 2010; Kitching et al. 2012; Mandelbaum et al. 2015), the last few years have seen an explosion of new strong (+weak) gravitational lensing analyses, particularly through the data obtained in the recent Hubble Frontier Fields (HFF) initiative (e.g. Medezinski et al. 2013, 2016; Zitrin et al. 2013, 2015; Jauzac et al. 2014, 2015a,b, 2016a,b; Johnson et al. 2014; Lam et al. 2014; Richard et al. 2014; Sendra et al. 2014; Diego et al. 2015a,b, 2016a,b; Grillo et al. 2015, 2016; Sharon & Johnson 2015; Wang et al. 2015; Caminha et al. 2016a,b, 2017; Mo et al. 2016; Mohammed et al. 2016; Sebesta et al. 2016; Treu et al. 2016; Lotz et al. 2017). Although highly effective where possible, SL is limited to regions where the lens projected mass approaches or exceeds a critical density, which extends only a few tens of arcseconds from the cluster centre. Conversely, WL techniques can generally be applied outside of this region, but must smooth over a large-enough region to encompass the number of galaxies required for reasonable signal, making detailed shape information of detectable substructure difficult to obtain. In general, it is very difficult to measure precise shape information for galaxy cluster substructure or even be able to detect substructure in the (1012–1013) M⊙ regime at all, particularly so when sufficiently removed from the centre and/or in complex or merging systems. Simulations of ΛCDM cluster formation suggest that substructure in even well-behaved clusters can be significantly extended away from the central core. Numerous N-body results indicate that the radial distribution of DM subhaloes in galaxy clusters is more extended than that of luminous galaxies (Okamoto & Nagashima 2001; Springel et al. 2001; Gao et al. 2004; Kravtsov, Gnedin & Klypin 2004b), while Nagai & Kravtsov (2005) in agreement with previous studies find that the radial distribution of subhaloes is significantly less concentrated than the overall cluster DM distribution. The subhaloes near the cluster centre lose more than 70 per cent of their initial mass due to tidal stripping, while haloes at much farther distances (approaching the virial radius) lose only ∼ 30 per cent on average, resulting in a subhalo population that is biased away from the central cluster regions. On scales relevant to lensing, Cain, Bradač & Levinson (2016) indicates that at least half of these DM subhaloes would exist outside of the Einstein radius of the main cluster halo, significantly impacting the probability of a multiply-imaged system occurring. We propose that one of the most-promising untapped ways to uncover the detailed structure of dynamic clusters is to look at the higher-order lensing statistic known as ‘flexion’ (Goldberg & Natarajan 2002; Goldberg & Bacon 2005; Bacon et al. 2006). Flexion induces arcing effects in lensed galaxies – similar to but smaller in scale than the giant arcs of SL – forming ‘arclets’ near local mass overdensities. The advantage to using this higher-order information is that it provides a new avenue to detect local substructure independently of shear or SL measurements. Flexion signal has the potential to be much more sensitive to small-scale perturbations of the convergence field than shear while remaining theoretically viable outside the SL region, with utility directly proportional to galaxy-image-substructure proximity. Flexion as an indicator of local structure will provide significantly more information if there is a high density of sources, and thus the HFFs are in a unique position to calibrate flexion measurement and reconstruction techniques for other, wider surveys and applications. In this work, we present a gravitational lensing flexion analysis of the widely-studied galaxy cluster Abell 2744. Our paper is structured accordingly. Section 2 provides a background of the lensing formalism used in this work along with information concerning the data set, observational pipeline, and data reduction process. In Section 3, we summarize our analysis methodology, including details and implementation of flexion measurement and our mass estimators. Section 4 presents the results of our analysis on Abell 2744 and its associated parallel field. Our derived convergence maps, best-fitting parametrized lens models, and flexion aperture mass signal maps are presented with comparisons of results from other lensing groups, and we conclude with an extended summary and discussion. We adopt the standard conventions of ΛCDM cosmology with Ωm = 0.3, ΩΛ = 0.7, Hubble constant H0 = 100 h km s−1 Mpc−1 unless otherwise specified, and with all magnitudes listed in the AB system. Likewise, we adopt complex notation for denoting directional vectors and define the gradient operator as ∂ = ∂1 + i∂2. 2 BACKGROUND 2.1 Lensing and flexion formalism Gravitational flexion was originally developed as a way to quantify gradients in the dimensionless surface mass density field, κ, through highly distorted, gravitationally lensed ‘arclet’ images (Goldberg & Natarajan 2002; Goldberg & Bacon 2005; Bacon et al. 2006; Schneider & Er 2008). As the flexion field is the third derivative of the underlying gravitational potential, it is far more sensitive than shear to substructure in cluster haloes while not requiring the multiple images necessary for a full SL analysis. Accurately, sampling the flexion field at source positions can give critical information about the local convergence field and ultimately any nearby (potentially non-luminous) structural peaks. To reconstruct the lens-plane convergence, κ, with measurable galaxy-image properties, we can establish a foreground lensing potential ψ through the Poisson equation:   $$\mathrm{\partial} ^*\mathrm{\partial} \psi = 2\kappa$$ (1)with derivatives in angular units in the lens (galaxy cluster) plane. The physical mass of the lens system may be determined by relating the line-of-sight projected surface mass density Σ to the dimensionless convergence field:   $$\Sigma (\boldsymbol {\theta }) = \kappa (\boldsymbol {\theta }) \cdot \Sigma _{{\rm crit}}\,,$$ (2)and the critical mass density defined in the usual way. The deflection potential may now be related to a number of observables, including the complex shear field,   $$\gamma = \gamma _1 + i\gamma _2 = \frac{1}{2}\mathrm{\partial} \mathrm{\partial} \psi ,$$ (3)which induces tangential elongation in a source image. Spin-1 flexion, $${\mathcal {F}}$$, has a natural interpretation as the gradient of the convergence, inducing distortion that ‘skews’ source images (coma), while spin-3 flexion, $${\mathcal {G}}$$, has m = 3 rotational symmetry and acts as the gradient of the shear, producing ‘arcing’ effects (trefoil) in an image (see fig. 1 of Bacon et al. 2006):   \begin{eqnarray} {\mathcal {F}}= |{\mathcal {F}}|{\rm e}^{i\phi } = \frac{1}{2}\mathrm{\partial} \mathrm{\partial} ^*\mathrm{\partial} \psi = \mathrm{\partial} \kappa , \end{eqnarray} (4)  \begin{eqnarray} {\mathcal {G}}= |{\mathcal {G}}|{\rm e}^{i3\phi } = \frac{1}{2}\mathrm{\partial} \mathrm{\partial} \mathrm{\partial} \psi = \mathrm{\partial} \gamma . \end{eqnarray} (5)In the very WL regime (κ ≪ 1 ,  |γ| ≪ 1), a lensing coordinate transformation may be sufficiently approximated through linear combinations of these distortion terms. The regions where flexion signal is expected to provide the most benefit however typically have fairly large convergence values that are a significant fraction of unity, and therefore we must consider the effects that magnification and the mass-sheet degeneracy impose on our lensing measurements. As first recognized in Falco, Gorenstein & Shapiro (1985) and expanded by Gorenstein, Shapiro & Falco (1988) and Schneider & Seitz (1995), every lensing observable defined through differential shape distortions will remain unchanged under application of a fixed mass sheet: $$\kappa _{\lambda }^{\prime } = \lambda \kappa + (1-\lambda ),$$ where λ can take any real non-trivial value. The same is true with shear and Flexion:   \begin{eqnarray} \gamma ^{\prime } = \lambda \gamma \,, \quad {\mathcal {F}}{}^{\prime } = \lambda {\mathcal {F}}{}\,, \quad {\mathcal {G}}{}^{\prime } = \lambda {\mathcal {G}}{} \,. \end{eqnarray} (6)To avoid the difficulties associated with degenerate measurements, it is possible to reformulate the lensing fields into new quantities that are invariant to the mass-sheet degeneracy:   $$g = \frac{\gamma }{1-\kappa }\,, \quad \Psi _1 = \frac{{\mathcal {F}}{}}{1-\kappa }\,, \quad \Psi _3 = \frac{{\mathcal {G}}{}}{1-\kappa }\,.$$ (7)where g is the reduced shear, Ψ1 is the reduced spin-1 flexion, and Ψ3 is the reduced spin-3 flexion (Schneider & Er 2008). In this work, references to shear and flexion will refer to their reduced versions unless otherwise noted. 2.2 Observations and data With ∼140 Hubble Space Telescope(HST) orbits for each of six massive cluster lenses and their respective parallel fields, the goal of the Hubble Frontier Field programme is to probe the early universe at redshifts up to and including the z = (6–8) regime, with a limiting magnitude of magAB ∼ 28.5–29 in seven passbands for a 5σ point source detection in a 0.4 arcsec aperture. The first cluster to be fully imaged was Abell 2744 at a redshift of z = 0.308 with a field size of 3.5 arcmin × 3.5 arcmin. At this redshift, assuming the previously stated cosmology gives a scale of 3.175 h−1 kpc arcsec−1 (or 4.536 kpc arcsec−1 if assuming h = 0.7). The angular scale of the data (0.05 arcsec scale, 0.03 arcsec after drizzling), the abundance of colour information, the considerable mass (e.g. Jauzac et al. (2015b) measure 2.76 × 1014 M⊙ within 250 kpc of the core of Abell 2744), and line-of-sight orientations of the cluster lenses create favourable conditions for strong flexion signal, and it is to this first cluster that we apply our flexion pipeline. We take advantage of the self-calibrated v1.0 high-level science images1 released to the Mikulski Archive for Space Telescope (MAST) archives from the HST Frontier Fields Science Products Team. These mosaics are designed as a final product to be used directly for scientific analysis, having been distortion-corrected, combined, and registered to a common pixel scale for the highest possible quality. We use the ACS F435w, F606w, and F814w filters for colour magnitudes, and make shape measurements in the F814w filter, with total respective filter integration times corresponding to 24, 14, and 46 orbits, respectively. 3 METHODS 3.1 Data reduction Flexion analyses pose unique challenges as compared to WL shear. If we consider the simplified case of a singular isothermal sphere (SIS) lens profile, the magnitudes of the shear and flexion field strengths fall off as 1/θ and 1/θ2, respectively. The much sharper distance scale on which flexion image-distortions begin to become prominent results in a scenario where only a small handful of proximal source images are able to strongly contribute to the detection of each local substructure. The majority of galaxy source images that are not readily adjacent to apparent substructure will generally exhibit only a weak or noisy flexion signal. These weaker and more radially distant sources are able to contribute either statistically or by constraining or excluding the possibility of local substructure through their absence of significant flexion distortion, in a similar manner as used for predicting multiply-lensed image positions in SL analysis. This dual-regime-spanning behaviour of flexion signal is likely to lead to a much more individualistic treatment of flexion sources when compared to shear, as any one galaxy image may be a significant part of a substructure detection. While shear reconstructions may require averaging over hundreds of galaxies to produce a significant signal, a flexion reconstruction might only depend on a handful of significant sources near each local structure, and thus would require a much higher signal-to-noise (S/N) ratio per galaxy. To optimize the quality of our galaxy sources and the accuracy of their associated uncertainty maps in preparation for flexion measurement, further handling and manipulation of the Space Telescope Science Institute (STScI)-calibrated HST images is necessary. The process of identifying source galaxies down to the background level, deblending a dense cluster field, extracting individual galaxies while masking other contaminating sources, and targeting the correct sources for measurement is non-trivial. We create a post-processing pipeline to optimize the specifics of our data creation, combining it with the Source Extractor (Bertin & Arnouts 1996) utility which can be used for identifying background signatures that match specified criteria, namely involving pixel thresholding levels, aperture size, deblending, and noise parameters. Alongside the identification of galaxy objects in the field, there are selection cuts that can be made both before and after measurements are made to decrease possible contaminants in the flexion signal, and which depend on the inherent and apparent properties of each source. For example, flexion signal depends on the gradient of the lens convergence and is therefore very sensitive to the apparent size of the lensed source, causing larger objects to be preferentially selected. Additionally, the intrinsic shape of a galaxy can play a significant role in a measurement of its ‘intrinsic flexion’ – the flexion that would be measured if the object was not gravitationally lensed at all. Irregular galaxies in particular can mimic a strong external flexion signal. As the flexion signal can depend on far fewer objects than shear, data reduction and source selection is an extremely important component of this analysis. Our data reduction strategy involves the following general steps: We crop all images to exclude regions that have poor stacking coverage, decreasing background estimation bias in Source Extractor and increasing performance. Using redshift and magnitude catalogues (Owers et al. 2011; Merten et al. 2011), we position-match and identify any established objects belonging to either the cluster itself or the foreground. We identify the foreground-associated pixels down to background level and replace them with a randomized noise map set to the mean of the local background. We implement a ‘hot-cold’ strategy (Rix et al. 2004; Leonard et al. 2007) running Source Extractor in dual-image mode across all appropriate bands – first clearing the field of large bright objects not previously removed, then targeting the expected source population with a minimum of 15 pixels over a 2σ detection threshold. As spectroscopic and photometric redshifts are available for only a fraction of our detected objects, we also use a series of criteria to exclude probable source contaminants or low S/N objects. Namely, we reject bright or large galaxies with mag F814w < 24 or full width at half-maximum (FWHM) > 0.9 arcsec as probable member or foreground contaminants, establish a low signal bound by excluding galaxies with mag F814w > 28, FWHM < (2 × FWHM PSF), or S/N ratio <20. Additionally, we remove objects which are flagged as incomplete (FLAG ≥8) in the Source Extractor catalogues, indicating closely associated or fragmented objects, and mark objects which have FLAG totals <8 to use as a discriminant in later visual inspection (see Section 3.2.7). We generate square postage stamps of selected sources. These are centred at the galaxy centroids and have a windowed radial extent set to 4.5 times the calculated half-light radii, chosen to be large enough to ensure that the flexion-susceptible galaxy wings and background zero-constraints are included. Alongside the image stamp, we generate an associated 1σ noise stamp that includes both background sky noise and Poisson noise. 3.2 Flexion measurement – Analytic Image Modeling While research into the potential power of flexion measurement continues to remain popular (Er, Tereno & Mao 2012; Viola, Melchior & Bartelmann 2012; Cain, Bradač & Levinson 2016; Cardone et al. 2016; Rexroth, Natarajan & Kneib 2016), there have only been a handful of flexion analyses applied to real data to date (Leonard et al. 2007; Okura, Umetsu & Futamase 2008; Cain, Schechter & Bautz 2011; Leonard, King & Goldberg 2011), with several using the same data set and the majority in the widely studied rich cluster Abell 1689. Flexion measurements of real, individual sources have been achieved through a few techniques, including by measuring a combination of various third (and higher) order moments of the light distribution (as originally suggested by Goldberg & Natarajan 2002; Okura et al. 2008), decomposing the projected galaxy shape on to a polar orthonormal ‘shapelets’ basis set and truncating the series at a particular threshold (Goldberg & Leonard 2007; Massey et al. 2007b), and exploring the local potential field through parametrized ray-tracing (Cain et al. 2011), known as Analytic Image Modeling (AIM). We utilize a modified version of the AIM technique to determine the localized flexion field values for individual sources. AIM is distinct in that instead of measuring derived quantities (such as weighted surface brightness moments), it fits the lensed galaxy objects using a parametric model. By comparing the observed image to the uncertainty-weighted model image, we are able to optimize the parameters over reasonable bounds and thus constrain the flexion fields (and other shape information.) 3.2.1 Model parametrization Galaxy profiles are largely fit to parametric models with a radial luminosity distribution, and for this work we implement an elliptical Sérsic intensity profile. The Sérsic profile is particularly useful in that it encompasses a range of different models, including ones already well established through galaxy-luminosity profiling, including exponential, Gaussian, and de Vaucoulers profiles. Our model takes the form   $$I(\theta ) = I_e \exp {\left\lbrace -b_n\left[\left(\frac{\theta }{\theta _s}\right)^{1/n}-1\right]\right\rbrace },$$ (8)where θ is the circularized radius, θs is the radius of the isophote containing half the total flux of the galaxy (the half-light radius), Ie is the brightness at this effective radius, and n is the Sérsic index, which controls the steepness of the profile. We also include a complex ellipticity such that   $$|e| \,=\, \frac{q^2 - 1}{q^2 + 1}\,,$$ (9)and where q represents the axis ratio. Finally, bn is empirically determined to set θe to the half-light isophotal radius, and we adopt an approximated functional form   $$b_n = 1.992n - 0.3271, \quad 0.5 < n < 8.0$$ (10)valid within the given Sérsic index range (Capaccioli 1989). 3.2.2 Flexion, shear, and intrinsic ellipticity coupling As shown in Viola et al. (2012), measurements of the spin-1 and spin-3 distortions in a lensed image cannot be straightforwardly related to the $${\mathcal {F}}{}$$ and $${\mathcal {G}}{}$$ flexions, in general. In the intermediate regime between very weak and SL, coupling terms between the shear and flexion produce non-negligible spin-1 and spin-3 terms that contribute to the distortions originating from flexion alone. A significant inherent ellipticity in a source will also produce a similar effect through flexion–ellipticity coupling. In particular, moment-based approaches can introduce large biases in flexion measurements when failing to take these terms in into account, especially when accompanied by a centroid shift that depends on the ellipticity and galaxy morphology itself. Forward model-fitting techniques such as AIM avoid some of these issues by including an ellipticity/shear optimization term that accounts for some of the parameter coupling, but inherent source ellipticity can still introduce bias during parameter optimization – in general there is no appropriate application of a spin-1, spin-2, and spin-3 distortion to a circular source that is able to reproduce the same shape as when these distortions are applied to a source with inherent ellipticity. Attempting to fit both inherent ellipticity and shear simultaneously to an observed galaxy will likely result in non-convergence, as first identified in Cain et al. (2011) using simulated data. We approach this issue by fixing the lensing shear γ to zero and allowing the intrinsic ellipticity to absorb the two degenerate parameters into one. Although it is possible to set the lensing shear to that of a predetermined model, we instead aim to obtain an unbiased estimate of the shear field alongside the targeted flexion fields. As noted, it is possible that bias may be introduced by neglecting the flexion-coupling terms from inherent source ellipticity. We investigate this issue by rerunning our analysis with the source lensing shear values fixed to model data. We utilize an independent shear map drawn from public HFF SL cluster reconstructions (Lotz et al. 2017) and interpolate shear priors for our flexion sources before fitting with AIM. The resulting convergence and aperture mass reconstructions are very similar to those produced with null prior flexion, with peak locations and enclosed mass deviations within measurement uncertainties, |δθpeak| < 2 arcsec, |δMenclosed| < 1013 M⊙. The signs of these deviations are not consistent across all peak measurements and do not clearly show any evidence of systematic trends towards a particular direction. From this, we confirm the trend found in Cain et al. (2011) that shear priors on average do not drastically influence measured flexion values when using the AIM measurement technique. Removing this parameter degeneracy, there are seven model parameters for the unlensed galaxy light profile, and four effective parameters for the lensing fields. While we can reasonably expect this model to work well for simple source galaxies, it is also important to explore typical parameter spaces for existing data and to consider any known limitations. 3.2.3 Typical model suitability One of the biggest problems for flexion measurement, in general, is the inherent inhomogeneity of galaxy structure, particularly when combined with the emphasis on a select few sources (again as a result of strong local mass sensitivity). A complete analysis of ‘inherent flexion’ in unlensed populations has yet to be fully investigated at this point, though previous analyses have done so in cluster lensed populations (Leonard et al. 2007; Okura et al. 2008; Cain et al. 2011). Most galaxies are not simply circular or elliptical, and many can be expected to contain structural irregularities such as arms, bars, or star-forming regions. Elliptical models tend to have a large range of Sérsic indices, whose values depend strongly on luminosity (Blanton & Moustakas 2009). Furthermore, ellipticals can also be separated into ‘boxy’ and ‘discy’ subtypes. For the most part, however, a Sérsic model can be used to model the overall structure of most galaxy types. For elliptical galaxies, the Sérsic index n generally reflects what is apparently a single component galaxy. Spiral galaxies can also be described by a Sérsic index, but in this case n reflects a balance between the disc and the bulge, two clearly distinct components. Morphology typically focuses on the separation of the disc from the bulge, usually treating the disc as an exponential profile and fitting the bulge to a Sérsic profile, with n dependent on the type of the central component. Along with lenticular galaxies, these form the basis of our unlensed galaxy reconstruction efforts. The range of possible parameters create a broad distribution of shapes; drawing a good statistical likelihood of inherent shape is difficult. In this work, we aim to partially investigate inherent flexion by performing a simultaneous flexion analyses on the offset parallel Abell 2744 field. 3.2.4 Distortions and PSF corrections Our images consist of many offset, rotated, and stacked individual exposures over noncontiguous time periods (with a total average integration time of over 104 000 s, or ∼29 h in the F814w band) increasing the object S/N by a factor roughly equivalent to the number of orbits, at the cost of possible geometric/elliptical distortions and a more complicated effective Point Spread Function (PSF) that is non-trivial to model. The 90-min HST orbital period introduces thermal breathing fluctuations, while the telescope focus deviates over periods of weeks, and our data spans a range of these cycles. Furthermore, the off-axis position of the ACS introduces a spatially varying geometric distortion across the WFC chips. Despite the complexities inherent in the exact shape of the ACS PSF, Schrabback et al. (2010) finds that ∼ 97 per cent of the total PSF ellipticity variation in random pointings can be attributed to a single modelled parameter – the telescope focus. While there are techniques that have been successfully applied to model the PSF of stacked ACS images (e.g. Bacon et al. 2003; Rhodes et al. 2007; Harvey et al. 2015; Jauzac et al. 2016b), these have mostly found use in shear measurement, with accuracy requirements at the ∼1 per cent level. As the induced flexion in a source from the PSF is   $${\mathcal {F}}_{\mathrm{induced}} \sim {\mathcal {F}}_{\mathrm{PSF}} \frac{a_{\mathrm{PSF}}^4}{a_{\mathrm{source}}^4 + a_{\mathrm{PSF}}^4}$$ (11)as derived in Leonard et al. (2007), for typical WFC values ($${\mathcal {F}}_{{\rm PSF}} \sim 10^{-3} \, {\rm arcsec^{-1}}, a_{{\rm PSF}} \sim 0.1125{\rm \,arcse}{\rm c}$$) there is not a significant flexion contribution from the PSF, provided the source is sufficiently large. As one of our pre-measurement cuts requires Source Extractor sources to have an FWHM greater than twice that of the ACS PSF (2 × 0.1125 arcsec), this minimizes the amount of any possible PSF-induced flexion. To test the flexion anisotropy across the ACS WFC chips, we used the tinytim2 software to simulate an ACS PSF under typical instrument parameters. By varying the input pixel coordinates across the two chips, we created a grid of spatially varying PSFs representing a single exposure, and were then able to measure the flexion signal using our AIM implementation. The resultant $${\mathcal {F}}$$-flexion vectors vary smoothly in magnitude and direction across the chips with a maximum magnitude of $$|{\mathcal {F}}| = 0.005$$, much smaller than the lensing flexion signal expected in Abell 2744, particularly in sources larger than the PSF. In addition, stacking the offset and rotated individual exposures into our integrated data image reprojects and averages the directional biases of the flexion anisotropies, further reducing possible measurement biases. In practice, as an approximation we use a simple Gaussian convolution with FWHM $$= 2\sqrt{2\mathrm{ln}\,2}\sigma = 0.1125{\rm \,arcsec}$$ to account for the general smoothing of a PSF. We expect the PSF to have relatively small effects on the overall shape and shape parameters of our selected galaxies, but as a self-consistency check we perform the best fit again without any applied PSF. If the fit shape or model parameters are markedly different we exclude that object from our analyses. 3.2.5 Implementation To determine the best-fitting model parametrization, we use a Levenberg–Marquardt least-squares minimization scheme, the core of which is implemented through a modified version of the python translation3 of the idl code mpfit. Our pipeline combines this minimization technique with galaxy models simulated through the GalFlex4 module, which ray-trace second-order weak gravitational lensing effects (flexion), along with shear and convergence, on simulated source images in real space. While a full exploration of the 11-dimensional parameter space is possible through combined downhill and sampling techniques (e.g. simulated annealing), for this work we take advantage of initial galaxy parameter estimates to aid global minimization. 3.2.6 Goodness of fit and parameter uncertainties Already noted in Section 3.1, the source-selection process is a critical part of flexion analysis. In a lensing reconstruction of Abell 1689 using ACS data, Leonard et al. (2007) combines SL, shear, and flexion analysis to confirm substructure in the form of a second mass peak apart from the central cluster. However, the central mass peak is not recovered in their flexion analysis despite having a large mean number density of flexion sources ($$\bar{n}_g \sim 75$$ arcmin−2). A later flexion-only analysis of A1689 by Cain et al. (2011) is able to infer mass structure consistent with previous measurements with a mean flexion source density of $$\bar{n}_g \sim 26$$ arcmin−2. Okura et al. (2008) finds that they can recover significant structure in the core of A1689 with only nine galaxies using Subaru data, through a meticulous source-selection procedure resulting in a final source density of $$\bar{n}_g \sim 7.75$$ arcmin−2. The lack of a significantly stronger signal in the Leonard et al. (2007) analysis despite the larger source density indicates the necessity of selecting which galaxy sources effectively ‘make the cut.’ Our final flexion source catalogues is constructed by excluding the galaxies that have either poor or unphysical best-fitting parameters. A commonly used metric for determining a ‘good’ fit is the reduced chi-squared statistic $$\chi ^2_{\mathrm{red}}$$, defined as the weighted sum of squared deviations (the χ2) divided by the number of independent degrees of freedom of the system. Although we use χ2 as the optimized statistic in AIM, we do not use the $$\chi ^2_{\mathrm{red}}$$ as the main selector for cutting our post-fit sources as in Cain et al. (2011) – highly correlated pixel data combined with a highly non-linear model render the actual meaning of the $$\chi ^2_{\mathrm{red}}$$ uncertain. We point the reader to Andrae, Schulze-Hartung & Melchior (2010) for further details. Nevertheless, this statistic is useful to broadly separate groups, as large values do indicate a likely poor fit. Fig. 1 shows an example of accepted, uncertain, and rejected galaxy source fits and their residuals, separated by $$\chi ^2_{\mathrm{red}}$$ for convenience. Figure 1. View largeDownload slide Best-fitting Sérsic galaxy models in the Abell 2744 HST F814w band. Each row shows a typical reduced χ2 value used in part to group objects by a goodness-of-fit statistic. As systematic uncertainties in inherent model suitability, best-fitting model error analyses, and model-parameter-error-size correlation are not yet well-understood, this statistic is not the defining accepted/rejected source criterion, but rather a weighted component of an extended source-selection pipeline. Figure 1. View largeDownload slide Best-fitting Sérsic galaxy models in the Abell 2744 HST F814w band. Each row shows a typical reduced χ2 value used in part to group objects by a goodness-of-fit statistic. As systematic uncertainties in inherent model suitability, best-fitting model error analyses, and model-parameter-error-size correlation are not yet well-understood, this statistic is not the defining accepted/rejected source criterion, but rather a weighted component of an extended source-selection pipeline. Instead of the reduced chi-squared statistic, we employ an impartial stamp inspection Graphical Universe Interface (GUI) as described in Section 3.2.7, in combination with the following parameter cuts: We reject fits where the reduced flexions {|Ψ1|, |Ψ3|} ≥ 1.0 arcsec−1. These large values indicate either a poor or unphysical fit, or that the source is in a region where WL assumptions are not valid. We also impose the constraint that $$\lbrace \sigma _{\Psi _1}, \sigma _{\Psi _3}\rbrace > 0.001$$ arcsec−1 to prevent overfitting the associated parameters {Ψ1, Ψ3}, as in Cain et al. (2011). We reject object fits with axis ratio q = 1.0 in the same vein, as a best fit which is exactly symmetric is likely to either be a star or overfitted. Finally, for practical convenience, we also eliminate objects whose best-fit χred2 is large (≳5), in most cases representing either an improperly deblended or foreground galaxy. 3.2.7 Interactive-limited stamp inspection While the Source Extractor utility turns the non-trivial procedure of source identification, deblending, and shape estimation into something more routine for most purposes, the various strategy and parameter combinations – particularly for crowded fields with diverse redshift populations – can still prove a need for artful navigation. Parameters must be selected so as to identify both the small and faint sources, as well as the larger, brighter, and closer ones; as described in Section 3.1, a common technique in WL studies is to perform a ‘hot-cold’ routine. However, deblending issues still can persist throughout the pipeline. These issues can entail multiple galaxies masquerading as one, or vice versa, especially with extended galaxies containing luminous star-forming regions. Therefore, an additional layer of inspection can benefit a signal dependent on an accurate measure of a single galaxy's shape. While performing this kind of source inspection in a shear analysis would be time-prohibitive, this is not necessarily true for flexion analysis. As shown in Viola et al. (2012), accurate flexion measurements are particularly dependent on high S/N sources, and so by selecting objects according to suitable criteria the number of usable galaxies for flexion measurement is narrowed down to a much smaller and more feasible subset. Our stamp/fit inspection routine is implemented through a GUI as part of our flexion pipeline. Built with the python bindings for Qt (PyQt4), we compile the base code and all relevant modules into an executable to create a stand-alone application. Seen in Fig. 2, the GUI can be used to quickly and efficiently scroll through input stamp data sets for acceptance/rejection. We show a slightly extended region around the selected stamp in order to inspect its immediate neighbours for blending issues or galaxy fragmentation. We purposefully restrict the field of view in an attempt to minimize any observational bias. We also inspect the best-fitting parameters to check for problematic behaviours, including having values on the established parameter bounds or exactly null, or having particular combinations of values that clearly indicate a poor fit (e.g. a best-fit half-light-radius larger than the stamp itself). Figure 2. View largeDownload slide Inspecting our image and best-fitting stamp and fit parameters is critical for an analyses that can depend on just a few handful of strong signal images. We use a custom GUI to streamline the process of accepting/rejecting/flagging these sources for use in the mass reconstruction. Figure 2. View largeDownload slide Inspecting our image and best-fitting stamp and fit parameters is critical for an analyses that can depend on just a few handful of strong signal images. We use a custom GUI to streamline the process of accepting/rejecting/flagging these sources for use in the mass reconstruction. 3.3 Mass estimators A major component of flexion analysis is the non-trivial process of incorporating flexion measurements into an overall mass reconstruction. Ideally, both shear and flexion measurements could be used to constrain the mass profile more than either could alone. As we've noted, the shear probes large scales more effectively, while flexion is particularly suited to investigate smaller-scale structure. However, the flexion field alone can give valuable information about both small-scale structure in a lensing cluster and the overall mass profile. To do so, we employ two independent techniques to get a measure of a mass structure based solely on positional flexion measurements: a direct calculation of the convergence field, stemming directly from the fact that flexion is the gradient of the convergence; and the construction of a mass-signal map through the aperture mass technique (Schneider 1996), which relates the weighted integral of flexion measurements within an aperture to the respective weighted integral of the convergence within that aperture (Cain et al. 2011; Leonard et al. 2011). Each technique offers a different interpretation of the flexion measurements, and the utilization of both can provide a broader image of any substructure present within the galaxy cluster. 3.3.1 Direct flexion-convergence convolution Using the various directional derivative relations between lensing distortions and a gravitational potential, we are able to relate (up to the mass-sheet degeneracy) the convergence field κ with the flexion fields $${\mathcal {F}}, {\mathcal {G}}$$ through the convolution:   $$\kappa (\boldsymbol {\theta }) = \int _{\mathbb {R}^2} d^2 \theta ^{\prime } \, \mathcal {D}_M (\boldsymbol {\theta } - \boldsymbol {\theta }^{\prime }) \, {\mathcal {F}}(\boldsymbol {\theta }^{\prime }),$$ (12)where $$\mathcal {M}$$ is the $${\mathcal {F}}$$ or $${\mathcal {G}}$$ flexion field in the lens plane, and the respective convolution kernels are given by   \begin{eqnarray} \mathcal {D}_{{\mathcal {F}}} = \frac{\theta _1 + i\theta _2}{2\pi |\theta |^2}, \end{eqnarray} (13)  \begin{eqnarray} \mathcal {D}_{{\mathcal {G}}} = \frac{(\theta _1 + i\theta _2)^3}{2\pi |\theta |^4}. \end{eqnarray} (14)as outlined in Leonard et al. (2011). To avoid the effects of infinite sampling noise (Bartelmann & Schneider 2001), the irregularly spaced flexion data is first smoothed with an uncertainty-weighted Gaussian kernel before convolving the resulting finely spaced flexion grid with the appropriate flexion convolution kernel. The resulting convergence map is a soft function of the selected smoothing scale, which generally reflects the scale of structure that is emphasized. To obtain estimates of the uncertainty in our convergence reconstructions, we investigate the corresponding uncertainties in both the values of the measured flexions (or shear) and the shot noise as a result of the scattered spatial positions of selected flexion (shear) sources. We simulate the prior probability distributions for these observables and ultimately apply Monte Carlo resampling methods to determine the variance in the derived convergence map as well as observational derivatives such as our estimates of enclosed peak masses. In order to sample the probability distribution, we make an assumption that our flexion measurement uncertainties {$$\sigma _{{\mathcal {F}}{}_i}, \sigma _{{\mathcal {G}}{}_i}$$} are approximately normally distributed and follow the relation:   \begin{eqnarray} \sigma _{{\mathcal {F}}{}_i ({\mathcal {G}}{}_i)} = \sqrt{ \left(\frac{\sigma _{\mathcal {F}\cdot a}}{a}\right)^2 + (\sigma _{\mathcal {F}_i, {\rm AIM}})^2}, \end{eqnarray} (15)where the term $$\sigma _{\mathcal {F}\cdot a}$$ (or $$\sigma _{\mathcal {G}\cdot a}$$) indicates the 1σ variation in the measured values of intrinsic flexion inherent in the unlensed galaxy population, added in quadrature to $$\sigma _{\mathcal {F}_i, {\rm AIM}}$$ (or $$\sigma _{\mathcal {G}_i, {\rm AIM}}$$), the individual flexion measurement uncertainty derived from the best-fitting model's covariance matrix. For clarity, we note that the population dimensionless flexion dispersions $$\sigma _{\mathcal {F}\cdot a}$$ and $$\sigma _{\mathcal {G}\cdot a}$$ are calculated from the combined sets of real and imaginary components {$$\mathcal {F}_1 \cdot a\,, \, \mathcal {F}_2 \cdot a$$} , {$$\mathcal {G}_1 \cdot a\,, \mathcal {G}_2 \cdot a$$} of spin-1 and spin-3 flexions, respectively, which are assumed to be normally distributed with intrinsic mean value of zero. For each iteration in a large set of successive convergence reconstructions, the resampled flexion values are drawn from a normal distribution with mean and standard deviation equal to the measured flexion values and corresponding 1σ uncertainties (including both model-based and population uncertainty). By stacking the value-resampled reconstructions, we are able to construct a two-dimensional variance map along the stacked axis. Our peak mass uncertainties are obtained from the variance of the enclosed aperture mass as measured in each of these frames. While this convolution relation offers a direct linear inversion map with a straightforward interpretation, it suffers from the same difficulties as the shear convolution, namely that: The integral in equation (12) extends over $$\mathbb {R}^2$$, while usable flexion data are constrained to a few square arcminutes. Smoothing over the lensing fields is necessary to avoid infinite noise from discrete sampling, which is a particular concern as the strength of the flexion signal falls off much more quickly than shear further from mass structure. Source-position shot noise can significantly affect local field values depending on smoothing scale. A low source density and large variance of nearest-neighbour separation distance can strongly influence the relative importance of a single source in any particular smoothing aperture. Using multiple scales can mitigate potential local field extremes by allowing narrow-band peak outliers to be identified and removed. The cluster regions where flexion signal is expected to be strong are usually those where the convergence cannot be assumed to be small, and thus the observables are better approximated by the reduced flexion Ψ1 and Ψ3. To offer an alternative and complementary method to the direct convergence inversion, we also utilize our implementation of the aperture mass technique. 3.3.2 Aperture mass statistic The aperture mass statistic, $$\mathcal {F}_{{\rm aper}}$$, is able to relate measured flexion to the underlying lens convergence through more general weighting kernels than the direct flexion-convergence convolution. Apertures of characteristic radius R are laid down in a grid-like pattern over the field. Within each aperture the convergence is convolved with a filter function w(r). This convolution is then related to the measured flexion convolved with an appropriate filter function Q(r). The role of the aperture mass technique is to evaluate mass structure detections. By randomly rotating the flexion vectors in the image and running the reconstruction multiple times, it is possible to create a S/N map of detected mass signal, analogous to the flux S/N ratio of a galaxy, for example. The aperture statistic is defined by relating the spin-1 flexion and convergence as follows:   \begin{eqnarray} \mathcal {F}_{{\rm aper}}{}(\theta _0;R) = \int _{|\theta |\le R} {\rm d}^2\, \theta \, \kappa (\theta +\theta _0)w(|\theta |), \end{eqnarray} (16)  \begin{eqnarray} \mathcal {F}_{{\rm aper}}{}(\theta _0;R) = \int _{|\theta |\le R} {\rm d}^2\, \theta \, \mathcal {F}_E(\theta ;\theta _0)Q(|\theta |), \end{eqnarray} (17)where $$\mathcal {F}_E$$ is the component of the first flexion oriented towards the centre of the aperture. The weight and filter functions are free to be selected, while obeying the relations   $$Q(r) = -\frac{1}{r} \int _0^r w(r^{\prime }) r^{\prime } {\rm d}r^{\prime } \quad w(r) = -\frac{1}{r} Q(r) - \frac{{\rm d}Q}{{\rm d}r},$$ (18)where the variable r represents |θ − θ0|. The weight function is constrained to go to zero smoothly at the aperture boundary and have a mean value of zero,   $$\int _0^R w(r^{\prime }) r^{\prime } {\rm d}r^{\prime } = 0.$$ (19)In this work, we use a family of polynomial weighting functions as described and used in the literature (Schneider et al. 1998; Cain et al. 2011; Leonard et al. 2011). These filter functions are chosen to test the persistence of mass estimates across a variety of scales. The convergence weighting function is defined as   $$w(r) = A_l \frac{(2 + l)^2}{\pi } \left(1 - \frac{r^2}{R^2}\right)^l \left(\frac{1}{2+l} - \frac{r^2}{R^2}\right)$$ (20)along with the flexion kernel   $$Q(r) = -A_l \frac{2 + l}{2\pi } r \left(1 - \frac{r^2}{R^2}\right)^{1+l}$$ (21)and normalization   $$A_l = \frac{4}{\sqrt{\pi }}\frac{\Gamma (\frac{7}{2} + l)}{\Gamma (3 + l)}$$ (22)as a function of the polynomial index, ℓ, and the characteristic aperture radius, R. A larger polynomial index ℓ makes the kernel more sensitive to smaller scale structure within the aperture, at the expense of having larger noise fluctuations. Lower indices smooth over a larger area, reducing noise at the cost of resolution. The aperture radius R has a similar effect: as the radius increases the flexion signal associated with the centre of the aperture falls off quickly, smoothing over smaller scale structures. To implement this method, we turn the flexion convolution into a discrete sum of kernel-weighted flexion sources at each point θ0 over a regular coordinate grid. We include flexion uncertainty in the kernel Q(r) by introducing a normalized weighting term qn applied to each flexion source:   $$Q(r) \rightarrow Q_i(r) \, \,q_i ; \quad q_i = \frac{1}{\sigma _{{\mathcal {F}}{}_i}^2},$$ (23)where the flexion uncertainty, $$\sigma _{{\mathcal {F}}{}_i}$$, refers to the total combined uncertainty from both the AIM measurement and the inherent population dispersion (equation 15). In practical application, the methods used to implement both the flexion-convergence inversion and the aperture mass statistic techniques are quite similar, though the noise and smoothing properties as well as the physical meanings of the reconstructions are different. 3.3.3 Flexion profile fitting While the two techniques as described above are able to give more than sufficient information for high-quality mass reconstructions, they are both subject to similar types of restrictions and assumptions. As a third perspective we directly fit gravitational lens models to our flexion data to further predict and constrain the position, shape, and size of cluster substructure. Parametric lens modelling has the advantage of using observed data more directly than by using a smoothed data set along with simultaneously being capable of breaking the mass-sheet degeneracy by constraining the slope of the shear/flexion data. We fit a SIS profile to each of the identified peak locations and allow the centre coordinates to moderately vary. The best-fitting values of the parameter set {θ1, 0, θ2, 0, θE} are optimized through a comparison between our flexion data and the model flexion,   $${\mathcal {F}}{} = -\frac{\theta _E}{2(\theta -\theta _0)^2} \, e^{i\phi }\qquad {\mathcal {G}}{}(\theta ) = \frac{2\theta _E}{3(\theta -\theta _0)^2} \, e^{3i\phi },$$ (24)where the preceding constant θE is the Einstein radius of the lens, which is related to the lensing depth and velocity dispersion, $$\sigma _{\nu _{{\rm 1D}}}$$, as   $$\theta _E = 4\pi \left(\frac{\sigma _{\nu _{{\rm 1D}}}}{c}\right)^2 \left(\frac{D_{{\rm ls}}}{D_{{\rm s}}}\right).$$ (25) We utilize the non-reduced formulations of the flexion in this model fitting – this avoids the direction-reversing behaviour that the SIS model tends to predict for the measured flexion within appreciable distances of the origin. The source galaxies that would be affected by this behaviour would exhibit radial shear with highly distorted profiles according to this model, and would be selected out of our source sample before flexion measurement. A result of this approximation is that the best-fitting model estimates for the Einstein radii and velocity dispersions may be somewhat larger than if a more accurate flexion profile and mass-degeneracy treatment were utilized. While fitting to more complicated profiles is possible, the SIS lens is commonly referenced in both past and current lensing work as a simple and analytically tractable model with well-understood model parameters. As a complement to our non-parametric reconstruction maps, the SIS flexion model should prove sufficient for the purposes of this work. 4 RESULTS AND SUMMARY 4.1 Source selection results Our extraction and pre-measurement selection routine returned 1344 out of the 4249 detected galaxies in the F814w-observed 3.5 arcmin × 3.5 arcmin field of Abell 2744. The after-measurement selection cuts give us our final catalogue of 969 sources, corresponding to 82.1 sources arcmin−2 or about 23 per cent of the initially detected galaxy population. In the cluster's associated parallel field, our pre-selection routine returned 1001 out of 4142 detected galaxies. Post-measurement cuts reduce our final catalogues to 923 usable galaxies, corresponding to 78.2 sources arcmin−2. As the observational details of the parallel field are identical to those of the cluster (e.g. a field of the exact same size, with the same total integration time in each band, in the same number orbits), a source density comparable to the cluster field is in line with the expectation that our cluster selection routine is able to exclude contaminating cluster members. We measure the spin-1 and spin-3 dimensionless flexion dispersions in our cluster source population to be $$\sigma _{\mathcal {F}\cdot a} = 0.09$$ and $$\sigma _{\mathcal {G}\cdot a} = 0.12$$, with similar quantities calculated for the selected parallel field population. The mean dimensionless flexions $$\langle {\mathcal {F}}{}\cdot a\rangle , \langle {\mathcal {G}}{}\cdot a\rangle$$ were found to be consistent with zero in both the cluster and parallel populations as well. Previous flexion analyses (Leonard et al. 2007; Okura et al. 2008; Cain et al. 2011) have found population dimensionless flexion dispersions of $$\sigma _{\mathcal {F}\cdot a}\simeq (0.03 {\rm -}0.05)$$ and $$\sigma _{\mathcal {G}\cdot a} \simeq (0.04 {\rm -}0.06)$$, roughly a factor of 2–3 times smaller than our calculated values. It is unclear if this population difference is simply due to chance or if possible systematic differences exist between the different flexion procedures. As all other previous flexion analyses have focused on the galaxy cluster Abell 1689, it is possible that differences in observational details (redshifts, integration time, etc.) or properties of the selected source populations may account for some or all of the variance in measured flexion. In addition to our flexion statistics, by deriving shear estimators from measured ellipticity $$\widetilde{\gamma } \approx e$$, we also note the AIM-associated shear population statistics: $$\sigma _{\gamma _1} = 0.30, \sigma _{\gamma _2} = 0.32,$$ 〈γ1〉 = 0.0037, 〈γ2〉 = −0.0361. These shear dispersions are typical of WL studies for the most part, although there is a non-negligible mean shear component γ2. As the high-resolution HST data set used in this work extends only a few arcminutes from the cluster centre, this systematic mean shear bias is much less than the expected lensing shear signal (γ ∼ 0.1) at fields bounds. However, most wide field shear analyses commonly attempt to measure cluster lensing signal at extended distances much farther than those explored here, and a similar mean shear residual would need to be accounted for in order to observe an actual lensing signal in the range of a few percentages. We approximate our effective background source redshifts by averaging the lensing depths determined in two previous lensing studies using the same HFF images. The source sample used in Medezinski et al. (2016) was determined to have 〈zeff〉 = 1.27 ± 0.08 with a mean lensing depth of Dls/Ds = 0.68 ± 0.11 while the sample used in Jauzac et al. (2016b) had an effective source redshift of 〈zeff〉 = 1.586 with a mean lensing depth of Dls/Ds = 0.729. Taking a simple average yields an effective critical density of Σcrit = 3.6 × 1010 h−1M⊙ arcsec−2. 4.2 Mass reconstruction results Following the non-parametric methods as described in Section 3.3, we recreate the mass distribution of Abell 2744 over a 5 arcmin × 5 arcmin square plane on which our 3.5 arcmin × 3.5 arcmin data field is diagonally circumscribed. We use an equirectangular projection of the sky plane as a planar coordinate system, oriented such that declination, δ, increases along the positive y-axis to the north; and right ascension, α, decreases along the positive x-axis towards the west. The origin of our relative sky coordinate system is set to the equatorial coordinates: {α = 3.5890213°, δ = −30.397196°} (deg, J2000.0), as the centre position in our stacked HST pixel data image. The flexion measurements from the lensed images can be visualized in a number of ways. In Fig. 3, we plot a cluster-scale reconstruction of the convergence map of Abell 2744 as well as the parallel field noise reconstruction, using $${\mathcal {F}}$$-flexion data exclusively. At a smoothing scale of 15 arcesc, only two of three substructures are detected – a large central peak that follows the elliptical luminosity structure of the central cluster galaxies, declining smoothly towards the western edges, and a second peak centred just south of the bright spider-like foreground galaxy to the east, offset 1.43 arcmin (273 h−1 kpc) from the cluster core. Figure 3. View largeDownload slide Independent $${\mathcal {F}}$$-flexion data are used to create cluster-scale mass field reconstructions of Abell 2744 and noise in its associated parallel field. The top left figure shows the cluster convergence (up to a mass-sheet degeneracy scaling) using a 15 arcsec smoothing kernal, while the top right shows a parallel field noise reconstruction under the same parameters. The bottom left and right figures are the flexion aperture mass S/N reconstructions for the cluster and parallel fields, respectively, at aperture radius 60 arcsec and index ℓ = 5. The axes are in units of arcseconds. Figure 3. View largeDownload slide Independent $${\mathcal {F}}$$-flexion data are used to create cluster-scale mass field reconstructions of Abell 2744 and noise in its associated parallel field. The top left figure shows the cluster convergence (up to a mass-sheet degeneracy scaling) using a 15 arcsec smoothing kernal, while the top right shows a parallel field noise reconstruction under the same parameters. The bottom left and right figures are the flexion aperture mass S/N reconstructions for the cluster and parallel fields, respectively, at aperture radius 60 arcsec and index ℓ = 5. The axes are in units of arcseconds. While the 15 arcsec smoothing scale presented in Fig. 3 only hints at an offshoot of increased convergence extending from the central region towards the north-east, and a second towards the south-east, finer smoothing scales in addition to our mass aperture reconstructions reveal apparent mass peaks within both tendrils, with a more substantial north-eastern peak and a lesser south-eastern detection. Using a 60 arcsec smoothing kernel, on the other hand, reveals a third structure, towards the north-east. We can get a better sense of the detailed shape of the substructures – and in particular the central substructure – by using a complementary map, the aperture mass reconstruction, in Fig. 4. Finally, in Fig. 5 presents the results of our shear-only convergence reconstructions, with noise maps and two different scaling kernels. We describe more detailed interpretations of these maps below. Figure 4. View largeDownload slide Magnified section of flexion aperture mass reconstructions centred around the detected core substructure. The aperture characteristic truncation scale $$\mathcal {R} = \lbrace 30{\rm },45,60,75\, {\rm arcsec}\rbrace$$ varies along the vertical axis while the flexion weighting function's polynomial index ℓ = {3, 5, 7} spans the horizontal. The plots are arranged such that the effective resolution of the aperture maps is the broadest towards the lower left corner of the figure, and the finest towards the top right corner, with approximate resolution scaling $$\Delta \sim \mathcal {R}_{{\rm aper}} / \ell$$. At higher effective resolutions, the detected core substructure appears to begin to deblend into multiple components – the most significant peak focuses towards the region associated with the southern BCG pair, while mass signature begins to show both slightly north-west of the BCG pair and also around the slightly easterly spiral galaxy. Our peak-moment analysis at $$\mathcal {R} = 60{\rm \,arcsec}, \ell = 5$$ shows a 13 arcsec offset between the maximum value peak and the centroid, indicating either a very skewed or dimodal structure. Figure 4. View largeDownload slide Magnified section of flexion aperture mass reconstructions centred around the detected core substructure. The aperture characteristic truncation scale $$\mathcal {R} = \lbrace 30{\rm },45,60,75\, {\rm arcsec}\rbrace$$ varies along the vertical axis while the flexion weighting function's polynomial index ℓ = {3, 5, 7} spans the horizontal. The plots are arranged such that the effective resolution of the aperture maps is the broadest towards the lower left corner of the figure, and the finest towards the top right corner, with approximate resolution scaling $$\Delta \sim \mathcal {R}_{{\rm aper}} / \ell$$. At higher effective resolutions, the detected core substructure appears to begin to deblend into multiple components – the most significant peak focuses towards the region associated with the southern BCG pair, while mass signature begins to show both slightly north-west of the BCG pair and also around the slightly easterly spiral galaxy. Our peak-moment analysis at $$\mathcal {R} = 60{\rm \,arcsec}, \ell = 5$$ shows a 13 arcsec offset between the maximum value peak and the centroid, indicating either a very skewed or dimodal structure. Figure 5. View largeDownload slide Independent AIM-derived shear data are used in a cluster-scale mass field reconstruction of Abell 2744 and a noise reconstruction in its associated parallel field (both up to a mass-sheet degeneracy scaling). The top left and right figures show the convergence field resulting from a 15 arcsec smoothing kernal applied to the cluster and parallel fields, respectively, while the bottom left and right figures display the results while using a 30 arcsec kernal instead. The axes are in units of arcseconds. Figure 5. View largeDownload slide Independent AIM-derived shear data are used in a cluster-scale mass field reconstruction of Abell 2744 and a noise reconstruction in its associated parallel field (both up to a mass-sheet degeneracy scaling). The top left and right figures show the convergence field resulting from a 15 arcsec smoothing kernal applied to the cluster and parallel fields, respectively, while the bottom left and right figures display the results while using a 30 arcsec kernal instead. The axes are in units of arcseconds. 4.2.1 Convergence-flexion inversion reconstructions There is a degree to which the smoothing scale for a mass reconstruction is somewhat arbitrary, though we typically set it so that the noise does not dominate. Further, when turning a reconstruction into concrete mass estimates, there is bound to be some ambiguity. For instance, the mass sheet degeneracy allows us to adjust the cluster up and down, and is dependent on the boundary conditions(e.g. Bradač, Lombardi & Schneider 2004). However, the effect on our measured peak mass within an aperture is limited – for example, increasing or reducing the assumed boundary convergence by κ ± 0.1 results in a contained mass difference of ∓0.04 × 1014 h−1 M⊙ for our primary structure, or 2 per cent of our calculated mass, with similar results (∼3 per cent) for the secondary and tertiary peaks. We focus on a smoothing scale of 15 arcsec (except as discussed below), and focus on comparisons to the literature within particular apertures. Our reconstructed peak positions and ellipticities, compared with those found previous researchers, can be found in Table 1, while the estimated masses within apertures (selected for comparison with previous researchers) can be found in Table 2. Table 1. Detected substructure geometries. The listed structural properties are derived from the moments of each peak's associated pixels. Included pixels are determined through nearest-neighbour proximity and have a minimum threshold of 2.5σ above the peak-subtracted background level. The detection geometries from three different reconstruction methods are compared: flexion aperture mass inversion, $$\mathcal {F}_{{\rm aper}}{}$$; direct flexion-convergence inversion, $$\kappa _{{\mathcal {F}}{}}$$; and optimization of an SIS model to measured flexion data. Coordinates are specified in arcseconds and relative to the centre of our HST pixel image, at α = 3.589 0213, δ = −30.397 196 (deg, J2000.0). Method/  Xpeak  Ypeak  a  b  e  ϕ  Peak  structure  (arcsec)  (arcsec)  (arcsec)  (arcsec)    (deg)  value  $$\quad {\mathcal {F}}{}_{{\rm aper}}$$            S/N  (1) Core  11.47 ± 13.3  −9.98 ± 17.1  31.36 ± 3.6  16.13 ± 1.3  0.58 ± 0.1  67.05 ± 8.5  3.50 ± 0.84  (2) NE  −55.15 ± 6.7  72.26 ± 14.0  18.03 ± 1.8  9.78 ± 0.6  0.55 ± 0.1  105.04 ± 6.3  2.68 ± 0.85  (3) E  −68.26 ± 9.9  −0.59 ± 10.0  17.16 ± 1.9  10.47 ± 0.8  0.46 ± 0.1  136.00 ± 8.0  2.34 ± 0.73  $$\quad \kappa _{{\mathcal {F}}{}}$$            $$\kappa \qquad$$  (1) Core  3.82 ± 4.0  −1.76 ± 4.4  30.1 ± 0.3  19.19 ± 0.1  0.42 ± 0.01  54.02 ± 0.6  2.50 ± 0.29  (2) NE  –  –  –  –  –  –  –  (3) E  −61.7 ± 2.3  −5.28 ± 2.4  26.0 ± 0.2  23.86 ± 0.2  0.08 ± 0.02  71.02 ± 2.9  2.65 ± 0.30  Lensmodel (SIS)            σν (km s−1)  (1) Core  3.14 ± 4.3  −11.97 ± 6.5  –  –  –  –  1904 ± 346  (2a) NE  −47.16 ± 2.2  41.41 ± 3.4  –  –  –  –  747 ± 380  (2b) NE  −35.99 ± 2.7  68.83 ± 4.2  –  –  –  –  594 ± 690  (3a) E  −52.79 ± 2.5  −21.78 ± 2.2  –  –  –  –  1161 ± 293  (3b) E  −65.44 ± 3.2  0.21 ± 6.2  –  –  –  –  881 ± 552   κγ            $$\kappa \qquad$$  (1) Core  3.82 ± 1.99  −8.81 ± 4.57  50.22 ± 2.58  37.51 ± 1.29  0.28 ± 0.06  0.50 ± 81.6  0.32 ± 0.03  (1) NE  −40.95 ± 3.33  55.81 ± 3.91  6.71 ± 0.24  3.06 ± 0.11  0.65 ± 0.03  128.74 ± 2.20  0.036 ± 0.024  Method/  Xpeak  Ypeak  a  b  e  ϕ  Peak  structure  (arcsec)  (arcsec)  (arcsec)  (arcsec)    (deg)  value  $$\quad {\mathcal {F}}{}_{{\rm aper}}$$            S/N  (1) Core  11.47 ± 13.3  −9.98 ± 17.1  31.36 ± 3.6  16.13 ± 1.3  0.58 ± 0.1  67.05 ± 8.5  3.50 ± 0.84  (2) NE  −55.15 ± 6.7  72.26 ± 14.0  18.03 ± 1.8  9.78 ± 0.6  0.55 ± 0.1  105.04 ± 6.3  2.68 ± 0.85  (3) E  −68.26 ± 9.9  −0.59 ± 10.0  17.16 ± 1.9  10.47 ± 0.8  0.46 ± 0.1  136.00 ± 8.0  2.34 ± 0.73  $$\quad \kappa _{{\mathcal {F}}{}}$$            $$\kappa \qquad$$  (1) Core  3.82 ± 4.0  −1.76 ± 4.4  30.1 ± 0.3  19.19 ± 0.1  0.42 ± 0.01  54.02 ± 0.6  2.50 ± 0.29  (2) NE  –  –  –  –  –  –  –  (3) E  −61.7 ± 2.3  −5.28 ± 2.4  26.0 ± 0.2  23.86 ± 0.2  0.08 ± 0.02  71.02 ± 2.9  2.65 ± 0.30  Lensmodel (SIS)            σν (km s−1)  (1) Core  3.14 ± 4.3  −11.97 ± 6.5  –  –  –  –  1904 ± 346  (2a) NE  −47.16 ± 2.2  41.41 ± 3.4  –  –  –  –  747 ± 380  (2b) NE  −35.99 ± 2.7  68.83 ± 4.2  –  –  –  –  594 ± 690  (3a) E  −52.79 ± 2.5  −21.78 ± 2.2  –  –  –  –  1161 ± 293  (3b) E  −65.44 ± 3.2  0.21 ± 6.2  –  –  –  –  881 ± 552   κγ            $$\kappa \qquad$$  (1) Core  3.82 ± 1.99  −8.81 ± 4.57  50.22 ± 2.58  37.51 ± 1.29  0.28 ± 0.06  0.50 ± 81.6  0.32 ± 0.03  (1) NE  −40.95 ± 3.33  55.81 ± 3.91  6.71 ± 0.24  3.06 ± 0.11  0.65 ± 0.03  128.74 ± 2.20  0.036 ± 0.024  View Large Table 2. Sky position and enclosed mass for substructure identified in this work, Merten et al. (2011), Medezinski et al. (2016), and Jauzac et al. (2016b). The coordinates Xpeak and Ypeak are specified in arcseconds and are relative to the centre of our HST pixel image, at α = 3.589 0213, δ = −30.397 196 (deg, J2000.0). These refer to the positions of either the peak maximum value (as in this work), or the peak centroid (as in Merten et al. 2011), as reported in each respective analysis. In order to directly compare our mass and size estimates with those in the relevant literature we assume here the dimensionless Hubble parameter h = 0.7. Uncertainties in coordinates and enclosed masses for the different techniques used in this work, $$\mathcal {F}_{{\rm aper}}{}, \kappa _{\mathcal {F}}{}, \kappa _{\gamma }, {\rm SIS}_{(1)}, {\rm SIS}_{(2)}$$, are derived from the variance of peak-associated pixel values throughout independent sets of bootstrap realizations as discussed in the text. Errors are derived from the standard deviation of parameter values as calculated for the selected pixels applied to each bootstrap realization.   Source  Xpeak (arcsec)  Ypeak (arcsec)  Renclosing (kpc)  M( < Renclosing) (1013 M⊙)  Core  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  11.47 ± 13.3  −9.98 ± 17.1   n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  3.82 ± 4.0  −1.76 ± 4.4  150  16.2 ± 1.2    κγ (this work)  3.81 ± 2.0  −8.81 ± 4.6  150 (250)  6.27 ± 0.35 (13.3 ± 0.8)    SIS (this work)  3.14 ± 4.3  −11.97 ± 6.5  150  19.9 ± 7.2    Jauzac et al. (2016b)  8.58  −10.72  150  13.55 ± 0.09    Medezinski et al. (2016)  12.91 ± 10.8  −15.84 ± 9.6  200c  77 ± 34    Merten et al. (2011)  $$-2.24^{+11}_{-8}$$  $$-16.71^{+6}_{-13}$$  250  22.4 ± 5.5  E  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −68.26 ± 9.9  −0.59 ± 10.0  $$\varnothing$$  $$\varnothing$$    $$\kappa _{\mathcal {F}}{}$$ (this work)  −61.7 ± 2.3  −5.28 ± 2.4  110  8.81 ± 0.52    SIS(1) (this work)  −52.79 ± 2.5  −21.78 ± 2.2  150  7.4 ± 3.7    SIS(2) (this work)  −65.44 ± 3.2  0.21 ± 6.2  150  4.3 ± 5.3    Jauzac et al. (2016b)  ..  ..  ..  ..    Medezinski et al. (2016)  ..  ..  ..  ..    Merten et al. (2011)  ..  ..  ..  ..  NE  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −55.15 ± 6.7  72.26 ± 14.0  n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  −  −  68  2.92 ± 0.26b    SIS(1) (this work)  −47.16 ± 2.2  41.41 ± 3.4  150  3.1 ± 3.1    SIS(2) (this work)  −35.99 ± 2.7  68.83 ± 4.2  150  1.9 ± 4.5    Jauzac et al. (2016b) a   −46.91  81.16  150  5.00 ± 0.40    Medezinski et al. (2016)  −52.02 ± 19.8  78.25 ± 31.8  200c  28 ± 16    Merten et al. (2011)  ×  ×  ×  ×    Source  Xpeak (arcsec)  Ypeak (arcsec)  Renclosing (kpc)  M( < Renclosing) (1013 M⊙)  Core  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  11.47 ± 13.3  −9.98 ± 17.1   n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  3.82 ± 4.0  −1.76 ± 4.4  150  16.2 ± 1.2    κγ (this work)  3.81 ± 2.0  −8.81 ± 4.6  150 (250)  6.27 ± 0.35 (13.3 ± 0.8)    SIS (this work)  3.14 ± 4.3  −11.97 ± 6.5  150  19.9 ± 7.2    Jauzac et al. (2016b)  8.58  −10.72  150  13.55 ± 0.09    Medezinski et al. (2016)  12.91 ± 10.8  −15.84 ± 9.6  200c  77 ± 34    Merten et al. (2011)  $$-2.24^{+11}_{-8}$$  $$-16.71^{+6}_{-13}$$  250  22.4 ± 5.5  E  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −68.26 ± 9.9  −0.59 ± 10.0  $$\varnothing$$  $$\varnothing$$    $$\kappa _{\mathcal {F}}{}$$ (this work)  −61.7 ± 2.3  −5.28 ± 2.4  110  8.81 ± 0.52    SIS(1) (this work)  −52.79 ± 2.5  −21.78 ± 2.2  150  7.4 ± 3.7    SIS(2) (this work)  −65.44 ± 3.2  0.21 ± 6.2  150  4.3 ± 5.3    Jauzac et al. (2016b)  ..  ..  ..  ..    Medezinski et al. (2016)  ..  ..  ..  ..    Merten et al. (2011)  ..  ..  ..  ..  NE  $$\mathcal {F}_{{\rm aper}}{}$$ (this work)  −55.15 ± 6.7  72.26 ± 14.0  n/a  n/a    $$\kappa _{\mathcal {F}}{}$$ (this work)  −  −  68  2.92 ± 0.26b    SIS(1) (this work)  −47.16 ± 2.2  41.41 ± 3.4  150  3.1 ± 3.1    SIS(2) (this work)  −35.99 ± 2.7  68.83 ± 4.2  150  1.9 ± 4.5    Jauzac et al. (2016b) a   −46.91  81.16  150  5.00 ± 0.40    Medezinski et al. (2016)  −52.02 ± 19.8  78.25 ± 31.8  200c  28 ± 16    Merten et al. (2011)  ×  ×  ×  ×  Notes.aThe NE structure is consistent with the Jauzac et al. (2016b) S1 substructure. bTo directly compare our derived masses and distances in physical units (M⊙ and kpc), we temporarily adopt the shared common cosmological model with Hubble constant h = 0.7. View Large As the NE peak is not clearly detected at this scale, the NE aperture is set to the centroid of the $$\mathcal {R} = 60{\rm \,arcsec}$$$$\mathcal {F}_{{\rm aper}}{}$$ reconstruction. The aperture sizes were selected such that the they approximately encircled the identified peaks down to a value around the calculated 1σ background level for both the $$\mathcal {F}_{{\rm aper}}{}$$ and κ reconstructions, with consideration for a representative size that could be easily standardized or replicated and that worked reasonably well for each reconstruction type. 4.2.2 Aperture mass reconstructions While the convergence–flexion relation allows us to calculate physical properties of the cluster or peaks inside a given aperture (provided the convergence is normalized to the mass-sheet degeneracy), the aperture mass method gives a more qualitative result of where mass peaks lie. The lower left and lower right graphs in Fig. 3 show the cluster/parallel aperture mass/noise reconstructions using aperture R = 60 arcsec and polynomial index ℓ = 5. Mass S/N contours in the cluster indicate a strong central peak aligned and centred on the Brightest Cluster Galaxy (BCG) and extending along the cluster's luminosity semimajor axis. A corresponding lesser peak is identified to the east along with another to the north-east. The signal tapers to zero near the edge of the data, with the exception of a few structures along the edge. It is important to note, however, that it is not straightforward to compare the significance values of these peak detections between works directly. Our significance values are derived through different means (our aperture mass method) than those of our comparison papers, and therefore the statistical meaning of these sets of values is not necessarily the same. While direct summation techniques (e.g. Fourier estimates of κ) are essential in providing mass estimates of individual structures, we stress the importance of the aperture mass statistic in this work because of its usefulness in identifying individual structures. Some care must be taken with the choice of smoothing radius and order of the polynomial index, however, in that they yield an effective smoothing scale of $$\sim \mathcal {R}_{{\rm aper}} / \ell$$. In Fig. 3, this gives an approximate smoothing scale of ∼12 arcsec. However, in Fig. 4, we see that the location and number of central peak(s) is(are) highly sensitive to these parameters. It is significant to note that a north-eastern local maxima is not clearly identified in our convergence map at the depicted scale size, but is well-pronounced within all explored $$\mathcal {F}_{{\rm aper}}{}$$ parametrizations. This may suggest that the aperture mass statistic might be better suited for substructure detection than a flexion-convergence inversion, by using the data more directly or through more suitable weighting profiles. Further study would be necessary in order to draw valid conclusions in this area, however. 4.2.3 Shear reconstructions We also compare the scales of structural predictions between flexion and shear by producing convergence reconstructions through direct shear inversion. Although our treatment of the PSF (using a simple Gaussian) would likely introduce ellipticity bias in an analysis on much wider scales, we expect a much stronger shear signal in the inner 3.5 arcmin × 3.5 arcmin cluster field. Fig. 5 shows our convergence reconstructions using just the averaged galaxy ellipticities measured through AIM fitting, again up to a mass-sheet degeneracy. Following a similar format, the top two graphs show the cluster and parallel field κ maps using a 15 arcsec shear smoothing kernel, while the bottom graphs utilize a 30 arcsec shear kernel. The clear broad central structure of the cluster peak at each scale is apparent, and both estimate the cluster centroid to be very close to what the flexion maps predict. However, both shear smoothing scales do not show any strong indications of local substructure. In the 15 arcsec map, the convergence extends slightly towards and hints at the eastern and north-eastern structures highlighted by our flexion analysis, but does not reproduce the flexion map's well-defined structural peaks, while the 30 arcsec reconstruction shows no indication of possible local structure at all. Employing our peak-finding algorithm to the 15 arcsec shear-convergence reconstruction with a σdetect = 2.5 detection threshold and a 1σbg background pixel-association threshold, we found the resulting region to have semimajor axis a = 50.22 ± 2.58 and semiminor axis b = 37.51 ± 1.29. As noted, the geometry values and uncertainties here refer to the reconstructed peak as detected in the κ map as reconstructed from 15 arcsec-smoothed shear data. Through the SL+WL analysis in Jauzac et al. (2016b), a WL model centred on the cluster core and SL region was found to have semimajor axis a = 60 arcsec, semiminor axis b = 42 arcsec, and orientation angle θ = 60°. This result is quite similar to our own reconstruction that finds an axis ratio of approximately 2:1, with an orientation angle of approximately 67°. To test the PSF assumptions and reconstruction method, we produced noise reconstructions of the shear in the parallel field in the same manner as the flexion procedure. Neither parallel field shear smoothing scale shows any appreciable structural significance, particularly when compared with the strong central profile of the cluster field, and the overall noise signal is centred at 0.0 similar to the flexion parallel field reconstructions. 4.2.4 Parallel reconstruction The parallel field convergence noise reconstruction under the same parameters is shown in the top right graph of Fig. 3. Compared to the cluster mass reconstruction, the derived magnitudes are larger but the noise peaks are much broader and do not correlate with any luminous matter. In addition, the magnitude of the highest value is approximately the same as that of the lowest, unlike our derived cluster convergence map. As we expect the parallel field to have 〈κ〉 ∼ 0, the observed |κ| parity around 0.0 is an indicator of the magnitude of systematic noise in an unconstrained flexion–convergence reconstruction. The parallel field map does not appear to have any significant peaks, beyond apparent edge effects (specifically at the north and northeast edges), nor do the minor peaks detected within the field (S/N ∼ 1) correlate with any luminous structure present. 4.3 Comparing flexion reconstructions to previous results One of the reasons that the galaxy cluster Abell 2744 was chosen to be observed as part of the HFF programme was that there already existed a wealth of observational information and a substantial amount of work detailed in the literature. One of the advantages of this strong existing framework is that it allows new methods and science to rest on, build from, and compare to studies and techniques that have already withstood academic rigour. While flexion signal would ideally be incorporated into larger and more encompassing combined approaches, it is also important to understand at a practical level how well flexion alone can reproduce the results of previous work on detecting complex and extended mass distributions. We compare our structure analysis with three works in particular: Merten et al. (2011), whose SL+WL reconstruction first detected four substructures arrayed around the core of the cluster; a later study by Medezinski et al. (2016) that independently detects four substructure components as well in a purely WL analysis, although not all structures detected by Merten were detected by Medezinski, and vice versa; and Jauzac et al. (2016b) who combine the techniques and results of the previous studies of Jauzac et al. (2015b) and Eckert et al. (2015) with additional information to form a broad view of the currently-known mass structure of Abell 2744, the summary of which can be found in Tables 1 and 2. Comparing with other estimates of Abell 2744, in a SL analysis Jauzac et al. (2015b) found a contained mass of M(<200 kpc) = (2.162 ± 0.005) × 1014 M⊙ within a 200 kpc aperture, which is similar to our estimate of the central potential. Although our relatively simple and straightforward mass calculation does not lie within the error bars of the precise SL estimate, Jauzac et al. stress that the precision of the mass models from gravitational-lensing studies depends strongly on the mass modelling technique, and mass estimates from different groups using different SL algorithms will find different results. As there are only a handful of flexion analyses of real galaxy data up to this point, it is difficult to determine if flexion reconstructions consistently over- or underestimate the mass of detected structure. Okura et al. (2008) obtain a convergence map of A1689 but do not estimate any enclosed masses. Cain et al. (2011) use aperture mass signal exclusively for reconstruction. Leonard et al. (2011) models aperture peaks with NFW and SIS models and obtain M200 estimates for these, but do not directly compare these masses with masses for those peaks attained in other works. As higher-order distortion measurements become more standard in gravitational lensing analyses we expect a clearer view of any biases to emerge. Reconstructions differ not only in terms of mass but also in terms of number of substructures. Using a joint SL+WL+X-ray analysis, Jauzac et al. (2016b) detect eight separate substructure peaks over the range of both the cluster core centre as well as the extended field that is outside of the observed HFF range, while our analysis only identifies three. Likewise, Merten et al. (2011) detect four distinct substructures within a field size of 600 arcsec × 600 arcsec, though the statistical significance of each is not made clear. While considering the question of matching the presence or absence of our detected substructure with other work, we identify a third major peak in our analysis that we previously ignored for its proximity to our field edge. However, the peak is over 25 arcsec from the edge, and the structure's had an elongated radial extent of over 50 arcsec, which brings a large portion of it well within an acceptable field of view. A particular motivation for the inclusion is the spatial correlation between our detected peak and an SL peak as identified in both Medezinski et al. (2016) and Jauzac et al. (2016b). 5 DISCUSSION We close with a few thoughts related to the differences between our work and previous analyses of A2744 (and other clusters). There are a few reasons why flexion analyses may not be able to pick up previously known mass peaks at any particular position, with the most prominent being that there are simply no (or very few) galaxy sources that both are close to the mass substructure and are suitable sources (e.g. large apparent size and good and realistic model fit) from which to measure a strong flexion signal. Because the flexion field is most prominent close to mass peaks and is buried by inherent shape noise farther away, if there are no chance source detections within a few radial arcseconds of the peak then it is likely that smaller substructure (on the 1012–1013 M⊙ scale) will remain undetected, without a substantial improvement in flexion measurement capability. Besides the structure correlated with potential DM peaks, our mass reconstructions also appear to capture a large portion of the gas distribution and X-ray surface brightness as examined and investigated independently in M11, M16, J16 using Chandra and XMM–Newton maps. The three analyses show generally consistent X-ray luminosity structure, with the global X-ray peak located just SW of the core BCG that expands broadly to the north-west. This structural orientation gives evidence towards the conclusion that the cluster had undergone a merger along the SE-NW axis that is supported by all three works along with other types of supporting evidence. In our own analysis, our derived mass distribution contours run parallel with J16's X-ray contours on a slanted line running from the core BCG NE towards the group of cluster galaxies located at approximately α = 3.596 5217°, δ = −30.388 257° (deg, J2000). While the gravitational shear field can be used effectively to determine the overall mass structure of galaxy clusters, its extended nonlocal effects as well as its inherent ellipticity-degeneracy limit its use to broad mass distributions, and thus it is not a viable candidate for higher resolution substructure detections. Additionally, while SL can and has led to precise characterizations of inner cluster cores, multiply-imaged source systems are not guaranteed to be located near substructure and quickly become sparse outside the dense core. As an intermediary signal spanning the strong and WL regimes, gravitational flexion signal has the ability to effectively probe significant cluster substructure on scales and at angular extents that cannot be practically detected through other means. As inherent flexion noise and systematic bias become more well-understood, flexion signal has the potential to be a key component in both exploring the behaviour of galaxy cluster formation and evolution as well as understanding the nature of DM structural dynamics. In this work, we have used an AIM implementation to measure flexion signal in the Abell 2744 galaxy cluster and inherent flexion in Abell 2744's associated parallel field. We show the efficacy of using flexion alone as an indicator of structure, exploring a much deeper view into the inner core of the cluster than shear would allow, and investigate the role of different mass estimators in both the cluster and parallel field. We identify and obtain mass estimates for both the central core of the cluster and a detected substructure offset 1.43 arcsec to the east of the core. Finally, we demonstrate that we are able to make simultaneous measurements of the shear field while measuring flexion through the AIM technique, and reconstruct the broader cluster mass structure while finding no such signal in the parallel field. Acknowledgements J.P.B. would like to thank both Austen Groener and Markus Rexroth for valuable discussions and dialog throughout this project. He would also like to recognize Mike Jarvis for his early coordination and support of the conceptual development and integration of the GalFlex flexion module as an affiliate into the GalSim software package. Finally, we thank the anonymous referee and editors for their constructive critiques, beneficial suggestions, and patience. The image data used in this paper are based on observations obtained with the NASA/ESA Hubble Space Telescope, retrieved from the MAST at the STScI. STScI is operated by the Association of Universities for Research in Astronomy, Inc. under NASA contract NAS 5-26555. This work also utilizes a subset of the gravitational lensing models produced by PIs Bradac, Natarajan & Kneib (CATS), Merten & Zitrin, Sharon, Williams, and the GLAFIC and Diego groups. These lens models have been partially funded by the HST Frontier Fields programmes conducted by STScI and obtained through MAST, as described above. Footnotes 1 https://archive.stsci.edu/pub/hlsp/frontier/abell2744/images/hst/ 2 http://www.stsci.edu/hst/observatory/focus/TinyTim 3 https://code.google.com/archive/p/astrolibpy/ 4 http://physics.drexel.edu/jbird/galflex/index.html REFERENCES Andrae R., Schulze-Hartung T., Melchior P., 2010, preprint (arXiv:1012.3754) Bacon D. J., Massey R. J., Refregier A. R., Ellis R. S., 2003, MNRAS , 344, 673 https://doi.org/10.1046/j.1365-8711.2003.06877.x CrossRef Search ADS   Bacon D. J., Goldberg D. M., Rowe B. T. P., Taylor A. N., 2006, MNRAS , 365, 414 https://doi.org/10.1111/j.1365-2966.2005.09624.x CrossRef Search ADS   Bartelmann M., Schneider P., 2001, Phys. Rep. , 340, 291 https://doi.org/10.1016/S0370-1573(00)00082-X CrossRef Search ADS   Bertin E., Arnouts S., 1996, A&AS , 117, 393 CrossRef Search ADS   Blandford R. D., Narayan R., 1992, ARA&A , 30, 311 https://doi.org/10.1146/annurev.aa.30.090192.001523 CrossRef Search ADS   Blanton M. R., Moustakas J., 2009, ARA&A , 47, 159 https://doi.org/10.1146/annurev-astro-082708-101734 CrossRef Search ADS   Bradac M., Schneider P., Lombardi M., Erben T., 2005, A&A , 437, 39 Bradač M., Lombardi M., Schneider P., 2004, A&A , 424, 13 https://doi.org/10.1051/0004-6361:20035744 CrossRef Search ADS   Bridle S. et al.  , 2010, MNRAS , 405, 2044 Cain B., Schechter P. L., Bautz M. W., 2011, ApJ , 736, 43 https://doi.org/10.1088/0004-637X/736/1/43 CrossRef Search ADS   Cain B., Bradač M., Levinson R., 2016, MNRAS  Caminha G. B. et al.  , 2016a, A&A , 587, A80 https://doi.org/10.1051/0004-6361/201527670 CrossRef Search ADS   Caminha G. B. et al.  , 2016b, A&A , 595, A100 https://doi.org/10.1051/0004-6361/201527995 CrossRef Search ADS   Caminha G. B. et al.  , 2017, A&A , 600, A90 https://doi.org/10.1051/0004-6361/201629297 CrossRef Search ADS   Capaccioli M., 1989, in Corwin H. G. Jr, Bottinelli L., eds, World of Galaxies (Le Monde des Galaxies) . Springer-Verlag, Berlin, p. 208 Google Scholar CrossRef Search ADS   Cardone V. F., Vicinanza M., Er X., Maoli R., Scaramella R., 2016, MNRAS , 462, 4028 https://doi.org/10.1093/mnras/stw1803 CrossRef Search ADS   Diego J. M., Broadhurst T., Molnar S. M., Lam D., Lim J., 2015a, MNRAS , 447, 3130 https://doi.org/10.1093/mnras/stu2660 CrossRef Search ADS   Diego J. M., Broadhurst T., Zitrin A., Lam D., Lim J., Ford H. C., Zheng W., 2015b, MNRAS , 451, 3920 https://doi.org/10.1093/mnras/stv1168 CrossRef Search ADS   Diego J. M. et al.  , 2016a, MNRAS , 473, 4279 https://doi.org/10.1093/mnras/stx2609 CrossRef Search ADS   Diego J. M., Broadhurst T. W J., Silk J., Lim J., Zheng W., Lam D., Ford H., 2016b, MNRAS , 459, 3447 https://doi.org/10.1093/mnras/stw865 CrossRef Search ADS   Diemand J., Kuhlen M., Madau P., Zemp M., Moore B., Potter D., Stadel J., 2008, Nature , 454, 735 https://doi.org/10.1038/nature07153 CrossRef Search ADS PubMed  Eckert D. et al.  , 2015, Nature , 528, 105 https://doi.org/10.1038/nature16058 CrossRef Search ADS PubMed  Er X., Tereno I., Mao S., 2012, MNRAS , 421, 1443 https://doi.org/10.1111/j.1365-2966.2012.20408.x CrossRef Search ADS   Falco E. E., Gorenstein M. V., Shapiro I. I., 1985, ApJ , 289, L1 https://doi.org/10.1086/184422 CrossRef Search ADS   Gao L., De Lucia G., White S. D. M., Jenkins A., 2004, MNRAS , 352, L1 https://doi.org/10.1111/j.1365-2966.2004.08098.x CrossRef Search ADS   Goldberg D. M., Bacon D. J., 2005, ApJ , 619, 741 https://doi.org/10.1086/426782 CrossRef Search ADS   Goldberg D. M., Leonard A., 2007, ApJ , 660, 1003 https://doi.org/10.1086/513137 CrossRef Search ADS   Goldberg D. M., Natarajan P., 2002, ApJ , 564, 65 https://doi.org/10.1086/324202 CrossRef Search ADS   Gorenstein M. V., Shapiro I. I., Falco E. E., 1988, ApJ , 327, 693 https://doi.org/10.1086/166226 CrossRef Search ADS   Grillo C. et al.  , 2015, ApJ , 800, 38 https://doi.org/10.1088/0004-637X/800/1/38 CrossRef Search ADS   Grillo C. et al.  , 2016, ApJ , 822, 78 https://doi.org/10.3847/0004-637X/822/2/78 CrossRef Search ADS   Harvey D., Massey R., Kitching T., Taylor A., Tittley E., 2015, Science , 347, 1462 https://doi.org/10.1126/science.1261381 CrossRef Search ADS PubMed  Heymans C. et al.  , 2006, MNRAS , 368, 1323 https://doi.org/10.1111/j.1365-2966.2006.10198.x CrossRef Search ADS   Hoekstra H., Bartelmann M., Dahle H., Israel H., Limousin M., Meneghetti M., 2013, Space Sci. Rev. , 177, 75 https://doi.org/10.1007/s11214-013-9978-5 CrossRef Search ADS   Jauzac M. et al.  , 2014, MNRAS , 443, 1549 https://doi.org/10.1093/mnras/stu1355 CrossRef Search ADS   Jauzac M. et al.  , 2015a, MNRAS , 446, 4132 https://doi.org/10.1093/mnras/stu2425 CrossRef Search ADS   Jauzac M. et al.  , 2015b, MNRAS , 452, 1437 https://doi.org/10.1093/mnras/stv1402 CrossRef Search ADS   Jauzac M. et al.  , 2016a, MNRAS , 457, 2029 https://doi.org/10.1093/mnras/stw069 CrossRef Search ADS   Jauzac M. et al.  , 2016b, MNRAS , 463, 3876 https://doi.org/10.1093/mnras/stw2251 CrossRef Search ADS   Johnson T. L., Sharon K., Bayliss M. B., Gladders M. D., Coe D., Ebeling H., 2014, ApJ , 797, 48 https://doi.org/10.1088/0004-637X/797/1/48 CrossRef Search ADS   Kaiser N., Squires G., 1993, ApJ , 404, 441 https://doi.org/10.1086/172297 CrossRef Search ADS   Kaiser N., Squires G., Broadhurst T., 1995, ApJ , 449, 460 https://doi.org/10.1086/176071 CrossRef Search ADS   Kitching T. D. et al.  , 2012, MNRAS , 423, 3163 https://doi.org/10.1111/j.1365-2966.2012.21095.x CrossRef Search ADS   Klypin A., Kravtsov A. V., Valenzuela O., Prada F., 1999, ApJ , 522, 82 https://doi.org/10.1086/307643 CrossRef Search ADS   Kneib J.-P., Natarajan P., 2011, A&AR , 19, 47 CrossRef Search ADS   Kravtsov A. V., Berlind A. A., Wechsler R. H., Klypin A. A., Gottlöber S., Allgood B., Primack J. R., 2004a, ApJ , 609, 35 https://doi.org/10.1086/420959 CrossRef Search ADS   Kravtsov A. V., Gnedin O. Y., Klypin A. A., 2004b, ApJ , 609, 482 https://doi.org/10.1086/421322 CrossRef Search ADS   Lam D., Broadhurst T., Diego J. M., Lim J., Coe D., Ford H. C., Zheng W., 2014, ApJ , 797, 98 https://doi.org/10.1088/0004-637X/797/2/98 CrossRef Search ADS   Leonard A., Goldberg D. M., Haaga J. L., Massey R., 2007, ApJ , 666, 51 https://doi.org/10.1086/520109 CrossRef Search ADS   Leonard A., King L. J., Goldberg D. M., 2011, MNRAS , 413, 789 https://doi.org/10.1111/j.1365-2966.2010.18171.x CrossRef Search ADS   Lewis A., Challinor A., 2006, Phys. Rep. , 429, 1 https://doi.org/10.1016/j.physrep.2006.03.002 CrossRef Search ADS   Lotz J. M. et al.  , 2017, ApJ , 837, 97 https://doi.org/10.3847/1538-4357/837/1/97 CrossRef Search ADS   Mandelbaum R., Seljak U., Kauffmann G., Hirata C. M., Brinkmann J., 2006, MNRAS , 368, 715 https://doi.org/10.1111/j.1365-2966.2006.10156.x CrossRef Search ADS   Mandelbaum R. et al.  , 2015, MNRAS , 450, 2963 https://doi.org/10.1093/mnras/stv781 CrossRef Search ADS   Massey R. et al.  , 2007a, MNRAS , 376, 13 https://doi.org/10.1111/j.1365-2966.2006.11315.x CrossRef Search ADS   Massey R., Rowe B., Refregier A., Bacon D. J., Bergé J., 2007b, MNRAS , 380, 229 https://doi.org/10.1111/j.1365-2966.2007.12072.x CrossRef Search ADS   Massey R., Kitching T., Richard J., 2010, Rep. Prog. Phys. , 73, 086901 https://doi.org/10.1088/0034-4885/73/8/086901 CrossRef Search ADS   Medezinski E. et al.  , 2013, ApJ , 777, 43 https://doi.org/10.1088/0004-637X/777/1/43 CrossRef Search ADS   Medezinski E., Umetsu K., Okabe N., Nonino M., Molnar S., Massey R., Dupke R., Merten J., 2016, ApJ , 817, 24 https://doi.org/10.3847/0004-637X/817/1/24 CrossRef Search ADS   Merten J. et al.  , 2011, MNRAS , 417, 333 https://doi.org/10.1111/j.1365-2966.2011.19266.x CrossRef Search ADS   Mo W. et al.  , 2016, ApJ , 818, L25 https://doi.org/10.3847/2041-8205/818/2/L25 CrossRef Search ADS   Mohammed I., Saha P., Williams L. L. R., Liesenborgs J., Sebesta K., 2016, MNRAS , 459, 1698 https://doi.org/10.1093/mnras/stw727 CrossRef Search ADS   Munshi D., Valageas P., van Waerbeke L., Heavens A., 2008, Phys. Rep. , 462, 67 https://doi.org/10.1016/j.physrep.2008.02.003 CrossRef Search ADS   Nagai D., Kravtsov A. V., 2005, ApJ , 618, 557 https://doi.org/10.1086/426016 CrossRef Search ADS   Okamoto T., Nagashima M., 2001, ApJ , 547, 109 https://doi.org/10.1086/318375 CrossRef Search ADS   Okura Y., Umetsu K., Futamase T., 2008, ApJ , 680, 1 https://doi.org/10.1086/587676 CrossRef Search ADS   Owers M. S., Randall S. W., Nulsen P. E. J., Couch W. J., David L. P., Kempner J. C., 2011, ApJ , 728, 27 https://doi.org/10.1088/0004-637X/728/1/27 CrossRef Search ADS   Press W. H., Schechter P., 1974, ApJ , 187, 425 https://doi.org/10.1086/152650 CrossRef Search ADS   Rexroth M., Natarajan P., Kneib J.-P., 2016, MNRAS , 460, 2505 https://doi.org/10.1093/mnras/stw1017 CrossRef Search ADS   Rhodes J. D. et al.  , 2007, ApJS , 172, 203 https://doi.org/10.1086/516592 CrossRef Search ADS   Richard J. et al.  , 2014, MNRAS , 444, 268 https://doi.org/10.1093/mnras/stu1395 CrossRef Search ADS   Rix H.-W. et al.  , 2004, ApJS , 152, 163 https://doi.org/10.1086/420885 CrossRef Search ADS   Schneider P., 1996, MNRAS , 283, 837 https://doi.org/10.1093/mnras/283.3.837 CrossRef Search ADS   Schneider P., Er X., 2008, A&A , 485, 363 https://doi.org/10.1051/0004-6361:20078631 CrossRef Search ADS   Schneider P., Seitz C., 1995, A&A , 294, 411 Schneider P., van Waerbeke L., Jain B., Kruse G., 1998, MNRAS , 296, 873 https://doi.org/10.1046/j.1365-8711.1998.01422.x CrossRef Search ADS   Schrabback T. et al.  , 2010, A&A , 516, A63 https://doi.org/10.1051/0004-6361/200913577 CrossRef Search ADS   Sebesta K., Williams L. L. R., Mohammed I., Saha P., Liesenborgs J., 2016, MNRAS , 461, 2126 https://doi.org/10.1093/mnras/stw1433 CrossRef Search ADS   Sendra I., Diego J. M., Broadhurst T., Lazkoz R., 2014, MNRAS , 437, 2642 https://doi.org/10.1093/mnras/stt2076 CrossRef Search ADS   Sharon K., Johnson T. L., 2015, ApJ , 800, L26 https://doi.org/10.1088/2041-8205/800/2/L26 CrossRef Search ADS   Sheth R. K., Tormen G., 1999, MNRAS , 308, 119 https://doi.org/10.1046/j.1365-8711.1999.02692.x CrossRef Search ADS   Springel V., White S. D. M., Tormen G., Kauffmann G., 2001, MNRAS , 328, 726 https://doi.org/10.1046/j.1365-8711.2001.04912.x CrossRef Search ADS   Tinker J., Kravtsov A. V., Klypin A., Abazajian K., Warren M., Yepes G., Gottlöber S., Holz D. E., 2008, ApJ , 688, 709 https://doi.org/10.1086/591439 CrossRef Search ADS   Treu T. et al.  , 2016, ApJ , 817, 60 https://doi.org/10.3847/0004-637X/817/1/60 CrossRef Search ADS   Van Waerbeke L. et al.  , 2000, A&A , 358, 30 Viola M., Melchior P., Bartelmann M., 2012, MNRAS , 419, 2215 https://doi.org/10.1111/j.1365-2966.2011.19872.x CrossRef Search ADS   Wang X. et al.  , 2015, ApJ , 811, 29 https://doi.org/10.1088/0004-637X/811/1/29 CrossRef Search ADS   Wittman D. M., Tyson J. A., Kirkman D., Dell'Antonio I., Bernstein G., 2000, Nature , 405, 143 https://doi.org/10.1038/35012001 CrossRef Search ADS PubMed  Zitrin A., Menanteau F., Hughes J. P., Coe D., Barrientos L. F., Infante L., Mandelbaum R., 2013, ApJ , 770, L15 https://doi.org/10.1088/2041-8205/770/1/L15 CrossRef Search ADS   Zitrin A. et al.  , 2015, ApJ , 801, 44 https://doi.org/10.1088/0004-637X/801/1/44 CrossRef Search ADS   © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society

### Journal

Monthly Notices of the Royal Astronomical SocietyOxford University Press

Published: May 1, 2018

## You’re reading a free preview. Subscribe to read the entire article.

### DeepDyve is your personal research library

It’s your single place to instantly
that matters to you.

over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month ### Explore the DeepDyve Library ### Search Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly ### Organize Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place. ### Access Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals. ### Your journals are on DeepDyve Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more. All the latest content is available, no embargo periods. DeepDyve ### Freelancer DeepDyve ### Pro Price FREE$49/month
\$360/year

Save searches from
PubMed

Create lists to

Export lists, citations