A test for skewed distributions of dark matter and a possible detection in galaxy cluster Abell 3827

TL;DR

This paper introduces a new gravitational lensing method to detect asymmetries in dark matter distributions, applied to galaxy cluster Abell 3827, providing tentative evidence of skewness consistent with self-interacting dark matter models.

Contribution

The paper develops a gravitational lensing technique incorporating a skewness parameter to identify asymmetric dark matter distributions in galaxy clusters.

Findings

Tentative detection of skewness in dark matter distribution in Abell 3827

Skewness measurement aligns with previous offset findings

Method shows potential for future dark matter collision studies

Abstract

Simulations of self-interacting dark matter (SIDM) predict that dark matter should lag behind galaxies during a collision. If the interaction is mediated by a high-mass force carrier, the distribution of dark matter can also develop asymmetric dark matter tails. To search for this asymmetry, we compute the gravitational lensing properties of a mass distribution with a free {\em skewness} parameter. We apply this to the dark matter around the four central galaxies in cluster Abell~3827. In the galaxy whose dark matter peak has previously been found to be offset, we tentatively measure a skewness $s=0.23^{+0.05}_{-0.22}$ in the same direction as the peak offset. Our method may be useful in future gravitational lensing analyses of colliding galaxy clusters and merging galaxies.

Tables

Table 1. Table 1: Parameters of the fiducial mass model fitted by Lenstool . Quantities in square brackets are fixed. Errors on other quantities show 68% statistical confidence limits, marginalising over uncertainty in all other parameters. Stellar mass components are modelled as Hernquist profiles, with a mass (computed from flux in the F606W band), scale radius and ellipticity (fitted using Galfit ; galaxy N4 is contaminated by a nearby star). Dark matter components are modelled as PISPs, with a 1D velocity dispersion, core and cut radii, ellipticity and skewness. Positions are given in arcseconds relative to (R.A.: 4330.47515 4330.47515 4330.47515 , Dec.: − 59.945996 59.945996 -59.945996 ), except galaxies’ dark matter components, which are relative to the position of their stars. Angles are anticlockwise from East.

$x$[″]

$y$[″]

Mass [$M_{\odot }$]

$r_{\mathrm{sc}}$[″]

$ϵ$

${\varphi }_{ϵ}$[∘]

$s$

${\varphi }_s$[∘]

$\Delta x$[″]

$\Delta y$[″]

${\sigma }_v$[km/s]

$r_{\mathrm{core}}\left[″\right]$

​​$r_{\mathrm{cut}}\left[″\right]$​​

N1

stars

[$-0.06$

]

[$0.04$

]

$\left[1.00\times {10}^{11}\right]$

[0.53]

[$0.12$

]

[$61$

]

dark matter

$-0.29$

${}_{-0.14}{}^{}{}_{}{}^{+0.25}$

$-0.71$

${}_{-0.16}{}^{}{}_{}{}^{+0.30}$

$149$

${}_{-12}{}^{}{}_{}{}^{+8}$

[$0.1$]

[$40$]

$0.02$

${}_{-0.01}{}^{}{}_{}{}^{+0.33}$

$151$

${}_{-116}{}^{}{}_{}{}^{+19}$

$0.21$

${}_{-0.22}{}^{}{}_{}{}^{+0.06}$

$86$

${}_{-44}{}^{}{}_{}{}^{+44}$

N2

stars

[$5.07$

]

[$2.05$

]

$\left[2.46\times {10}^{11}\right]$

[0.79]

[$0.17$

]

[$39$

]

dark matter

$-0.23$

${}_{-0.16}{}^{}{}_{}{}^{+0.30}$

$0.00$

${}_{-0.30}{}^{}{}_{}{}^{+0.30}$

$182$

${}_{-22}{}^{}{}_{}{}^{+29}$

[$0.1$]

[$40$]

$0.42$

${}_{-0.22}{}^{}{}_{}{}^{+0.05}$

$23$

${}_{-12}{}^{}{}_{}{}^{+32}$

$0.03$

${}_{-0.14}{}^{}{}_{}{}^{+0.11}$

$117$

${}_{-80}{}^{}{}_{}{}^{+41}$

N3

stars

[$9.69$

]

[$3.98$

]

$\left[2.77\times {10}^{11}\right]$

[0.33]

[$0.05$

]

[$31$

]

dark matter

$-0.05$

${}_{-0.25}{}^{}{}_{}{}^{+0.25}$

$-0.06$

${}_{-0.29}{}^{}{}_{}{}^{+0.18}$

$213$

${}_{-10}{}^{}{}_{}{}^{+8}$

[$0.1$]

[$40$]

$0.49$

${}_{-0.16}{}^{}{}_{}{}^{+0.01}$

$15$

${}_{-8}{}^{}{}_{}{}^{+14}$

$-0.02$

${}_{-0.11}{}^{}{}_{}{}^{+0.08}$

$169$

${}_{-109}{}^{}{}_{}{}^{+7}$

N4

stars

[$9.26$

]

[$-1.08$

]

$\left[2.08\times {10}^{11}\right]$

[1.37]

[$0.39$

]

[$127$

]

dark matter

$-1.35$

${}_{-0.34}{}^{}{}_{}{}^{+0.39}$

$0.51$

${}_{-0.27}{}^{}{}_{}{}^{+0.35}$

$255$

${}_{-10}{}^{}{}_{}{}^{+8}$

[$0.1$]

[$40$]

$0.02$

${}_{-0.01}{}^{}{}_{}{}^{+0.25}$

$136$

${}_{-28}{}^{}{}_{}{}^{+17}$

$0.08$

${}_{-0.09}{}^{}{}_{}{}^{+0.08}$

$147$

${}_{-80}{}^{}{}_{}{}^{+21}$

N6

stars

[$18.54$

]

[$2.47$

]

[0]

dark matter

[$0$

]

[0]

$38$

${}_{-25}{}^{}{}_{}{}^{+26}$

[$0.1$]

[$40$]

[$0$

]

[$0$

]

[$0$]

[$0$]

Cluster   dm

$5.53$

${}_{-1.61}{}^{}{}_{}{}^{+1.46}$

$2.33$

${}_{-1.59}{}^{}{}_{}{}^{+1.97}$

$683$

${}_{-75}{}^{}{}_{}{}^{+139}$

​​$30.12$

${}_{-6.43}{}^{}{}_{}{}^{+9.23}$​​

[$1000$]

$0.56$

${}_{-0.10}{}^{}{}_{}{}^{+0.13}$

$63$

${}_{-3}{}^{}{}_{}{}^{+2}$

[$0$]

[$0$]

Equations

$$\Sigma (\mathbf{r})\equiv \frac{{\sigma }_0^2{r}_{\text{cut}}}{2G\left(r_{\text{cut}}-r_{\text{core}}\right)}\left(\frac{1}{\sqrt{r_{\text{core}}^2+r^2}}-\frac{1}{\sqrt{r_{\text{cut}}^2+r^2}}\right),$$

$$x\to x_{\text{em}}=\frac{x}{1+ϵ},y\to y_{\text{em}}=\frac{y}{1-ϵ}.$$

$$r^{2}\rightarrow r_{\text{em}}^{2}=\frac{r^{2}}{\left(1-\epsilon^{2}\right)^{2}}\left[1+\epsilon^{2}-2\epsilon\cos{\big{(}2(\theta-\phi_{\epsilon})\big{)}}\right].$$

$$\alpha \left(\mathbf{r}\right)=$$

$$\psi \left(\mathbf{r}\right)=$$

$$\alpha \left(r\right)=$$

$$\Sigma \left(r\right)=$$

$$\psi \left(x,y\right)\to {\psi }^{\prime }\left(x,y\right)\equiv \psi \left(x^{\prime },y^{\prime }\right).$$

$$\psi^{\prime}_{x}=\left(\psi_{x^{\prime}}\left(x^{\prime},y^{\prime}\right)x^{\prime}_{x}+\psi_{y^{\prime}}\left(x^{\prime},y^{\prime}\right)y^{\prime}_{x}\right)\big{|}_{\left(x,y\right)},$$

$$r^2\to r^{\prime 2}=\frac{r^2{\left(1-s^2\right)}^{3/2}}{{\left(1+s\cos [\theta -{\varphi }_s]\right)}^2}$$

$$w=\frac{|s|}{\sqrt{0.3^2-s^2{\cos }^2({\varphi }_s-\varphi )}}.$$

$${\rho }_{\text{star}}(r)=\frac{{\rho }_s}{\left(r/r_s\right){\left(1+r/r_s\right)}^3},$$

$$(ma{g}_{F814W}-ma{g}_{F606W})=0.022\times ma{g}_{F606W}-\mathrm{1.129.}$$

$$\psi (r)\to {\psi }^{\prime }(r)\equiv w(r;\{a_i\})\psi (r),$$

$$w(r;\{a_i\})=1+sf(r,\theta )$$

$$w(r;s,r_s,\beta ,{\varphi }_s)=1+s{\text{ tan}}^{-1}\left(\beta \frac{r}{r_s}\right)\text{cos}(\theta -{\varphi }_s),$$

$$\frac{c^2}{8\pi G}\frac{D_s}{D_l{D}_{ls}}{\int }_{\mid r\mid <R}{\nabla }^2\psi \text{ d}A=\frac{c^2}{8\pi G}\frac{D_s}{D_l{D}_{ls}}{\int }_{\mid r\mid <R}{\nabla }^2\left(w\psi \right)\text{ d}A,$$

$$\Sigma \left(r\right)={\Sigma }_{{ϵ}_1}\left(r\right)w_1\left(r\right)+{\Sigma }_{{ϵ}_2}\left(r\right)w_2\left(r\right)$$

$$w_1\left(r\right)=$$

$$w_2\left(r\right)=$$

Full text

††pagerange: A test for skewed distributions of dark matter, and a possible detection in galaxy cluster Abell 3827–A.2††pubyear: 2014

1 Introduction

Most of the mass in the Universe is dark matter (e.g. Planck Collaboration, 2016). Dark matter appears invisible, because it does not interact (or interacts very weakly) with Standard Model particles including photons.

As the nature of dark matter remains unknown, there is no reason to a priori assume a particular theory of its origin. The wide range of proposed dark matter models predict different spatial distributions, particularly on small scales. Dark matter particles that interact with each other (SIDM) were proposed in (Spergel & Steinhardt, 2000) to explain small scale discrepancies between observations and simulations of collisionless dark matter. In the SIDM paradigm, energy transfer between particles makes the centre of galaxies (Vogelsberger et al., 2012) and galaxy clusters (Rocha et al., 2013) more circular and less dense, potentially resolving the core/cusp problem. Small substructures can also be erased – leading to the observed underabundance of galaxies in the Local Group, relative to simulations. During mergers between galaxies or galaxy clusters, dark matter interactions transfer momentum between the colliding dark matter haloes (Randall et al., 2008; Kahlhoefer et al., 2014; Robertson et al., 2017; Kim et al., 2016; Robertson et al., 2016). These scatterings can temporarily separate dark matter from its associated galaxies. Such dark matter lags behind the galaxies, toward the position of diffuse gas that is slowed by ram pressure (Clowe et al., 2004; Lage & Farrar, 2014; Harvey et al., 2015). Scattering processes during collisions can be seperated into two types: frequent low momentum transfer and infrequent high momentum transfer. These will have different qualitative behaviours.

Frequent low momentum transfer scattering will cause an effective drag force, which if greater than the gravitational restoring force, will seperate the entire DM halo from the galaxy during collisions. Crucially there will be no tail of scattered DM particles escaping the potential well (Kahlhoefer et al., 2014). Numerous studies have placed constraints on the cross-section of DM in this regime. Measuring dark matter galaxy offsets on a sample of 72 merging clusters, (Harvey et al., 2015) found $\tilde{\sigma}/{m_{\mathrm{DM}}}<0.5$ cm2g*-1*. Constraints from the Bullet Cluster place $\tilde{\sigma}/{m_{\mathrm{DM}}}\gtrsim 1.2$ cm2g*-1* (Kahlhoefer et al., 2014), while constraints from an offset galaxy in Abell 3827 yields $\tilde{\sigma}/{m_{\mathrm{DM}}}\gtrsim 2.0$ cm2g*-1* (Kahlhoefer et al., 2015).

In contrast infrequent high momentum transfer scattering (mediated by a high-mass force carrier, for example) will cause a small fraction of scattered particle to leave the potential well on the trailing side. Shortly after collision, this will appear as a tail of scattered DM particles (see Figure 5 in Kahlhoefer et al., 2014). Although the peak of the DM distribution will remain conincident with the galaxy, the tail of DM particles will lead to an apparant shift in the centre (Kahlhoefer et al., 2014).

Gravitational lensing offers the most direct way to map the spatial distribution of dark matter, and hence to infer its particle properties. Gravitational lensing refers to the deflection of light rays passing near any mass, including dark matter. Thanks to this deflection, (unrelated) objects behind dark matter appear characteristically distorted, or even visible along more than one (curved) line-of-sight. Even though dark matter is invisible, it is possible to invert this process and infer where it must be, by undistorting the observed images, or ray-tracing multiple images back onto each other.

Galaxy cluster Abell 3827 ($22$h $01\arcmin$ $49\aas@@fstack{\prime\prime}1$ $-59^{\circ}$ $57\arcmin$ $15\arcsec$, redshift $z$=0.099) is particularly well suited for this kind of study. It gravitationally lenses a $z$=1.24 galaxy with spiral arms and several knots of star formation that can be treated as independent background sources (Massey et al., 2015, hereafer M15). While most clusters contain only one brightest central galaxy, Abell 3827 contains four equally-bright galaxies within its central 10 kpc (Carrasco et al., 2010; Williams & Saha, 2011). This highly unusual configuration means that some of the galaxies appear close to gravitationally lensed images. Thus, under parametric model assumptions, the distribution of the dark matter can be measured. Because of the cluster’s relative proximity (in terms of gravitational lensing), it is possible to resolve small spatial offsets between the distribution of dark matter and stars in the foreground galaxies.

In this paper we present a new parametric lensing approach to search for the predicted asymmetry in the distribution of dark matter during mergers. A previous search looked for residuals after subtracting the symmetric component (Harvey et al., 2017), but that may be less sensitive because a tail of scattered particles shifts the best-fit position of the symmetric component backwards, thus removing some of the residual. We instead construct a single halo model with a free skewness parameter that qualitatively captures the asymmetry found in high momentum transfer scattering simulations. We implement and distribute this model in the publicly available Lenstool software111http://projets.lam.fr/projects/lenstool/wiki (Jullo et al., 2007). We test it on both mock data, where the skewness of the lens is known a priori, and on Abell 3827. Section 2 describes existing observations of Abell 3827. Section 3 introduces our new parametric lens model. Section 4 contains an analysis of Abell 3827. To be consistent with M15 we assume throughout this paper a flat $\Lambda$CDM cosmology with $\Omega_{M}$=0.3, $\Omega_{\Lambda}$=0.7 and $H_{0}$=70 km s*-1* Mpc -1. At the redshift of Abell 3827, $1\arcsec$ corresponds to $1.828$ kpc.

2 Data

Broad-band imaging of Abell 3827 has been obtained from the Gemini telescope at optical wavelengths (Carrasco et al., 2010) and from the Hubble Space Telescope (HST) in the UV, optical and near-infrared (M15). This revealed four similarly-bright elliptical galaxies within 10 kpc of each other, plus a background spiral galaxy, whose multiply-lensed images are threaded throughout the cluster core (figure 1).

Integral Field Unit spectroscopy has been obtained from the VLT. An initial 1 hour exposure with the Multi-Unit Spectroscopic Explorer (MUSE) identified four main groups of lensed images, and suggested two low S/N peaks as candidates for a demagnified central image (M15). A subsequent additional 4 hour exposure (programme 295.A-5018; Massey et al. in prep.) confirms both candidates (Ao7 at RA: $330.47047$, Dec: $-59.945183$, Ao8 at RA: $330.47079$, Dec: $-59.946112$). Indeed, Ao7 is also visible in HST imaging, after using the MuSCADeT multiwavelength method (Joseph et al., 2016) to estimate and subtract bright foreground emission (Figure 1). We therefore use all the images identified by M15, plus the two new ones.

3 Lens modelling method

We shall model the distribution of mass in the galaxy cluster as a sum of cluster-scale plus galaxy-scale halos (following e.g. Limousin et al., 2007a; Jauzac et al., 2014), each a perturbation around the Pseudo Isothermal Elliptical Mass Distribution (PIEMD; Kassiola & Kovner, 1993).

3.1 Pseudo-Isothermal Skewed Potential

The 2D surface mass density $\Sigma$ of a circularly symmetric pseudo-isothermal mass distribution, projected along a line of sight, is:

[TABLE]

where $\sigma_{0}$ is the 1D velocity dispersion, and where $r_{\text{core}}$ ($r_{\text{cut}}$) is an inner (outer) radius. To convert this into a PIEMD with ellipticity $\epsilon=\frac{a-b}{a+b}\geq 0$, where $a$ and $b$ are the semi-major and semi-minor axes respectively, Kassiola & Kovner (1993) apply their coordinate transformation (2.3.6):

[TABLE]

This maps a circle onto an ellipse centered at the origin, with its major axis along the $x$ axis. In general, including a rotation to set the major axis at angle $\phi_{\epsilon}$, this can be expressed in polar coordinates as:

[TABLE]

The angle $\alpha$ by which a light ray is deflected as it passes near the lens, and the equivalent 2D gravitational potential $\psi$ can be computed by integrating the density distribution:

[TABLE]

where $D_{l}$, $D_{s}$ and $D_{ls}$ are the angular diameter distance from the observer to the lens, observer to the source, and lens to the source respectively. For general mass distributions, these integrals are difficult to solve – but closed forms have been found for the PIEMD, using techniques from complex analysis that exploit its elliptical symmetry (Bourassa & Kantowski, 1975).

A related halo model is the Pseudo Isothermal Elliptical Potential (PIEP; Kassiola & Kovner, 1993). In this, the coordinate transformation is applied to a circular potential $\psi$ (rather than the density). It is then mathematically easier to obtain the deflection angle and density via differentiation:

[TABLE]

In detail, the PIEP potential $\psi$ is transformed so that:

[TABLE]

The first and second derivatives can then be computed with applications of the chain rule. For example, the first $x$-derivative of the potential is:

[TABLE]

where the subscript denotes partial differentiation. The resulting mass distribution is not the same as a PIEMD, because of the way the coordinate transformation propagates through the chain rule (or back up the integrals in equation 4). For large $\epsilon$, the mass distribution corresponding to a PIEP has undesirable features including concave (peanut-shaped) isodensity contours (Kassiola & Kovner, 1993).

3.2 Pseudo-Isothermal Skewed Mass Distribution

To perturb the mass distribution in a way that resembles the behaviour of SIDM in numerical simulations (see figure 5 of Kahlhoefer et al., 2014), we apply a further coordinate transformation that maps a circle onto an ellipse with its focus (rather than centre) at the origin:

[TABLE]

with $s$ being the third eccentricity such that $s=\sqrt{1-b^{2}/a^{2}}$, and the power $3/2$ being introduced to preserve area. Note the asymmetric $\cos{(\theta)}$ terms rather than the $\cos{(2\theta)}$ terms in the mapping described by equation (3).

We apply this transformation to the 2D gravitational potential corresponding to the PIEMD.222We would ideally apply this transformation to the PIEMD mass distribution, but the relevant integrals (equation 4) do not contract to a simple form. A skewed mass distribution could also be derived from the potential corresponding to a PIEP. We choose to perturb the PIEMD so that we recover this widely-used mass distribution in the $s\rightarrow 0$ limit, and to minimise undesired convex curvature in density isophotes. Analytic (albeit cumbersome) expressions for deflection angle and density can be readily calculated via differentiation (equation 5). We denote this the Pseudo Isothermal Skewed Potential (PISP); its isodensity contours are shown for various values of $\epsilon$ and $s$ in Figure 2.

For any skew, the peak of the mass distribution always lies at the same position, so it will be possible to use the same metric as M15 to measure any offset between the most gravitationally bound stars and dark matter. The total mass changes slightly with increasing skew, but this can be recalculated after a fit.

Like the PIEP, the density distribution of the PISP exhibits undesired behaviour with large skews, because the coordinate transformation was applied to $\psi$, not $\Sigma$. Isophotes of the density distribution become concave, and the skew ellipticity can overwhelm the underlying ellipticity. We avoid these effects by restricting $|s|<0.3$.333See the Appendix for an alternative model that does not suffer from this effect, but has other disadvantages. Since the PISP is invariant under transformations $s\rightarrow-s$ and $\phi_{s}\rightarrow\phi_{s}+\pi$, we fit parameters within the domain $s\in[-0.3,0.3]$ and $\phi_{s}$ in some interval of length $\pi$. This ensures that the parameter space can be explored symmetrically about $s=0$, allowing the case of zero skew to be recovered without bias. We set the edge of the domain of $\phi_{s}$ well away from any preferred direction (in practice, having explored parameter space via a quick search), to make sure an MCMC sampler operates efficiently near regions of interest.

3.3 Testing an implementation in Lenstool

We have implemented the PISP as potential 813 in the publicly-available software Lenstool (Jullo et al., 2007). Given a parameterised mass distribution, and the location of background sources, Lenstool can compute the position of observed multiple images. Given the position of observed images, it can also use Markov-Chain Monte Carlo (MCMC) optimisation to fit parameters of the lensing mass distribution.

To test whether Lenstool can accurately recover a known input skew, we run two sets of tests. We first consider an isolated lens, with three background sources at different redshifts: the example_with_images configuration that is packaged with Lenstool. As a null test, we adopt the input mass distribution with skew $s_{\mathrm{true}}=0$. From the observed positions of multiple images, Lenstool successfully recovered a best fit (maximum likelihood) value $s=-0.0008^{+0.02}_{-0.02}.$

We then set skew $s_{\mathrm{true}}=0.2$ and $\phi_{s_{\mathrm{true}}}=1.6$. We set source positions by projecting one image of each source back to its source plane, then create a mock set of multiply lensed images by re-projecting this source forward through the lens. When fitting this mock data, Lenstool successfully recovers best-fit values $s=0.2^{+0.001}_{-0.001}$ and $\phi_{s}=1.6^{+0.04}_{-0.05}$.

Second, we test the recovery of input skews in a complex cluster lens with a mass distribution based on the Abell 3827. Choosing one of the quadruply lensed background galaxy images, we repeat the procedure outlined above: projecting the light backwards and then forwards through a cluster lens with known mass distribution. The cluster is given the same parameters as our fiducial model for Abell 3827 (see §4.1), with the exception of the skew parameters. In this test, the dark matter associated with galaxy N1 is given skewness $s_{\mathrm{true}}=0.25$ and $\phi_{s_{\mathrm{true}}}=1.6$. As a null test, galaxies N2–N4 are set to have no skew, $s_{\mathrm{true}}=0$.

We run Lenstool with the same free parameters and priors as in §4.1). Within such a highly dimensional parameter space, we find that the best-fit values are sometimes noisy, for parameters that make only a small difference to the overall goodness of fit. However, the full posterior probability distribution function (PDF) is smooth and well-sampled. Hence, for the rest of this paper, we shall quote the modal peak and 68% width of the posterior PDF, which Lenstool also returns. This makes no difference for the simple model above, and successfully recovers $s=0.24^{+0.04}_{-0.31}$ and $\phi_{s}=1.6^{0.92}_{-0.99}$ for galaxy N1, and $s=0.01^{+0.14}_{-0.13}$, $0.07^{+0.10}_{-0.15}$, $0.11^{+0.11}_{-0.16}$ for galaxies N2, N3, N4 (with very weakly constrained $\phi_{s}$).

3.4 Prior bias for polar parameters

A skew is a two-component vector, and can be expressed in polar form as a magnitude $|s|$ and direction $\phi_{s}$, or in Cartesian form as an amount in orthogonal directions $\{s_{x},s_{y}\}$. We implemented the polar option, so that Lenstool’s MCMC optimiser can explore a circularly symmetric region, with no preferred direction that could bias the inferred skew. Lenstool also defines ellipticities in this way, for the same reason.

Nonetheless, it may often be desirable to know the posterior probability distribution of skewness along e.g. a direction of motion, and perpendicular to that (i.e. the Cartesian form). The posterior probability distributions of skew and skew angle are returned by Lenstool (in runmode=3) by the sampling density of the MCMC chain. This can be converted to the posterior of the skew in some direction $\phi$ by projecting and then weighting each sample by:

[TABLE]

The numerator is the Jacobian to convert the area of parameter space from polar to Cartesian coordinates. The denominator corrects for prior bias, because the restriction $|s|\in[-0.3,0.3]$ leads to a (semi-)circular prior on the projected skew.

4 Strong lens analysis of Abell 3827

We use the observed positions of lensed multiple images to fit a mass model of the cluster. Our choice of model parameters and their priors is based on those of M15, with some additional degrees of freedom. We assume $0.8\arcsec$ uncertainty on the position of lensed image Ao8, which has only been detected from the ground. We assume $0.2\arcsec$ uncertainty on the position of all other lensed images, which are identified by HST.

4.1 Fiducial mass model

The cluster’s large-scale mass distribution is modelled as a single PIEMD. Based on a comprehensive (but slow) initial exploration of parameter space, its position is given by a broad Gaussian prior with $\sigma=2\arcsec=3.66\,\text{kpc}$, centred on the position of galaxy N2. Flat priors are imposed on its ellipticity ($\epsilon<0.75$), core size ($r_{\mathrm{core}}<40\arcsec$) and velocity dispersion (300<$$\sigma_{\mathrm{v}}$$<1000 km/s). Its cut radius is fixed at $r_{\mathrm{cut}}=1000\arcsec$, well outside the strong lensing region, i.e. away from any multiple image constraints.

Central galaxies N1–N4 are each modelled as a stellar component (which was not included in the fiducial model of M15), plus a dark matter one. Following Giocoli et al. (2012), the stellar components are modelled with Hernquist (1990) profiles:

[TABLE]

where the scale radius $r_{s}$ is related to the half mass radius $R_{e}$, such that $R_{e}=r_{s}/0.551$, and the scale density $\rho_{s}=M_{\mathrm{total}}/\left(2\pi r_{s}^{3}\right)$. We fix the mass of the stellar component, and its half-mass radius, using the optical magnitudes and profiles measured by M15. These parameters are listed in Table 1.

The four central galaxies’ dark matter components are now modelled as PISPs. We impose flat priors on their positions, in $4\arcsec\times 4\arcsec$ boxes centred on their luminosity peaks, plus flat priors on their ellipticity ($\epsilon<0.5$) and velocity dispersion ($v_{\mathrm{disp}}<600$ km/s). We fix $r_{\mathrm{cut}}=40\arcsec=73$ kpc (Limousin et al., 2007a).

Galaxy N6 is much fainter than the others, so we approximate its total mass distribution as a single PIEMD. This has a fixed position and ellipticity to match the light distribution, and only its velocity dispersion is optimised (with a flat prior $v_{\mathrm{disp}}<500$ km/s).

We optimise the free parameters using Lenstool, with runmode=3. This runmode is used to fully explore the posterior (Jullo et al., 2007). (Modal) maximum likelihood parameters are shown in Table 1, and the corresponding mass distribution is shown in Figure 3. The best fit model achieves a $RMS$ offset between the observed and predicted positions of multiple images of $\left<\mathrm{rms}\right>_{i}$=$0.26\arcsec$. There are $54$ constraints and $35$ free parameters in our model. The modal $\chi^{2}/\text{dof}$=$67.1/19$ with a log likelihood of 8.18. The full posterior probability distribution for the dark matter associated with galaxies N1–N4 is shown in Figures 4 and 5.

4.2 Sensitivity to model choices

4.2.1 Stellar mass components

Galaxies definitely contain stars, and those stars have mass. Not accounting for this mass could bias the skew measurement. In an offset DM halo, the stellar mass will lead the DM peak, while any tail of DM particles will trail behind in the opposite direction. Not accounting for this stellar mass could weaken the skew measurement, or in the worst case scenario, if the stellar mass is greater than that of scattered particles, reverse the direction of the measured skew. We have explicitly modelled the stellar mass seperately to avoid any bias in the inferred skew. In practice, as in M15, we find that including the stellar mass component (or even multiplying/dividing its mass by a factor 2) does not change any other results, within their statistical errors.

4.2.2 Identification of new lensed images

Adding constraints from the two new lensed images Ao7 and Ao8 tightens constraints on nearby galaxies N3 and N4. These (demagnified) images are unresolved, and any of the features in the background spiral could be assigned to them. We have tried relabelling one or both of the demagnified images as either the bulge, Ao$n$, or one of the two brightest knots of star formation, Aa$n$ or Ab$n$. Lenstool’s outputs are statistically indistinguishable. In all cases, the entire background spiral galaxy is predicted to be lensed onto both the northern and the southern demagnified images.

4.2.3 Mass in other cluster member galaxies

We also tested the impact of adding more cluster member galaxies to the mass model. These galaxies were identified using a colour-magnitude selection using the F814W and F606W2 HST/ACS band imaging. Source detection was done using sextractor (Bertin & Arnouts, 1996) in dual mode, with reference taken in the F814W-band. We then identified as cluster members all galaxies brighter than $mag_{\mathrm{F814W}}<23$ and within 1$\sigma$ of the red sequence best-fit:

[TABLE]

Our final cluster member catalogue contains 147 galaxies.

These galaxies are added to the mass model as small scale perturbers. We assume fixed cut radius and velocity dispersion, scaled by their luminosities in the F814W-band. This methodology has been successfully validated by Harvey et al. (2016), and adopted widely in previous work (e.g. Jauzac et al., 2014; Richard et al., 2014; Limousin et al., 2007b). To derive $L^{\ast}$, we use the $K^{\ast}$ magnitudes obtained by Lin et al. (2006) as a funciton of cluster redshift. lenstool is then scaling the cut radius and velocity dispersion of each galaxies in our catalogue relative to a $K^{\ast}=16.6$ galaxy with velocity dispersion $\sigma^{\ast}=108.4\pm 27.5$ km/s, and cut radius $r_{\mathrm{cut}}^{\ast}=48.5\pm 16.0$ kpc.

Including all cluster memeber galaxies in a re-optimised mass model significantly affects neither offset nor skewness measurements of dark matter associated with central galaxies N1–N4. By far the most affected measurement is the skewness of galaxy N1, which increases to $s=0.28^{+0.01}_{-0.31}$. All other quantities remain consistent within random noise.

5 Discussion

5.1 Galaxy N1

The previously-detected offset between galaxy N1’s stars and dark matter persists at $>3\sigma$ in our new analysis. Adding two free parameters for the asymmetry of its dark matter slightly increases uncertainty in its position. The modal offset is $\left(\Delta x,\Delta y\right)=\left(-0.22^{+0.25}_{-0.14},-0.81^{+0.16}_{-0.17}\right)$ for an unskewed model, and $\left(\Delta x,\Delta y\right)=\left(-0.29^{+0.25}_{-0.14},-0.71^{+0.30}_{-0.16}\right)$ if skewness is allowed. The consistency between these suggests that the measured position of the density peak is robust against the skew parameter probed here.

If the offset is entirely due to an effective drag force through frequent dark matter self-interactions, it implies a momentum-transfer interaction cross-section $\tilde{\sigma}/{m_{\mathrm{DM}}}\gtrsim 1$ cm2g*-1*, assuming galaxy N1 is falling into the cluster for the first time (Kahlhoefer et al., 2015). In general, We agree with this interpretation444We have also repeated Kahlhoefer et al. (2015)’s calculation of $\tilde{\sigma}/{m_{\mathrm{DM}}}$ but integrating the effect of the restoring force on the entire distribution rather than just the peak. The difference is not significant., but note that the cross-section can be lower if the galaxy has completed multiple orbits; its current direction of motion is unknown.

We also find (at much lower statistical significance) that galaxy N1’s dark matter is skewed in a direction consistent with the SIDM interpretation of its offset. This could be signs of a tail of scattered DM particles and would favour high momentum transfer scattering. However, the weak statistical significance of our result makes it impossible to rule out the low momentum transfer case. Figure 6 shows the posterior PDF of the skew onto the vector pointing from the DM peak to the stellar luminosity in the fiducial model, such that a positive skew corresponds to the direction predicted by SIDM. The peak of the posterior and $1\sigma$ errors are $s=0.23^{+0.05}_{-0.22}$. If we individually project the skewness onto the offset direction individually in all MCMC samples, we find $s=0.26^{+0.03}_{-0.22}$. Finally in our model that contained the additional 147 cluster galaxies we find $s=0.28^{+0.01}_{-0.31}.$ In all cases $\sim 70\%$ of the posterior PDF lies at $s>0$. As the posterior peak is near the edge of the prior, which is chosen coservatively (see §3.2), a parametric halo model that does not break down for large skew parameters could result in a stronger detection.

5.2 Galaxies N2 and N3

The dark matter associated with galaxies N2 and N3 appears symmetric, and coincident with the stars. This result does not preclude offsets from existing along the line of sight. Furthermore, even if the galaxies are moving in the plane of the sky, they could be behind or in front of the most dense regions of the cluster core, and therefore passing through a less dense medium, experiencing less drag.

5.3 Galaxy N4

The measured position of N4’s dark matter is intriguing. Accounting purely for statistical error bars, thanks to the confirmed positions of demagnified images, we find that galaxy N4 is offset from the galaxy’s stars at the $3\sigma$ level. However, the offset position is mildly degenerate with the position of the cluster-scale dark matter (Figure 5), thus a flat prior on the cluster scale halo could lead to a different offset measurement. Furthermore, the measured ellipticity of the galaxy light is contaminated by light from an adjacent Milky Way star, and its position may also be.

Despite the measured offset of dark matter from galaxy N4, it shows no sign of skewness. If the offset is spurious, as discussed above, then a tail is not expected. If the offset is real, the lack of skew favours low momentum transfer scattering and (at very low S/N) is in mild tension with the skew detected in galaxy N1. Unknown systematics in the modelling of DM around either galaxy could be responsible. Nonetheless, there is also a possible physical explanation for this discrepancy. Galaxy N4 is in a higher density environment than galaxy N1, closer to the cluster core. It is possible that any tail of scatttered N4 particles has been tidally destroyed by the steeper gradient in gravitational potential.

6 Conclusion

We have developed a parametric lens models for asymmetrically skewed mass distributions. This can be used to search for scattered (self-interacting) dark matter in colliding systems. More generally, it will also be useful to investigate claims of dynamically-induced asymmetry (Prasad & Jog, 2016; Chemin et al., 2016), or tidal tails (which are asymmetric if the size of a body is large compared to its distance from the centre of potential).

We have also presented a new model for the distribution of mass in galaxy cluster Abell 3827. Our choice of flat priors for the position of all galaxies’ dark matter leads to a detected offset between a second galaxy’s stars and its dark matter. New VLT/MUSE observations tighten the constraints on that offset. Neither measured offset changes significantly if the models are allowed extra freedom to become skewed.

We find tantalising, but low significance evidence that the galaxy closest to multiply lensed images (and therefore the best constrained) has an asymmetric distribution of dark matter, skewed in the same direction as its offset from stars. We emphasise that our skew model, which captures the qualitative behaviour of scattered DM particles, is primarily motivated by mathematical convenience and that all skew measurements here are model dependent. More work will be needed to determine the significance of this result: whether it is physical, or an artefact of systematics in parametric lens modelling. Even in mock data where the true skew is known, skewness cannot be measured to high precision in a system as complex as Abell 3827. This is probably because the effect of skewness on lensed image positions is smaller than the effects of other free parameters.

A promising direction for future investigation may be provided by pairs of field galaxy in the SLACS survey, one of which has already been found to have an offset between dark and luminous matter (Shu et al., 2016). Whilst the SIDM model predicts a largest (most easily observable) offset in galaxies moving through a dense intracluster medium, it may be possible to more tightly constrain any asymmetry of dark matter in these simpler systems. If the directions of their dark matter tails correlate with the directions of their offsets, this evidence would support the hypothesis of SIDM with a high-mass mediator particle.

Acknowledgements

The authors are grateful for constructive conversations with Benjamin Clément, Eric Jullo, Felix Kahlhoefer, Thomas Kitching, Marceau Limousin and Subir Sarkar. We would also like to thank the anonymous referee whose comments have significantly improved the paper. Figure 4 was produced using the python module corner (Foreman-Mackey, 2016). PT is supported by the STFC. RM is supported by the Royal Society. AR and MJ are supported by the UK Science and Technology Facilities Council (grant numbers ST/H005234/1 and ST/L00075X/1).

Facilities: This paper uses data from observations GO-12817 (PI: R. Massey) with the NASA/ESA Hubble Space Telescope, obtained at the Space Telescope Science Institute, which is operated by AURA Inc, under NASA contract NAS 5-26555. This paper also uses data from observations made with ESO Telescopes at the La Silla Paranal Observatory under programmes 294.A-5014 and 295.A-5018 (PI: R. Massey). We thank the Paranal Science Operations team for carrying out those observations. Our Lenstool optimisation used the DiRAC Data Centric system at Durham University, operated by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). This equipment was funded by BIS National E-infrastructure capital grant ST/K00042X/1, STFC capital grant ST/H008519/1, STFC DiRAC Operations grant ST/K003267/1 and Durham University. DiRAC is part of the UK National e-Infrastructure.

Appendix A Alternative method to generalise lens mass distributions

Another way to introduce asymmetry is to apply a weighting function $w\left(\mathbb{r};\{a_{i}\}\right)$ to an elliptical lensing potential

[TABLE]

where $\{a_{i}\}$ are a set of parameters. The deflection and surface mass density are readily computed by differentiating.

We consider weighting functions of the form

[TABLE]

where $s$ is a skew parameter and $f(r,\theta)$ is written in polar coordinates. The second term acts as a perturbation away from elliptical symmetry of $O(s)$, with $s=0$ corresponding to the elliptically symmetric case. We chose $f(r,\theta)$ to meet the following criteria:

• •

To ensure that the space about $s=0$ is explored symmetrically in Lenstool, so that a non-zero skew is not artificially recovered, we require that $sf(r,\theta)=-sf(r,\theta+\pi)$.

• •

To avoid difficulties near the origin, we require $f(r,\theta)$ to be an increasing function of $r$. This is also physically motivated, as it is difficult to scatter particles from the centre of the potential well at $r=0$.

• •

To ensure that the surface mass density remains positive (or becomes negative only for large $r$ well outside any region of interest), we require $f\left(r,\theta\right)$ to be bounded.

A.1 The Weighting Function Pseudo-Isothermal Skewed Potential

Meeting the above conditions we form the weighted Pseudo-Isothermal Skewed Potential (wPISP) by applying the weighting function

[TABLE]

to the PIEMD potential, where $r_{s}$ is a new scale radius, $\beta$ sets the radial dependence of the skew, and $\phi_{s}$ is the skew angle.555The inverse tangent form of the radial dependence is not physically motivated, and other functional forms may also work. While it is mathematically unnecessary to have two degenerate parameters $r_{s}$ and $\beta$, to avoid computational divisions by zero, the distributed Lenstool implementation uses hardcoded $r_{s}=0.1\arcsec$ and allows $\beta$ to vary. This is now available as potential 812 in Lenstool.

The resulting surface mass densities are shown in Figure 7. The qualitative shape of the isodensity contours changes inside or outside the scale radius (owing to the sign change in second derivative of the inverse tangent). This feature could be used to isolate a behaviour that best matches numerical simulations, by fixing very large or very small $r_{s}$, or to capture more complex behaviour.

The total mass of a wPISP is identical to that of of a PIEMD. Since the weighting function is normalised by construction, the integrated mass density of a PIEMD and wPISP over a circular region are the same:

[TABLE]

where $\theta$-dependence cancels. Taking the limit as $\mid\mathbb{r}\mid\rightarrow\infty$, the left hand side will converge to the total mass of a PIEMD with ellipticity $\epsilon$, and the right to the mass of an equivalent wPISP. However, the position of the density peak varies slightly as a function of $s$. Care would need to be taken if using a wPISP to measure offsets of dark matter.

As with the PISP, this model breaks down for large values of $s$, since the weighting function has been applied to the potential and not the density. We have found that the value of $s$ where this occurs is sensitive to $\beta$ and the cut and core radii. For this reason, we recommend testing the boundaries of the parameter space for a breakdown of the desired skewed behaviour before substantial future work. Nonetheless, we tested the wPISP against the null hypothesis of the unskewed example_with_images system distributed with Lenstool (see §3.3). Fixing $\beta=0.01$ and starting with a flat prior $s\in[0.3,0.3]$, Lenstool successfully recovers skewness $s=0.002^{+0.002}_{-0.002}$.

A.2 Pseudo Isothermal Varying Ellipticity Mass Distribution

Despali et al. (2016) predict that, even with standard cold dark matter, the ellipticity of a cluster scale halo should change as a function of radius, becoming more elongated further from the centre. This prediction can be tested by using the weighting function formalism to add an extra parameter to halo models that mimics this behaviour. To achieve this, we combine a weighted sum of two different mass densities with different ellipticities into a Pseudo Isothermal Varying Ellipticity Mass Distribution (PIVEMD)

We suggest a mass density of the following form:

[TABLE]

where $\Sigma_{\epsilon_{1}}\left(\mathbb{r}\right)$ and $\Sigma_{\epsilon_{2}}\left(\mathbb{r}\right)$ are two elliptical profiles with ellipticity $\epsilon_{1}$ and $\epsilon_{2}$. All the other parameters for these two densities should be shared. In this case, we find it most effective (and possible) to apply the weighting function directly to the mass density. To be computationally efficient within an MCMC loop, deflection angles must be computed once, using numerical integration, and stored in a look-up table.

The weighting functions $w_{i}$ should meet the following criteria:

• •

To normalise the total mass, $w_{1}\left(\mathbb{r}\right)+w_{2}\left(\mathbb{r}\right)\equiv 1$, $\forall r$.

• •

So that one ellipticity dominates at small $r$ and the other at large $r$, let $w_{1}\left(\mathbb{r}\right)\rightarrow 1\text{ as }r\rightarrow\infty$, $w_{1}\left(\mathbb{r}\right)\rightarrow 0\text{ as }r\rightarrow 0$, $w_{2}\left(\mathbb{r}\right)\rightarrow 0\text{ as }r\rightarrow\infty$ and $w_{2}\left(\mathbb{r}\right)\rightarrow 1\text{ as }r\rightarrow 0$.

Although this is quite a general set of conditions, we can take, for example

[TABLE]

where $\beta$ controls the radial dependence. The resulting mass distribution for this weighted sum is illustrated in Figure 8.