# The Origin of Lopsided Satellite Galaxy Distribution in Galaxy Pairs

Chen Chris Gong1,2, Noam I Libeskind2,3, Elmo Tempel4,2,Quan Guo5,2, Stefan
Gottlöber2, Gustavo Yepes6, Peng Wang2, Jenny Sorce7,2, Marcel Pawlowski2
1University of Potsdam, Institute of Physics and Astronomy, Karl-Liebknecht-
Straße 32, 14476 Potsdam
2Leibniz-Institut für Astrophysik Potsdam (AIP), An der Sternwarte 16,
D-14482, Potsdam, Germany
3University of Lyon; UCB Lyon 1/CNRS/IN2P3; IPN Lyon, France
4Tartu Observatory, University of Tartu, Observatooriumi 1, 61602
T$\tilde{o}$ravere, Estonia
5Shanghai Astronomical Observatory (SHAO), Nandan Road 80, Shanghai 200030,
China
6 Departamento de Física Teórica and Centro de Investigación Avanzada en
Física Fundamental (CIAFF), Facultad de Ciencias,
Universidad Autónoma de Madrid, 28049 Madrid, Spain
7 University of Lyon, University of Lyon 1, Ens de Lyon, CNRS, Centre de
Recherche Astrophysique de Lyon UMR5574, F-69230, Saint-Genis-Laval, France
E-mail<EMAIL_ADDRESS>
###### Abstract

It is well known that satellite galaxies are not isotropically distributed
among their host galaxies as suggested by most interpretations of the
$\Lambda$CDM model. One type of anisotropy recently detected in the SDSS (and
seen when examining the distribution of satellites in the Local Group and in
the Centaurus group) is a tendency to be so-called “lopsided”. Namely, in
pairs of galaxies (like Andromeda and the Milky Way) the satellites are more
likely to inhabit the region in between the pair, rather than on opposing
sides. Although recent studies found a similar set up when comparing pairs of
galaxies in $\Lambda$CDM simulations indicating that such a set up is not
inconsistent with $\Lambda$CDM, the origin has yet to be explained. Here we
examine the origin of such lopsided setups by first identifying such
distributions in pairs of galaxies in numerical cosmological simulations, and
then tracking back the orbital trajectories of satellites (which at $z=0$
display the effect). We report two main results: first, the lopsided
distribution was stronger in the past and weakens towards $z=0$. Second, the
weakening of the signal is due to the interaction of satellite galaxies with
the pair. Finally, we show that the $z=0$ signal is driven primarily by
satellites that are on first approach, who have yet to experience a “flyby”.
This suggests that the signal seen in the observations is also dominated by
dynamically young accretion events.

###### keywords:

Local Group – galaxies: evolution – galaxies: formation – Galaxy: kinematics
and dynamics – cosmology: theory – dark matter

††pubyear: 2018

##  1 Introduction

The favoured scenario for how structures in the Universe form is known as the
$\Lambda$CDM paradigm, wherein the universe is composed primarily of dark
matter ($\sim 26\%$), dark energy ($\sim 70\%$) and baryons ($\sim 4\%$)
(Planck Collaboration et al., 2014). Accordingly, small perturbations in the
otherwise nearly homogeneous density field, decouple at high red-shift from
the expansion of the Universe and collapse to form the first dark matter
haloes. Dark matter haloes proceed to grow hierarchically via gravitational
instability (Zeldovich, 1980) in a “bottom-up” fashion: small haloes merge to
form larger ones and so on (White & Rees, 1978).

One of the defining features of this scenario is that the dark matter haloes
are continually merging with one another. Because of the nature of the power
spectrum of fluctuations measured in the CMB, the majority of the merger
events will involve small, low mass haloes (Lacey & Cole, 1993). Indeed
studies such as Genel et al. (2010) have shown that Milky Way mass haloes gain
just $20-30$% of their mass via discrete resolved mergers (with a mass ratio
of less than 10:1) and 70-80% via smooth ambient accretion (but see also:
Fakhouri & Ma, 2008; Stewart et al., 2008; Madau et al., 2008; Angulo & White,
2010; L’Huillier et al., 2012). Often the dynamical time needed for such small
haloes to lose enough angular momentum so that tidal forces can rip them apart
them is quite long (e.g. equation 2 in Wetzel & White, 2010). As a result
small “satellite” galaxies exist in abundance in the nearby environment of
larger galaxies (Tollerud et al., 2008; Koposov et al., 2008; Newton et al.,
2018).

Such dwarf galaxies are observed to cluster in the vicinity of larger galaxies
Pawlowski et al. (2012); Muñoz et al. (2015); Ordenes-Briceño et al. (2018);
Taylor et al. (2018). The geometry of this clustering has been a subject of
investigation at least since Holmberg (1969) identified an anisotropy in the
distribution of satellites surrounding disk galaxies. Indeed, there is an open
debate in the literature on if satellite galaxies are distributed
isotropically based on observation (e.g. Kroupa et al., 2005; Ibata et al.,
2013; Ibata et al., 2014b; Lüghausen et al., 2014; Pawlowski, 2018; Tully et
al., 2015; Cautun et al., 2015b; Müller et al., 2016; Bullock & Boylan-
Kolchin, 2017; Arakelyan et al., 2018, among others); and if the satellite
galaxies are distributed anisotropically, whether or not there is a persistent
sense of orientation or rotation among them (Metz et al., 2008; Pawlowski et
al., 2013; Ibata et al., 2014a; Müller et al., 2018); and when there is,
whether such a structure could maintain itself (Bowden et al. 2014). Most of
these studies have explicitly focused on the anisotropy of satellite
distributions of field galaxies (as opposed to galaxies in clusters or pair of
galaxies). Yet recently, the anisotropic distribution of the satellites of
galaxies pairs has been brought to attention, which has been seen in the Local
Group as well as in the Centaurus group (Müller et al. (2015), Müller et al.
(2017)).

The existence of galaxy pairs (such as the two main galaxies of the Local
Group and those of the Centaurus group) presents the cosmologist with a unique
testing ground for questions related to structure formation and the
cosmological model. This has been recognised at least since the timing
argument of Kahn & Woltjer (1959) who famously used the dynamics of the Local
Group to compute its mass according to a cosmological paradigm (see also
extensions such as Partridge et al., 2013). Indeed in recent years there has
been an abundance of work that examines in detail the properties of pairs in
both large sky surveys (e.g. Ellison et al., 2008; Shen et al., 2016; Behroozi
et al., 2015; Tempel & Tamm, 2015; Gordon et al., 2017, to name a few) or in
numerical simulations (e.g. Moreno et al., 2013; Perez et al., 2009; Kelley et
al., 2018, and references therein). The vast majority of these studies have
focused on pairs for the purpose of understanding the physics of galaxy
formation and galaxy conformity (Kauffmann et al., 2013), namely how living in
a pair (and the eventual merger) affects a galaxy’s star formation efficiency,
gas reservoir and morphology.

A number of recent projects have extended the research on satellite galaxy
distribution from those around isolated field galaxies towards those which
inhabit galaxy pairs, and have compared the distribution of small satellites
in these two types of environment when their hosts are of similar mass and
other characteristics. Since only a handful of such pairs have been both
spectroscopically confirmed and imaged with sensitivity down to the smallest
dwarfs (for example the Centaurus A - M83 pair), most people wishing to
examine satellites of galaxy pairs must turn to large photometric catalogues
(e.g. Wang et al., 2019, among many others) and employ stacking.

Libeskind et al. (2016) examined the distribution of satellite galaxies near
pairs of galaxies that resemble the “Local Group” (in a very broad sense)
identified in the SDSS. Since each pair in such a survey has just a few
satellites, any trend in the distribution of satellites is revealed by
stacking tens of thousands of such pairs. They recovered the well known result
that satellites are distributed anisotropically or ellipsoidally (and not
spherically) among their hosts (see also Sales & Lambas, 2004; Libeskind et
al., 2005; Agustsson & Brainerd, 2010; Libeskind et al., 2015a, b; Cautun et
al., 2015a). The analysis further revealed that satellites of galaxy pairs
identified in the SDSS show a statistically significant tendency to bulge
towards the partner galaxy, in a lopsided manner. This effect was coined “the
Lopsided distribution of satellite galaxies” owing to its ovular or egg-like
configuration. While the small satellites are found to preferentially inhabit
the region between the pair, this striking effect was shown not to be a result
of projection effects or simply the overlapping of two independent satellite
distributions. It should be noted that the satellites of Andromeda display
exactly the same lopsidedness with the majority on the near side of M31.
(Bowden et al. (2014) addressed this issue analytically ruling out tides as
being responsible, but suggesting that recent group accretion may explain the
effect.) Interestingly after the publication of Libeskind et al. (2016), Epps
& Hudson (2017) found a very similar signal by stacking lensing images of
pairs, albeit on slightly larger scales (see also De Graaff et al., 2017).

Since then, Pawlowski et al. (2017) has found such a lopsided effect of
satellite galaxy distribution in different $\Lambda$CDM simulations confirming
that such an arrangement is a feature of the cosmological paradigm, possibly
even a generic one. Although the strength of the signal found in simulations
is not identical to that found in observations, the discrepancy can easily be
attributed to the slight differences in the samples being compared, and to
fore- and background contaminations in the observed sample. Pawlowski et al.
(2017) noted that although successful in finding such a lopsided distribution,
the evolutionary history and causal mechanism for the lopsidedness was left
unexplained.

To some degree, in this study, we pick up where Pawlowski et al. (2017) left
off. We conduct an investigation into the origin of the lopsided satellite
distribution in galaxy pairs like the Local Group which are found in the
cosmological simulations. Due to the nature of our study we are then able to
track these lopsided distributions back in time and investigate their origin.
In doing so we characterize the accretion history of the satellites and their
host haloes.

##  2 Methods

In this section we describe how our sample of halo pairs (and their
satellites) are identified.

###  2.1 Simulation and the identification of haloes and subhaloes

For our analysis, we use the “Extremely Small Multidark” simulation (ESMD).
ESMD is an extension of the Multidark suite of pure $N$-body simulations111
see https://www.cosmosim.org for more detailed information. (Klypin et al.,
2016). This is a dark-matter-only $N$-body simulation assuming a $\Lambda$CDM
power spectrum of fluctuations according to the Planck cosmological parameters
(Rodríguez-Puebla et al. 2016; $\Omega_{\Lambda}=0.69$, $\Omega_{M}=0.33$,
$\sigma_{8}=0.83$ and $H_{0}=67.77$ km s-1 Mpc-1). The simulation consists of
$N_{\rm part}=2048^{3}$ particles in a periodic box of side length $L_{\rm
box}=64$ $h^{-1}$Mpc. Such a simulation achieves a mass resolution of
$2.6\times 10^{6}~{}h^{-1}$M⊙ per particle and a spatial softening length of
$1h^{-1}$kpc.

The haloes and subhaloes in the simulation are identified using the publicly
available Rockstar halo finder (Robust Over-density Calculation using K-Space
Topologically Adaptive Refinement) which is described in detail in Behroozi et
al. (2012). We briefly summarise how Rockstar works but refer the reader to
Behroozi et al. (2012) for details. Rockstar is a massively parallel method of
identifying haloes and their substructures based on adaptive hierarchical
refinement of friends-of-friends (FOF) groups in six phase-space dimensions
and one time dimension.

The halo finder begins by identifying 3DFOF groups based only on positional
information. Then with these 3DFOF groups, Rockstar adaptively chooses a
phase-space linking length based on the standard deviations of the particle
distribution in position and velocity space. In this way subhaloes can be
identified and are then assigned to their hosts similar to the way particles
are assigned to their 6DFOF subgroups. Additionally, particle-based merger
trees are created. For a given halo, its descendant is assigned as the halo in
the next time step which share the maximum number of common particles. We
consider only haloes resolved with more than 20 particles, resulting in a
minimum halo mass of $M>5.2\times 10^{7}h^{-1}M_{\odot}$. Although uncertain
due to the stochasticity of star formation processes at these low masses, such
dwarfs could host galaxies with a stellar mass below $10^{4}$ (Cautun et al.
2014; Reed et al. 2017) or could have star formation entirely suppressed
(Sales et al. 2017).

Figure 1: A sketch of a typical binary halo system. The pair is separated by a
distance $d_{\text{sep}}$. Haloes or subhaloes that are within 1/2 the
separation distance, i.e. $d_{\rm sat}<0.5d_{\text{sep}}$ from a host are
identified and referred to as satellites of that host (colored red and blue
here). For illustration, the position angle $\theta$ of a satellite (in
magenta) of its host (in blue) is measured between the line connecting the
pair and a satellite’s position vector. $\theta_{\text{rv}}$ indicates the
angle between the velocity of the satellite ($\vec{v}_{\rm sat}$) and its
position angle and is a measure of how radial the infall is. The bigger circle
indicates the exclusion range (see criteria 3), beyond which a third halo
which is more massive than the less massive member (red) of the pair cannot
trespass (e.g. the halo in green).

###  2.2 Identification of halo pairs and satellites

Loosely based on the mass and geometrical set up of the Local Group, we
identify halo pairs at $z=0$ using the following criteria:

  1. 1.

halo mass range: $10^{11}<M/(h^{-1}M_{\odot})<5\times 10^{12}$

  2. 2.

halo separation: $0.3<d_{\text{sep}}/(h^{-1}{\rm Mpc})<1.5$

  3. 3.

Exclusion: There is no third halo more massive than 1/2 of the mass of the
less massive host halo within a radius, centered on the binary’s midpoint, and
equal to the separation distance from the binary center.

We note that the values for the halo mass, halo pair separation and exclusion
mass are allowed to be within chosen ranges, which are loosely motivated by
the values seen in the Local Group. They are deliberately chosen, where we try
to both maintain a physical justification (approximating, even loosely, the
situation in the Local Group) and to obtain a maximal sample size. At $z=0$,
our criteria returns $2~{}252$ halo pairs with $625~{}033$ satellites among
them. The geometry of satellites in such a halo pair is shown in Fig. 1.

For each halo pair, the smaller haloes and subhaloes around them are
identified as satellite haloes within $0.5d_{\text{sep}}$ from its nearest
host halo. These include both subhaloes located within their host’s virial
radius and “field dwarfs" located outside the virial radius; we collectively
term all of these “satellite haloes”. Satellite position may be characterised
by the angle, $\theta$, made between the satellites position vector from the
closest primary and the line connecting the binary. The viewing angle is
rotated such that $\theta$ is maximum, and the satellite – host – partner
triangle is viewed “face on”. We treat all satellites equally regardless of
mass. In effect our analysis proceeds by considering the dynamics of each
satellite halo as it is accreted by its host halo pair.

###  2.3 Control (Overlapping) sample at $z=0$ constructed from two isolated
haloes

Following similar reasoning as in Libeskind et al. (2016), when the
distribution of $\theta$ is examined, any non-uniformity must be separated
from the possible effect of overlap. Namely, two haloes with spherically
isotropic satellite distributions, when brought near to each other, would
naturally show an increased number density of satellites in regions where
their satellite distributions overlap, which happens purely geometrically,
without the influence from the force of gravity. This follows works such as
Prada et al. (2006) who have shown that haloes of the mass range we are
considering are extended up to $2\sim 3$ virial radius.

In order to quantify the effect of overlap, halo pairs are artificially
constructed from isolated haloes and their satellites. To do this, we must
match each member of a pair to an isolated halo. The two isolated haloes and
their satellites are then artificially brought together and placed at the same
separation distance as the real halo pair. In detail, for each halo pair, we
note the pair separation $d_{\rm sep}$. We then identify two haloes in the
same mass range (i.e. $10^{11}-5\times 10^{12}h^{-1}M_{\odot}$) that are
isolated and are the closest in mass to the members of the real halo pair,
respectively. Here, “isolated” means that the nearest halo, which is either
more massive than the less massive member of the pair, or more massive than
the “isolated” halo itself, is at least $1.5\times d_{\rm sep}$ away. This
ensures that when we bring the two haloes together and place them at distance
$d_{\rm sep}$ apart, there is no third halo bigger than either member of the
real halo pair and no third halo bigger than either member of the artificial
halo pair. Once two isolated haloes are matched to a real pair, their
satellites are chosen in the same way as in the main sample, namely all haloes
within $d_{\rm sep}/2$ of each halo’s center are identified. The isolated
haloes (and their satellites) are then placed at a distance of $d_{\rm sep}$
from each other, therefore creating an artificial pair composed of two
isolated haloes.

Figure 2: The lopsided signal is characterised by the probability distribution
of $\cos\theta$ shown in top panel. In bottom panel, in order to emphasise the
lopsided signal we divide the distribution into two regions: the “facing”
(i.e. $\cos\theta>0$, solid) and “facing away” region from the other primary
halo (i.e. $\cos\theta<0$, dashed), and plot the two distributions of absolute
value of $\cos\theta$ for direct comparison. It is obvious that the solid line
dominates for $|\cos\theta|>0.8$, indicating there are more satellites close
to the line connecting the pair in the region in between the pair than 180
degrees away. The thin vertical lines denote the medians of the distributions
which are although close, are statistically inconsistent with a symmetric
(non-lopsided) distribution. The lopsidedness metric $q$ and its statistical
significance $\sigma$ are defined at the end of section 3.1. In short, $q$ is
the ratio between the number of satellites in the region with $\cos>0.8$ and
the number of those in $\cos<0.8$. Its $\sigma$ is computed from the variance
in this quantity from random trials.

##  3 Results

###  3.1 Lopsided $z=0$ signal

We begin by re-examining the findings by Pawlowski et al. (2017), who confirm
the observational result of Libeskind et al. (2016) with similar numerical
simulations to those used here. In Fig. 2 we show the normalized probability
distribution of $\cos\theta$ for all satellites of all host pairs in our
sample. The peaks of the distribution at $|\cos\theta|\approx 1$ in Fig. 2 top
panel indicates that the satellites have a strong tendency to be aligned with
the line connecting the halo pair in both regions.

A direct comparison between the two regions “facing” the other primary halo
and “facing away” from the other primary halo show lopsidedness exists for
this signal. An excess probability is seen for $|\cos\theta|\mathrel{\lower
2.58334pt\hbox{$\buildrel{\textstyle>}\over{\scriptstyle\sim}$}}0.8$ for the
$\cos\theta>0$ interval with respect to $\cos\theta<0$. Furthermore the
asymmetry can be quantified by the median of the absolute value:

$\displaystyle\langle|\cos\theta|\rangle=0.529$ $\displaystyle{\rm
for}\cos\theta>0$ $\displaystyle\langle|\cos\theta|\rangle=0.523$
$\displaystyle{\rm for}\cos\theta<0$

Although this does not seem like a large difference, given our sample size it
is strongly statistically significant, as described below.

To consider the statistical significance of the median values stated above, we
compare to 10,000 samples of the same size. Namely, we generate $625~{}033$
random, uniformly distributed numbers in the interval [0,1] (representing a
uniform $|\cos\theta|$ distribution). The median of these $625~{}033$ is then
computed. This is then repeated 10,000 times. In this way we compute the
distribution (standard deviation) of the median value of $625~{}033$ random
numbers drawn from [0,1] as $\sigma_{\rm m}=0.001.$ Therefore both the facing
($\cos\theta>0$) and facing away ($\cos\theta<0$) samples are inconsistent
with a spherical uniform distribution at roughly $\sim 10\sigma_{\rm m}$
level.

We now use this value of $\sigma_{\rm m}$ to determine how much more
asymmetric the “facing” signal is than the “facing away” signal, determining
that the facing sample is more aligned than the facing away sample at roughly
the $\sim 2.1\sigma_{\rm m}$ level. This is, in effect, one aspect of the
lopsided signal found by Libeskind et al. (2016): satellites facing the
opposite halo are more aligned with the line connecting the pair than
satellites facing away from the opposite halo. The fact that the satellites in
the “facing region” are just 2.1$\sigma_{\rm m}$ (and not more as found by
Libeskind et al., 2016), than their counterparts in the “facing away” region,
is surprising and demonstrates that both regions include flattened satellite
distributions, aligned with the geometry of the pair.

We can also define a lopsidedness metric $q$ that is used in the rest of the
paper to quantify the lopsidedness. For a distribution of $\cos\theta$ of size
$N$, the number of satellite haloes in the most aligned subsample, namely with
$0.8<\cos\theta<1.0$ is termed $N_{+}$. The number of satellite haloes in the
most counter-aligned subsample, namely $-1.0<\cos\theta<-0.8$ is termed
$N_{-}$. We then define $q=N_{+}/N_{-}$: when $q>1$ we have a bulging
satellite distribution facing the other pair member. When $q<1$ we have a
bulging distribution facing away from the other pair member. In order to asses
the statistical significance of a given value of $q$, we compute the expected
fluctuations in $q$ by generating 10 000 random data sets of the same size,
drawn from uniform distribution between $-1$ and $1$, and getting the expected
value of $q$ from these random sets: $\sigma=\operatorname{std}(q)$.

Figure 3: We quantify the lopsided signal by showing it as both a probability
distribution of $\cos\theta$ (x axis and top histograms) as well as a function
of satellite radial position (bottom 2D histogram). Satellite distance to the
host is normalized by the pair separation, $d_{\text{sep}}$ (y axis). We show:
(a) the signal and (b) control sample of artificial pairs constructed from
isolated haloes.

###  3.2 Quantifying the lopsided signal at $z=0$.

Figure 4: Probability distribution of $\cos\theta$ plotted for 4 subsamples of
satellites sorted by their distance to hosts. Although the satellite
distributions at all distances are anisotropic and flattened, it is only at
the greatest distances that they are lopsided.

In this section we wish both to study what features of the halo pair might
correlate with the lopsided signal and we wish to compare it with the control
sample of artificial pairs constructed from isolated haloes (see subsection
2.3). In Fig. 3 the cosine angle and radial distribution of satellites are
shown in normalized 2D histograms for the sample of pairs (Fig. 3a) and the
control sample of artificial pairs (Fig. 3b). Note that the upper panels in
this plot show histograms of the probability distribution of the angle (i.e.
irrespective of distance of satellite).

Fig. 3a shows that satellites are more likely to be situated along or near the
line connecting the binary (namely $\cos\theta\approx\pm 1$). Furthermore the
lopsided signal is driven by those satellite haloes that are far from their
host, near the binary center ($d_{\text{sat}}\approx 0.5d_{\text{sep}}$, more
on this effect in section 3.2.1). Note that the particle resolution inhibits
and artificially suppresses the ability for subhaloes to be identified in the
densest central regions of a halo (here clearly seen by the regions
$d_{\text{sat}}<0.1d_{\text{sep}}$). Regardless, Fig. 3a clearly shows that
the lopsided signal is driven by the fact that the satellite haloes that are
near the binary center also tend to extend along or near the line connecting
the two binary members.

Fig. 3b shows the cosine angle and radial distributions of the satellites for
the control sample of artificially generated pairs, meant to qualify the
effect of overlap. The overlapping effect clearly results in a gradual and
slight over-density towards the binary center at $\cos\theta\approx 1$ by
construction. However, the distribution of satellites of the artificially
constructed halo pairs is significantly different from the that of the real
sample. (A KS test reveals that the two distributions are highly unlikely to
be drawn from the same parent sample.) While the real sample is characterised
both by a flattening of the satellite distribution aligned with the partner
halo as well as by a lopsidedness bulging towards partner, the control sample
lacks the flattening but does exhibit lopsidedness due to the overlap effect
that we wish to capture. Is this overlap enough to explain the lopsided
signal? To some extent it does appear to be of the correct magnitude and the
$z=0$ effect could be at least partially ascribed to overlapping satellite
distributions. However, such an explanation would suggest that the lopsided
signal should decrease with increasing pair separation (i.e. with increasing
$z$), which, as we will show later in Section 4, is not the case.

####  3.2.1 Distance of satellite from host

In this section we examine how the lopsided signal changes as function of the
satellite’s distance to the host. We start by examining Fig. 3, where in the
top panels we plot the lopsided signal averaged over all distances and in the
bottom panels we show a “heat map” for the cosine of alignment angle $\theta$
as a function of satellite distance from host. Namely, each satellite provides
a $\cos\theta$ and a $d_{\text{sat}}/d_{\text{sep}}$ value, and we bin the
satellites in both values.

Fig. 3a indicates that for the satellites which are more distant from their
hosts, their distribution is more anisotropic. This is more explicitly shown
in Figure 4, where we examine the distribution of satellites in concentric
shells of increasing distance from their host.

Fig. 4 shows that the shape of the satellite distribution tends to be oblate
yet symmetric (i.e. not lopsided) for the regions closest to the host, and
only when the more distant satellites are considered: those at $d_{\rm
sat}>0.3~{}d_{\rm sep}$, does the distribution start to become more lopsided.
Such a finding suggests that it is the opposite binary partner which may be
responsible for destroying the aspherical but symmetric satellite
distributions, making them more oblate and lopsided.

(a) Less massive hosts (b) More massive hosts

Figure 5: We divide the lopsided signal into two subsamples, according to
whether the satellite belongs to the less massive member of the pair (panel
(a), left) or the more massive one (panel (b), right). All satellites are
defined as located within $d_{\rm sep}/2$ from their host.  Figure 6: The
lopsided signal is shown here for satellites of the less massive host partner,
split according to host mass ratio (defined as $M_{2}/M_{1}$ where
$M_{2}<M_{1}$). The lopsided signal is strongest for the most extreme mass
ratios. Systems of mass ratios larger than $0.4$ are not shown here as they do
not show substantial differences from those in the lower mass ratio range
(between 0.2 and 0.4).

####  3.2.2 Host mass and mass ratio

Here we examine how the mass ratio (defined as $M_{2}/M_{1}$ where
$M_{2}<M_{1}$) of the pair affects the lopsided satellite distribution. In
Fig. 5 we show the distribution of satellite position angle around the less
and more massive haloes independently. Fig. 5 reveals that the more massive
hosts show a weaker lopsided signal than the less massive partners.

We now examine the effect of the mass ratio on the probability distribution of
$\cos\theta$ for the less massive partner. The reader will recall from section
2.2 that the haloes in our sample may have masses anywhere in the range
$10^{11}\sim 5\times 10^{12}h^{-1}M_{\odot}$, i.e. may have a mass ratio from
0.02 to unity. In Figure 6 our sample is divided according to the mass ratio
of the pair, and the distribution of $\cos\theta$ for the satellites of the
less massive partner (only) is shown. Figure 6 shows that the lopsided
satellite signal is strongest for the smallest mass ratios (for example when
the more massive halo is 4 or 5 times larger than its smaller partner).

We conclude this section by stating that the lopsided signal is greatest in
the less massive halo pair, and for the most extreme pair halo mass ratios.
However, as in the previous section, we again ask the question if the signal
could be due to overlapping effects, especially suspect are the halo pairs
with large mass ratios, since the more massive halo could influence its
environment out to a few virial radii. Thus we construct artificial pairs that
match the mass ratio bins shown in Fig. 6 and examine the strength of the
lopsided effect due to overlap. In all cases the effect is indistinguishable
from the distributions shown in Fig. 3b (not shown). Given that the
distribution in the most extreme mass ratio is highly lopsided ($q=1.37$), it
is impossible that such an overlap effect can cause the lopsidedness shown in
Fig. 6.

##  4 The origin of the lopsided signal

In section 3.1 and 3.2 we established the existence of the lopsided satellite
distribution in pairs of haloes in a cosmological simulation at $z=0$. We now
attempt to uncover the origin of these distributions by following each
satellite’s orbital trajectory back in time. To do so, we use the Rockstar
Mergetree software described in section 2.1.

###  4.1 The role of infall

In order to examine the origin of the $z=0$ lopsided signal we wish to trace
back in time the most massive progenitor of all pairs and their satellites.
However, not all satellites can be traced back through previous snapshots,
since some of the smaller satellites will disappear as their mass fall below
the halo finder’s resolution. In order to keep the sample size at each
redshift the same, we are thus forced to cull from our sample those satellites
which during the back-tracking procedure fall below the resolution limit and
disappear from the merger tree. We arbitrarily set the redshift to which we
wish to track orbits to $z=2.5$ which maximises the sample size and the
redshift. Using $z=2.5$, we are forced to cull $\sim 35\%$ of the subhaloes
that are found at $z=0$ as these fall below the resolution limit at some point
between $z=0$ and 2.5.

In Fig. 7 we present five probability distribution of this sample as it is
traced back through cosmic time, at $z=0,0.5,1,1.5,2$ and 2.5. The reader will
first note that this culled sample still shows a lopsided signal at $z=0$. It
is immediately apparent that the lopsided signal was significantly stronger in
the past and that there is a considerable weakening as the system evolves
towards $z=0$. Namely, satellites of $z=0$ pairs, are more likely to be
accreted from (and therefore to have formed in) the regions between the two
hosts. Satellites are less likely to be accreted from (or to have formed in)
the regions diametrically on the other side, opposite the host’s partner. The
weakening of the lopsided signal as cosmic time advances may be quantified by
examining the lopsided metric $q$ as function $z$, shown in Fig. 8; a
monotonically decreasing trend is seen.

In Fig. 9 we examine $\cos\theta_{\rm rv}$, the angle formed between a
satellite’s position vector (with respect to the most massive progenitor of
its $z=0$ host), and its (relative) velocity vector. Namely at each redshift
under consideration, we examine how radial the satellite’s velocity vector is
with respect to its host. Note that a satellite with $\cos\theta_{\rm rv}=-1$
is moving radially towards its host with 0 angular momentum; a satellite with
$\cos\theta_{\rm rv}=1$ is moving radially away from its host, also with zero
angular momentum, while a satellite with $\cos\theta_{\rm rv}=0$ is moving
perpendicularly or tangentially to its position vector. A coherent picture
emerges: at higher redshifts, the distributions of $\cos\theta_{\rm rv}$ is
highly skewed, indicating a strong tendency for radial accretion. At low $z$,
as satellites approach their host, the velocity vector no longer has such a
strong tendency to be aligned with the position vector. This is likely due to
a combination of the gravitational pull of the host’s partner (as the pair
separation shrinks) as well as other “tidal” or nonlinear virial effects due
to the increase in number density of other satellite galaxies.

Figure 7: The evolution of the cosine of the angle $\theta$ as described in
Fig. 1. The lopsided signal distribution is stronger in the past than it is
today.
Figure 8: The evolution of the lopsided signals is described by the lopsided
metric $q$, shown here as a function of $z$ for all pairs considered in Fig.
7.
Figure 9: Evolution of cosine velocity incident angle distribution of
satellite haloes.

In Fig. 10 we show how a satellite’s velocity alignment angle
($\cos\theta_{\rm rv}$, see Fig. 1) depends on its distance (measured in units
of the pair separation at $z=0$, i.e. $d_{\rm sep0}$) for four different
redshifts, $z=0,0.5,1.5$ and 2.5. Fig. 10 shows that when satellites are very
distant, most come in on “radial-like” orbits i.e. with $\cos\theta_{\rm
rv}\approx-1$. As $z$ goes to 0, $\cos\theta_{\rm rv}$ begins to fill the
interval [-1,1], which might be an indication that many of the satellites gain
angular momentum. At $z=0.5$, where the region of satellite moving away from
their host (namely with $\cos\theta_{\rm rv}>0$) begins to be filled, we see
that it is almost entirely populated by satellites located at distances of
$<0.5d_{\rm sep0}$. Perhaps most pertinent to our study, we find that by $z=0$
up to $\sim$40% of the satellites are actually moving away from their hosts
(i.e. $\cos\theta_{\rm rv}>0$), after they probably have flown by their hosts
at some point in the past.

Since Fig. 10 indicates that when, for the most part, the lopsided signal is
strong, the angular momentum might be low and satellites are travelling on
infalling orbits (and vice versa), we hypothesize that an encounter between
the accreted satellite and its host may be responsible for both weakening the
lopsided effect and generating orbital angular momentum of the satellite. This
is a well established aspect of the virialization process which accompanies
halo relaxation, namely the stabilization of angular momentum from infalling
material (Stewart et al. 2013).

Figure 10: Normalized 2D histogram for the distribution of the radial
infalling angle $\theta_{\text{rv}}$ as a function of satellite distance from
host (in units of $d_{\rm sep0}$, the $z=0$ pair separation) at various
redshift. Colours are used to visualise which regions of these planes are more
populated by data points. Note that the sample can not extend beyond
$0.5d_{\rm sep0}$ at $z=0$, by construction. For higher redshift the red line
is used to mark the region that is searched at redshift 0.

###  4.2 Flybys and their impact on the lopsided signal

In order to test the hypothesis that accreted satellites that undergo “flyby”
event weaken the lopsided effect of satellite distribution, we split the
satellite sample up into those satellites that have “flown by” their host and
those that haven’t. This is done by following back in time the trajectory of
each satellite, and by examining the distance between the satellite and its
host. If this distance experiences a local minimum (that is not at $z=0$) then
the satellite has flown by its host.

In Fig. 11 we present the $z=0$ probability distribution of $\cos\theta$ for
those satellites labelled as flybys (a) and for those that have not
experienced an encounter (b). For completeness we show the same plot for the
full $z=0$ sample222The reader will note that this is identical to the black
curve in Fig. 7. It is not the same as the probability distribution shown in
Fig. 2 or 3a since it includes only those satellites that can be traced back
in time to $z=2.5$.. Fig. 11 strikingly demonstrates the difference in these
two samples. Those satellites that have incurred a flyby at one point in the
past, are both more centrally concentrated and therefore (see Fig. 4) more
uniformly distributed, and less lopsided. Those satellites that have not
experienced a close flyby are both more aligned with $\cos\theta\approx\pm 1$
and more lopsided. The lopsided signal appears to strongly be driven by
satellites that are infalling and who have not been gravitationally slung
about their host.

Note that the flybys have not fully erased the lopsided signal, they have
merely reduced the total lopsided signal to a very weak degree ($q=1.051$), so
weak in fact, that it is consistent with the effect of overlapping haloes.
Nevertheless, as we have seen so far, the overlapping effect alone cannot
explain the lopsidedness in the signal given that it weakens with decreasing
redshift.

Figure 11: Flyby signal (a) weakens the total lopsidedness of the signal (c)
at z = 0.

##  5 Discussion and Conclusions

The conclusion we draw regarding the origin of the lopsided signal is the
following. The satellites that cause the signal tend to be formed in between
the pair, in the filament connecting the pair. These satellites are then
attracted to each member of the pair due to gravity. Although those satellites
that are accreted and “swung” to the other side of the attracting pair member
may end up preferentially closer to the line connecting the binary, they are
outnumbered by the satellites still being accreted from the region between the
pair. The filamentary formation and accretion scenario may be related to the
3-point correlation function (Peebles, 1980) which has been measured in the
SDSS by McBride et al. (2011) but which shows only a very weak signal below
$\sim 6h^{-1}$ Mpc, (see also Ma & Fry, 2000; Jing & Börner, 2004), due mostly
to non-linear effects (such as flybys).

We examined how simply overlapping two satellite distributions would add to
the observed signal. By constructing artificial pairs we found that indeed at
z=0 overlapping satellite distributions have lopsidedness that is on the same
order of magnitude as that of the observed haloes. We conclude that, the
effects described above, namely that the signal is weakened by satellites that
have had a pericentric passage and that the lopsidedness is stronger at higher
z, rule out overlapping satellite distributions as being primarily
responsible. In summary: the lopsided signal originates at high redshift and
is weakened as the halo pair distance shrinks and satellites are accreted and
have pericentric passages. Such dynamical effects weaken the signal to the
degree that at z=0 it is consistent with artificially created overlapping
satellite distributions but not due to it.

As mentioned in Sec. 1, Libeskind et al. (2016) examined how satellites were
spatially distributed in spectroscopically identified SDSS pairs. Pawlowski et
al. (2017) went on to examine whether such satellite distributions could be
found in $\Lambda$CDM cosmological simulations, coming to the conclusion that
indeed such a lopsided distribution was consistent with what the cosmological
paradigm predicted. In order to follow up on the question of “why do
satellites tend to be in the in-between regions of galaxy pairs?”, in this
work we examined how pairs of galaxies identified in cosmological simulations
at $z=0$ accrete satellites. By employing a large, well resolved, cosmological
simulation we are able to numerically follow the formation of such pairs and
their satellites. We trace back each host progenitor and each satellite to
$z=2.5$. At each snapshot we construct the triangle composed of the
(progenitors of the) two members of the binary and the (progenitor of the)
satellite and examine the evolution of the angle satellite – host – host
partner. We find that the signal is strongly driven by satellites that are on
their first approach to the binary pair (by “first approach” we mean satellite
that have never been closer to their host). Satellites that were at some point
in the past closer to their hosts (a category we call “flyby” satellites)
exhibit a much weaker lopsided distribution. Lopsided distributions of
satellite galaxies among pairs of galaxies are thus characterised by
dynamically young accretion events that have either recently come to the pair,
or are still in the process of being accreted. This is similar to the picture
described by Wetzel (2011) where accretion is more radial at higher $z$. We
have also found that the lopsided distribution is more often found in the
lower mass member of the pair, where lower mass implies lower degree of
virialization and hence a stronger lopsided signal, especially when its
partner is much larger. A higher mass member of the pair has in general more
satellites which are accreted, and hence corresponds to a denser environment.

The pairs examined in this paper (or indeed in other works that focus on the
lopsided signal) are not identical to the Local Group or the Centaurus group,
although such a lopsided distribution is also seen in the distribution of
satellites along Andromeda where 19-23 out of 27 satellites are on the near
side of M31 (Ibata et al. 2013) and (albeit to a lesser extent) around the
Milky Way where the census is poorer due to obscuration from the galactic disk
and inhomogeneous coverage from sky surveys such as the SDSS (McConnachie,
2012), as well as in the Centaurus A - M83 pair (Müller et al. 2015). That
being said, the selection criteria employed here is loosely motivated such
that the Local Group would also be identified according to the definitions we
use. Our criteria are purposefully kept more relaxed such that we may get a
more general picture of these systems.

Our work confirms one of the assertions of Pawlowski et al. (2017) that
lopsided distributions are not in conflict with $\Lambda$CDM models of
structure formation, at least when comparing surveys like the SDSS to the
virialized haloes that form in dark matter only $N$-body simulations. However,
our discovery that the origin of these systems is primarily newly accreted
satellites is somewhat surprising, since its well known that satellites may
orbit their hosts for many gigayears before being accreted (for example:
Tormen et al., 1998; Yang et al., 2009; Bahé et al., 2019, among others).

We leave the reader with the following notion regarding satellites of “Local
Group like” pairs: if the existence of lopsided distributions in surveys and
simulations is driven by satellites on first approach, we infer that, either
these are a different population and unlike satellites found within the virial
radius of Milky Way like haloes that may survive for many giga years, or the
lopsided signal in Local Group is also driven by such mostly accreted
satellites. Indeed, recent work by Hammer et al. (2018) suggests that this may
be the case and that the satellites of the Milky Way may indeed be on their
first approach.

## Acknowledgments

NIL acknowledges financial support of the Project IDEXLYON at the University
of Lyon under the Investments for the Future Program (ANR-16-IDEX-0005). NIL,
PW and QG, acknowledge support from the joint Sino-German DFG research Project
“The Cosmic Web and its impact on galaxy formation and alignment” (DFG-LI
2015/5-1). ET was supported by ETAg grants (IUT40-2, IUT26-2) and by EU
through the ERDF CoE grants TK133 and MOBTP86. GY acknowledges support by
Ministerio de Economia y Competitividad and the Fondo Europeo de Desarrollo
Regional (MINECO/FEDER, UE) in Spain through grants AYA2015-63810-P and
PGC2018-094975-B-C21. JS acknowledges support from the “Centre National
d’études spatiales (CNES)” postdoctoral fellowship program. The ESMD
simulation has been performed at LRZ Munich within the project pr74no.

## References

  * Agustsson & Brainerd (2010) Agustsson I., Brainerd T. G., 2010, ApJ, 709, 1321
  * Angulo & White (2010) Angulo R. E., White S. D. M., 2010, MNRAS, 401, 1796
  * Arakelyan et al. (2018) Arakelyan N. R., Pilipenko S. V., Libeskind N. I., 2018, MNRAS, 481, 918
  * Bahé et al. (2019) Bahé Y. M., et al., 2019, MNRAS, p. 397
  * Behroozi et al. (2012) Behroozi P. S., Wechsler R. H., Wu H.-Y., 2012, ApJ, 762, 109
  * Behroozi et al. (2015) Behroozi P. S., et al., 2015, MNRAS, 450, 1546
  * Bowden et al. (2014) Bowden A., Evans N. W., Belokurov V., 2014, ApJ, 793, L42
  * Bullock & Boylan-Kolchin (2017) Bullock J. S., Boylan-Kolchin M., 2017, ARA&A, 55, 343
  * Cautun et al. (2014) Cautun M., Hellwing W. A., van de Weygaert R., Frenk C. S., Jones B. J. T., Sawala T., 2014, Monthly Notices of the Royal Astronomical Society, 445, 1820
  * Cautun et al. (2015a) Cautun M., Wang W., Frenk C. S., Sawala T., 2015a, MNRAS, 449, 2576
  * Cautun et al. (2015b) Cautun M., Bose S., Frenk C. S., Guo Q., Han J., Hellwing W. A., Sawala T., Wang W., 2015b, MNRAS, 452, 3838
  * De Graaff et al. (2017) De Graaff A., Cai Y.-C., Heymans C., Peacock J. A., 2017 (arXiv:1303.5076)
  * Ellison et al. (2008) Ellison S. L., Patton D. R., Simard L., McConnachie A. W., 2008, AJ, 135, 1877
  * Epps & Hudson (2017) Epps S. D., Hudson M. J., 2017, MNRAS, 468, 2605
  * Fakhouri & Ma (2008) Fakhouri O., Ma C.-P., 2008, MNRAS, 386, 577
  * Genel et al. (2010) Genel S., Bouché N., Naab T., Sternberg A., Genzel R., 2010, ApJ, 719, 229
  * Gordon et al. (2017) Gordon Y. A., et al., 2017, MNRAS, 465, 2671
  * Hammer et al. (2018) Hammer F., Yang Y., Wang J., Arenou F., Puech M., Flores H., 2018, arXiv e-prints, p. arXiv:1812.10714
  * Holmberg (1969) Holmberg E., 1969, Arkiv for Astronomi, 5, 305
  * Ibata et al. (2013) Ibata R. A., et al., 2013, Nature, 493, 62
  * Ibata et al. (2014a) Ibata N. G., Ibata R. A., Famaey B., Lewis G. F., 2014a, Nature, 511, 563
  * Ibata et al. (2014b) Ibata R. A., Ibata N. G., Lewis G. F., Martin N. F., Conn A., Elahi P., Arias V., Fernando N., 2014b, ApJ, 784, L6
  * Jing & Börner (2004) Jing Y. P., Börner G., 2004, ApJ, 607, 140
  * Kahn & Woltjer (1959) Kahn F. D., Woltjer L., 1959, ApJ, 130, 705
  * Kauffmann et al. (2013) Kauffmann G., Li C., Zhang W., Weinmann S., 2013, MNRAS, 430, 1447
  * Kelley et al. (2018) Kelley T., Bullock J., Garrison-Kimmel S., Boylan-Kolchin M., Pawlowski M., S. Graus A., 2018 (arXiv:1303.5076)
  * Klypin et al. (2016) Klypin A., Yepes G., Gottlöber S., Prada F., Heß S., 2016, MNRAS, 457, 4340
  * Koposov et al. (2008) Koposov S., et al., 2008, ApJ, 686, 279
  * Kroupa et al. (2005) Kroupa P., Theis C., Boily C. M., 2005, A&A, 431, 517
  * L’Huillier et al. (2012) L’Huillier B., Combes F., Semelin B., 2012, A&A, 544, A68
  * Lacey & Cole (1993) Lacey C., Cole S., 1993, MNRAS, 262, 627
  * Libeskind et al. (2005) Libeskind N. I., Frenk C. S., Cole S., Helly J. C., Jenkins A., Navarro J. F., Power C., 2005, MNRAS, 363, 146
  * Libeskind et al. (2015a) Libeskind N. I., Guo Q., Kipper R., Tempel E., 2015a, MNRAS, 450, 2727
  * Libeskind et al. (2015b) Libeskind N. I., Hoffman Y., Tully R. B., Courtois H. M., Pomarède D., Gottlöber S., Steinmetz M., 2015b, MNRAS, 452, 1052
  * Libeskind et al. (2016) Libeskind N. I., Guo Q., Tempel E., Ibata R., 2016, ApJ, 830, 121
  * Lüghausen et al. (2014) Lüghausen F., et al., 2014, MNRAS, 442, 2362
  * Ma & Fry (2000) Ma C.-P., Fry J. N., 2000, ApJ, 531, L87
  * Madau et al. (2008) Madau P., Diemand J., Kuhlen M., 2008, ApJ, 679, 1260
  * McBride et al. (2011) McBride C. K., Connolly A. J., Gardner J. P., Scranton R., Newman J. A., Scoccimarro R., Zehavi I., Schneider D. P., 2011, ApJ, 726, 13
  * McConnachie (2012) McConnachie A. W., 2012, AJ, 144, 4
  * Metz et al. (2008) Metz M., Kroupa P., Libeskind N. I., 2008, ApJ, 680, 287
  * Moreno et al. (2013) Moreno J., Bluck A. F. L., Ellison S. L., Patton D. R., Torrey P., Moster B. P., 2013, MNRAS, 436, 1765
  * Müller et al. (2015) Müller O., Jerjen H., Binggeli B., 2015, A&A, 583, A79
  * Müller et al. (2016) Müller O., Jerjen H., Pawlowski M. S., Binggeli B., 2016, A&A, 595, A119
  * Müller et al. (2017) Müller O., Scalera R., Binggeli B., Jerjen H., 2017, A&A, 602, A119
  * Müller et al. (2018) Müller O., Pawlowski M. S., Jerjen H., Lelli F., 2018, Science, 359, 534
  * Muñoz et al. (2015) Muñoz R. P., et al., 2015, The Astrophysical Journal, 813, L15
  * Newton et al. (2018) Newton O., Cautun M., Jenkins A., Frenk C. S., Helly J. C., 2018, MNRAS, 479, 2853
  * Ordenes-Briceño et al. (2018) Ordenes-Briceño Y., et al., 2018, The Astrophysical Journal, 859, 52
  * Partridge et al. (2013) Partridge C., Lahav O., Hoffman Y., 2013, MNRAS: Letters, 436, L45
  * Pawlowski (2018) Pawlowski M., 2018, Modern Physics Letters A, 33
  * Pawlowski et al. (2012) Pawlowski M. S., Pflamm-Altenburg J., Kroupa P., 2012, Monthly Notices of the Royal Astronomical Society, 423, 1109
  * Pawlowski et al. (2013) Pawlowski M. S., Kroupa P., Jerjen H., 2013, MNRAS, 435, 1928
  * Pawlowski et al. (2017) Pawlowski M. S., Ibata R. A., Bullock J. S., 2017, ApJ, 850, 132
  * Peebles (1980) Peebles P. J. E., 1980, The large-scale structure of the universe. Princeton University Press, Princeton, N.J.
  * Perez et al. (2009) Perez J., Tissera P., Blaizot J., 2009, MNRAS, 397, 748
  * Planck Collaboration et al. (2014) Planck Collaboration et al., 2014, A&A, 571, A16
  * Prada et al. (2006) Prada F., Klypin A. A., Simonneau E., Betancort-Rijo J., Patiri S., Gottlober S., Sanchez-Conde M. A., 2006, The Astrophysical Journal, 645, 1001
  * Reed et al. (2017) Reed D. S., Schneider A., Potter D., Stadel J., Moore B., Smith R. E., 2017\. WORLD SCIENTIFIC, pp 2453–2458
  * Rodríguez-Puebla et al. (2016) Rodríguez-Puebla A., Behroozi P., Primack J., Klypin A., Lee C., Hellinger D., 2016, MNRAS, 462, 893
  * Sales & Lambas (2004) Sales L., Lambas D. G., 2004, MNRAS, 348, 1236
  * Sales et al. (2017) Sales L. V., et al., 2017, MNRAS, 464, 2419
  * Shen et al. (2016) Shen S.-Y., et al., 2016, Research in Astronomy and Astrophysics, 16, 43
  * Stewart et al. (2008) Stewart K. R., Bullock J. S., Wechsler R. H., Maller A. H., Zentner A. R., 2008, ApJ, 683, 597
  * Stewart et al. (2013) Stewart K. R., Brooks A. M., Bullock J. S., Maller A. H., Diemand J., Wadsley J., Moustakas L. A., 2013, ApJ, 769, 74
  * Taylor et al. (2018) Taylor M. A., Eigenthaler P., Puzia T. H., Muñoz R. P., Ribbeck K. X., Zhang H.-X., Ordenes-Briceño Y., Bovill M. S., 2018, The Astrophysical Journal, 867, L15
  * Tempel & Tamm (2015) Tempel E., Tamm A., 2015, A&A, 576, L5
  * Tollerud et al. (2008) Tollerud E. J., Bullock J. S., Strigari L. E., Willman B., 2008, ApJ, 688, 277
  * Tormen et al. (1998) Tormen G., Diaferio A., Syer D., 1998, MNRAS, 299, 728
  * Tully et al. (2015) Tully R. B., Libeskind N. I., Karachentsev I. D., Karachentseva V. E., Rizzi L., Shaya E. J., 2015, ApJ, 802, L25
  * Wang et al. (2019) Wang P., Guo Q., Libeskind N. I., Tempel E., Wei C., Kang X., 2019, MNRAS, 484, 4325
  * Wetzel (2011) Wetzel A. R., 2011, MNRAS, 412, 49
  * Wetzel & White (2010) Wetzel A. R., White M., 2010, MNRAS, 403, 1072
  * White & Rees (1978) White S. D. M., Rees M. J., 1978, MNRAS, 183, 341
  * Yang et al. (2009) Yang X., Mo H. J., van den Bosch F. C., 2009, ApJ, 693, 830
  * Zeldovich (1980) Zeldovich Y. B., 1980, MNRAS, 192, 663


# Lipschitz regularity of energy-minimal mappings between doubly connected
Riemann surfaces

David Kalaj  Faculty of Natural Sciences and Mathematics, University of
Montenegro, Cetinjski put b.b. 81000 Podgorica, Montenegro<EMAIL_ADDRESS>
###### Abstract.

Let $M$ and $N$ be doubly connected Riemann surfaces with
$\mathscr{C}^{1,\alpha}$ boundaries and with nonvanishing conformal metrics
$\sigma$ and $\wp$ respectively, and assume that $\wp$ is a smooth metric with
bounded Gauss curvature $\mathcal{K}$ and finite area. Assume that
${\mathsf{H}}_{\rho}(M,N)$ is the class of all $\mathscr{W}^{1,2}$
bomeomorphisms between $M$ and $N$ and assume that
$\mathcal{E}^{\wp}:\overline{\mathsf{H}}_{\rho}(M,N)\to\mathbf{R}$ is the
Dirichlet-energy functional, where $\overline{\mathsf{H}}_{\rho}(M,N)$ is the
closure of ${\mathsf{H}}_{\rho}(M,N)$ in $\mathscr{W}^{1,2}(M,N)$. By using a
result of Iwaniec, Kovalev and Onninen in [8] that the minimizer, is locally
Lipschitz, we prove that the minimizer, of the energy functional
$\mathcal{E}^{\wp}$, which is not a diffeomorphism in general, is a globally
Lipschitz mapping of $M$ onto $N$.

###### Key words and phrases:

Dirichlet energy, Riemann surfaces, minimal energy, harmonic mapping,
conformal modulus

###### 2000 Mathematics Subject Classification:

Primary 58E20; Secondary 30C62, 31A05

###### Contents

  1. 1 Introduction
    1. 1.1 Harmonic mappings between Riemann surfaces
    2. 1.2 Statement of results
  2. 2 Background and preliminaries
    1. 2.1 Stationary mappings
  3. 3 Proof of the main result

##  1\. Introduction

The primary goal of this paper is to study Lipschitz behaviour of stationary
deformations of the Dirichlet energy of mappings between doubly-connected
Riemann surfaces. The existence part has been proved. More precisely, Iwaniec,
Koch, Kovalev and Onninen in [7] proved that there exists so-called
deformation $f$ that maps a double connected domain $\mathbb{X}$ onto a doubly
connected domain $\mathbb{Y}$ in the complex plane and which minimizes the
Diriclet energy integral throughout the class os deformations
$\mathcal{D}(\mathbb{X},\mathbb{Y})$ which contains the class of all Sobolev
$\mathscr{W}^{1,2}(\mathbb{X})$ homeomorphism. The minimizer is a harmonic
diffeomorphism, provided that the conformal moduli satisfies the relation
$\mathrm{Mod}(\mathbb{X})\leqslant\mathrm{Mod}(\mathbb{Y})$. Moreover, the
Dirichlet energy is invariant under conformal change of the original domain.
This is why the original domain can be chosen to be equal to
$\mathbb{X}=\mathbb{A}(r,R):=\\{z:r<|z|<R\\}$. The Hopf’s differential defined
by $\mathrm{Hopf}(f)=f_{z}\overline{f}_{z}$ has very special form for so-
called stationary deformations namely

$\mathrm{Hopf}(f)(z)=\frac{c}{z^{2}},\ \ z\in\mathbb{A}(r,R).$

Later this result has been generalized by the author in [13] for so-called
admissible metrics $\wp$ in $\mathbb{Y}$.

In this case Hopf’s differential is defined by

$\mathrm{Hopf}(f)=\wp^{2}(f(z))f_{z}\overline{f}_{z}$

and has a very special form for a so-called stationary deformation namely

$\mathrm{Hopf}(f)(z)=\frac{c}{z^{2}},\ \ z\in\mathbb{A}(r,R).$

Further in [8], Iwaniec, Kovalev and Onninen proved that every stationary
deformation is locally Lipschitz in the domain.

In [11] the author proved that if
$f:\mathbb{X}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\mathbb{Y}$
is a $\rho-$harmonic diffeomorphic minimizer, and
$\partial\mathbb{X},\partial\mathbb{Y}\in\mathscr{C}^{2}$ then $f$ is
Lipschitz continuous up to the boundary.

In [14] the author and Lamel proved that a minimizer of the Euclidean
Dirichlet energy which is a diffeomophic surjection of $\mathbb{X}$ onto
$\mathbb{Y}$ has smooth extension up to the boundary. More precisely they
proved that, they proved that, if
$f:\mathbb{X}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\mathbb{Y}$
is a Euclidean diffeomorphic minimizer, so that
$\partial\mathbb{X},\partial\mathbb{Y}\in\mathscr{C}^{1,\alpha}$, then the
$f\in\mathscr{C}^{1,\alpha^{\prime}}(\overline{\mathbb{X}})$, where
$\alpha^{\prime}=\alpha$, if
$\mathrm{Mod}(\mathbb{X})\geqslant\mathrm{Mod}(\mathbb{Y})$ and
$\alpha^{\prime}=\alpha/(2+\alpha)$ if
$\mathrm{Mod}(\mathbb{X})<\mathrm{Mod}(\mathbb{Y})$.

In [12] the author extended the main result in [14] and proved the following
extension of the Kellogg theorem. Every diffeomorphic minimiser of Dirichlet
energy of Sobolev mappings between doubly connected Riemanian surfaces
$(\mathbb{X},\sigma)$ and $(\mathbb{Y},\rho)$ having $\mathscr{C}^{n,\alpha}$
boundary, $0<\alpha<1$, is $\mathscr{C}^{n,\alpha}$ up to the boundary,
provided the metric $\rho$ is smooth enough. Here $n$ is a positive integer.
It is crucial that, every diffeomorphic minimizer of Dirichlet energy is a
harmonic mapping with a very special Hopf differential and this fact is used
in the proof.

###  1.1. Harmonic mappings between Riemann surfaces

Let $M=(\mathbb{X},\sigma)$ and $N=(\mathbb{Y},\wp)$ be Riemann surfaces with
metrics $\sigma$ and $\wp$, respectively. If a mapping $f:M\to N,$ is $C^{2}$,
then $f$ is said to be harmonic (to avoid the confusion we will sometimes say
$\wp$-harmonic) if

(1.1) $f_{z\overline{z}}+{(\log\wp^{2})}_{w}\circ ff_{z}\,f_{\bar{z}}=0,$

where $z$ and $w$ are the local parameters on $M$ and $N$ respectively. Also
$f$ satisfies (1.1) if and only if its Hopf differential

(1.2) $\mathrm{Hopf}(f)=\wp^{2}\circ ff_{z}\overline{f_{\bar{z}}}$

is a holomorphic quadratic differential on $M$. Let

$|\partial f|^{2}:=\frac{\wp^{2}(f(z))}{\sigma^{2}(z)}\left|\frac{\partial
f}{\partial z}\right|^{2}\text{ and
}|\bar{\partial}f|^{2}:=\frac{\wp^{2}(f(z))}{\sigma^{2}(z)}\left|\frac{\partial
f}{\partial\bar{z}}\right|^{2}$

where $\frac{\partial f}{\partial z}$ and $\frac{\partial f}{\partial\bar{z}}$
are standard complex partial derivatives. The $\wp-$Jacobian is defined by

$J(f):=|\partial f|^{2}-|\partial\bar{f}|^{2}.$

If $u$ is sense preserving, then the $\wp-$Jacobian is positive. The Hilbert-
Schmidt norm of differential $df$ is the square root of the energy $e(f)$ and
is defined by

(1.3) $|df|=\sqrt{2|\partial f|^{2}+2|\partial\bar{f}|^{2}}.$

For $g:M\mapsto N$ the $\wp-$ _Dirichlet energy_ is defined by

(1.4) $\mathcal{E}^{\wp}[g]=\int_{M}|dg|^{2}dV_{\sigma},$

where $\partial g$, and $\bar{\partial}g$ are the partial derivatives taken
with respect to the metrics $\varrho$ and $\sigma$, and $dV_{\sigma}$, is the
volume element on $(M,\sigma)$, which in local coordinates takes the form
$\sigma^{2}(w)du\wedge dv,$ $w=u+iv$. Assume that energy integral of $f$ is
bounded. Then a stationary point $f$ of the corresponding functional where the
homotopy class of $f$ is the range of this functional is a harmonic mapping.
The converse is not true. More precisely there exists a harmonic mapping which
is not stationary. The literature is not always clear on this point because
for many authors, a harmonic mapping is a stationary point of the energy
integral. For the last definition and some important properties of harmonic
maps see [9]. It follows from the definition that, if $a$ is conformal and $f$
is harmonic, then $f\circ a$ is harmonic. Moreover if $b$ is conformal,
$b\circ f$ is also harmonic but with respect to (possibly) an another metric
$\wp_{1}$.

Moreover if $N$ and $M$ are double connected Riemann surfaces with non
vanishing metrics $\sigma$ and $\wp$, then by [9, Theorem 3.1], there exist
conformal mappings
$X\colon\mathbb{X}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}M$
and
$X^{*}\colon\mathbb{Y}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}N$
between double connected plane domains $\mathbb{X}$ and $\mathbb{Y}$ and
Riemann surfaces $M$ and $N$ respectively.

Notice that the harmonicity neither Dirichlet energy do not depend on metric
$\sigma$ on domain so we will assume from now on $\sigma(z)\equiv 1$. This is
why throughout this paper $M=(\mathbb{X},\mathbf{1})$ and $N=(\mathbb{Y},\wp)$
will be doubly connected domains in the complex plane $\mathbb{C}$ (possibly
unbounded), where $\mathbf{1}$ is the Euclidean metric. Moreover $\wp$ is a
nonvanishing smooth metric defined in $\mathbb{Y}$ with bounded Gauss
curvature $\mathcal{K}$ where

(1.5) $\mathcal{K}(z)=-\frac{\Delta\log\wp(z)}{\wp^{2}(z)},$

(we put $\kappa:=\sup_{z\in\mathbb{Y}}|\mathcal{K}(z)|<\infty$) and with
finite area defined by

$\mathcal{A}(\wp)=\int_{\mathbb{Y}}\wp^{2}(w)dudv,\ \ w=u+iv.$

We call a metric $\rho$ an admissible one if there is a constant $C_{\wp}>0$
so that

(1.6) ${|\nabla\wp(w)|}\leqslant C_{\wp}{\wp(w)},\ \ \ w\in\mathbb{Y}\ \ i.e.\
\ \nabla\log\wp\in L^{\infty}(\mathbb{Y})$

which means that $\wp$ is so-called approximately analytic function (c.f.
[5]).

Assume that the domain of $\wp$ is the unit disk
$\mathbf{D}:=\\{z:|z|<1\\}\subset\mathbf{C}$. From (1.6) and boundedness of
$\wp$, it follows that it is Lipschitz, and so it is continuous up to the
boundary. Again by using (1.6), the function $\psi(t)=\wp(te^{i\alpha})$,
$0<t<1$, $\alpha\in[0,2\pi]$ satisfies the differential inequalities
$-C_{\wp}\leqslant\partial_{t}\log\psi(t)\leqslant C_{\wp}$, which by
integrating in $[0,t]$ imply that
$\psi(0)e^{-C_{\wp}t}\leqslant\psi(t)\leqslant\psi(0)e^{C_{\wp}t}$. Therefore
under the above conditions there holds the double inequality

(1.7)
$0<\wp(0)e^{-C_{\wp}}\leqslant{\wp(w)}\leqslant\wp(0)e^{C_{\wp}}<\infty,\ \
w\in\mathbf{D}.$

A similar inequality to (1.7) can be proved for $\mathbb{Y}$ instead of
$\mathbf{D}$. The Euclidean metric ($\wp\equiv 1$) is an admissible metric.
The Riemannian metric defined by $\wp(w)={1}/{(1+|w|^{2})^{2}}$ is admissible
as well. The Hyperbolic metric $h(w)={1}/{(1-|w|^{2})^{2}}$ is not an
admissible metric on the unit disk neither on the annuli
$\mathbb{A}(r,1):=\\{z:r<|z|<1\\}$, but it is admissible in
$\mathbb{A}(r,R):=\\{z:r<|z|<R\\}$, where $0<r<R<1$. We call such a metric
_allowable_ one (cf. [1, P. 11]). If $\wp$ is a given metric in $\mathbb{Y}$,
we conventionally extend it to be equal to $0$ in $\partial\mathbb{Y}$. As we
already pointed out, we will study the minimum of Dirichlet integral of
mappings between certain sets. We refer to introduction of [7] and references
therein for good setting of this problem and some connection with the theory
of nonlinear elasticity. Notice first that a change of variables $w=f(z)$ in
(1.4) yields

(1.8)
${\mathcal{E}^{\wp}}[f]=2\int_{\mathbb{X}}\wp^{2}(f(z))J_{f}(z)\,dz+4\int_{\mathbb{X}}\wp^{2}(f(z))\lvert
f_{\bar{z}}\rvert^{2}dz\geqslant 2\mathcal{A}(\wp)$

where $J_{f}$ is the Jacobian determinant and $\mathcal{A}(\wp)$ is the area
of $\mathbb{Y}$ and $dz:=dx\wedge dy$ is the area element w.r. to Lebesgue
measure on the complex plane. A conformal mapping of
$f:\mathbb{X}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\mathbb{Y}$;
that is, a homeomorphic solution of the Cauchy-Riemann system $f_{\bar{z}}=0$,
would be an obvious choice for the minimizer of (1.8). For arbitrary multiply
connected domains there is no such mapping.

Any energy-minimal diffeomorphism satisfies Euler-Laplace’s equation, since
one can perform first variations while preserving the diffeomorphism property.
However, in the case of Euclidean metric $\wp\equiv 1$, the existence of a
harmonic diffeomorphism between certain sets does not imply the existence of
an energy-minimal one, see [7, Example 9.1]. Example 9.1 in [7] has been
constructed with help of affine self-mappings of the complex plane. For a
general metric $\wp$, affine transformations are not harmonic, thus we cannot
produce a similar example.

###  1.2. Statement of results

The main results of this paper is the following extension of the main result
in [11].

######  Theorem 1.1.

Suppose that $\mathbb{X}$ and $\mathbb{Y}$ are doubly connected domains in
$\mathbb{C}$ with $\mathscr{C}^{1,\alpha}$ boundaries. Assume that
$\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$ is the class of deformations as in
definition (2.2) below. Let $\wp$ be an admissible metric in $\mathbb{Y}$.
Then there exists a deformation $w$ that minimizes $\wp-$energy throughout the
class of deformations $\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$, and it is
Lipschitz continuous up to the boundary of $\mathbb{X}$.

######  Remark 1.2.

In contrast to the main result obtained in [11], in Theorem 1.1 the minimizer
is not necessarily a diffeomorphism, and such a minimizer exists almost always
except in some degenerate cases (see Proposition 2.5 below).

The following example is taken from [13] for general radial metric. For a
corresponing example for Euclidean metric see [7].

######  Example 1.3.

Assume that $\wp(w)=\rho(|w|)$ is a radial metric defined in the annulus
$\Omega^{*}=\mathbb{A}(1,R)$. Assume that $\Omega=\mathbb{A}(1,r)$. Chose

$r<\exp\left(\int_{1}^{R}\frac{\rho(y)dy}{\sqrt{y^{2}\rho^{2}(y)-R^{2}\rho^{2}(R)}}\right).$

Let

$r_{\diamond}=\exp\left(\int_{1}^{R}\frac{\rho(y)dy}{\sqrt{y^{2}\rho^{2}(y)-R^{2}\rho^{2}(R)}}\right).$

Then $r<r_{\diamond}<1$.

Further for $\gamma=-\delta^{2}\rho^{2}(\delta)$, we have well defined
function

(1.9)
$q_{\diamond}(s)=\exp\left(\int_{\sigma}^{s}\frac{dy}{\sqrt{y^{2}-\delta^{2}\rho^{2}(\delta)\varrho^{2}}}\right),\delta\leqslant
s\leqslant\sigma.$

Then the infimum $\mathsf{E}^{\rho}(\Omega,\Omega^{*})$ is realized by a non-
injective deformation
$h\colon\Omega\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\Omega^{\ast}$

$h(z)=\begin{cases}R\frac{z}{\lvert z\rvert}&\mbox{for }r<\lvert
z\rvert\leqslant r_{\diamond}\\\ h_{\diamond}(z)&\mbox{ for
}r_{\diamond}\leqslant\lvert z\rvert<1\end{cases}$

where
$h_{\diamond}(z)=p_{\diamond}(|z|)e^{it}=(q_{\diamond})^{-1}(|z|)e^{it}$,
$z=|z|e^{it}$ and $q_{\diamond}$ is defined in (1.9) below

Here the radial projection $z\mapsto Rz/\lvert z\rvert$ hammers
$A(r,r_{\diamond})$ onto the circle $|z|=R$ while the critical $\rho-$Nitsche
mapping $h_{\diamond}$ takes $A(r_{\diamond},1)$ homeomorphically onto
$\Omega^{\ast}$.

##  2\. Background and preliminaries

A homeomorphism of a planar domain is either sense-preserving or sense-
reversing. For homeomorphisms of the Sobolev class $\mathscr{W}^{1,1}_{\rm
loc}(\mathbb{X},\mathbb{Y})$ this implies that the Jacobian determinant does
not change sign: it is either nonnegative or nonpositive at almost every point
[3, Theorem 3.3.4]. The homeomorphisms considered in this paper are sense-
preserving.

Let $\mathbb{X}$ and $\mathbb{Y}$ be domains in $\mathbb{C}$. Let
$a\in\mathbb{Y}$ and $b\in\partial\mathbb{Y}$. We define

$\operatorname{dist}_{\wp}(a,b):=\inf_{\gamma}\int_{\gamma}\wp(w)|dw|,$

where $\gamma$ ranges over all rectifiable Jordan arcs connecting $a$ and $b$
within $\mathbb{Y}$ if the set of such Jordan arcs is not empty (otherwise we
conventionally put $\operatorname{dist}_{\wp}(a,b)=\infty$). To every mapping
$f\colon\mathbb{X}\to\overline{\mathbb{Y}}$ we associate a boundary distance
function

$\delta^{\wp}_{f}(z)=\operatorname{dist}_{\wp}(f(z),\partial\mathbb{Y})=\inf_{b\in\partial\mathbb{Y}}\operatorname{dist}_{\wp}(f(z),b)$

which is set to $0$ on the boundary of $\mathbb{X}$.

We now adapt for our purpose the concepts of $c\delta-$uniform convergence and
of deformation defined for Euclidean metric and bounded domains in [7].

######  Definition 2.1.

A sequence of continuous mappings $h_{j}\colon\mathbb{X}\to{\mathbb{Y}}$ is
said to converge _$c\delta$ -uniformly_ to
$h\colon\mathbb{X}\to\overline{\mathbb{Y}}$ if

  1. (1)

$h_{j}\to h$ uniformly on compact subsets of $\mathbb{X}$ and

  2. (2)

$\delta^{\wp}_{h_{j}}\to\delta^{\wp}_{h}$ uniformly on
$\overline{\mathbb{X}}$.

We designate it as $h_{j}\xrightarrow[]{{}_{c\delta}}h$. Concerning the item
(1) we need to notice that the Euclidean metric and the nonvanishing smooth
metric $\wp$ are equivalent on compacts of $\mathbb{Y}$.

######  Definition 2.2.

A mapping $h\colon\mathbb{X}\to\overline{\mathbb{Y}}$ is called a _$\wp$
deformation_ if

  1. (1)

$h\in W_{loc}^{1,2}$ and $|dh|\in L^{2}$ (which we write shortly $h\in
W_{l}^{1,2}$);

  2. (2)

The Jacobian $J_{h}:=\det Dh$ is nonnegative a.e. in $\mathbb{X}$;

  3. (3)

$\int_{\mathbb{X}}\wp^{2}(h(z))J_{h}\leqslant\mathcal{A}(\wp)$;

  4. (4)

there exist sense-preserving diffeomorphisms
$h_{j}\colon\mathbb{X}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\mathbb{Y}$,
called an _approximating sequence_ , such that
$h_{j}\xrightarrow[]{{}_{c\delta}}h$ on $\mathbb{X}$.

The set of $\wp$ deformations $h\colon\mathbb{X}\to\overline{\mathbb{Y}}$ is
denoted by $\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$.

Notice first that the condition (1) of the previous definition in the case of
bounded domains $\mathbb{X}$ and $\mathbb{Y}$ can be replaced by $h\in
W^{1,2}$ see [13, Lemma 2.3.].

######  Lemma 2.3.

Let $\mathbb{X}$ and $\mathbb{Y}$ be bounded doubly connected planar domains.
Assume that the boundary components of $\mathbb{X}$ do not degenerate into
points. If a sequence
$\\{h_{j}\\}\subset\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$ converges weakly
in $W^{1,2}$, then its limit belongs to
$\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$. In particular, every such limit
function has continuous extension to the boundary and the boundary is mapped
into the boundary.

Now we formulate the following existence result proved in [13].

######  Proposition 2.4.

[13, Lemma 2.15] Let $\mathbb{X}$ and $\mathbb{Y}$ be bounded doubly connected
planar domains. Assume that the boundary components of $\mathbb{X}$ do not
degenerate into points. There exists
$h\in\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$ such that
${\mathcal{E}^{\wp}}[h]=\mathsf{E}^{\wp}(\mathbb{X},\mathbb{Y})$, where

$\mathsf{E}^{\wp}(\mathbb{X},\mathbb{Y}):=\inf\\{\mathcal{E}^{\wp}[h]:h\in\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})\\}.$

Assume that $\mathsf{H}_{\wp}(\mathbb{X},\mathbb{Y})$ is the class of
$\mathscr{W}^{1,2}$ homeomorphic mappings between $\mathbb{X}$ and
$\mathbb{Y}$ so that

$\int_{\mathbb{X}}\wp^{2}(f(z))(|f_{z}|^{2}+|f_{\bar{z}}|^{2})dxdy+\int_{\mathbb{X}}\wp^{2}(f(z))|f(z)|^{2}dxdy<\infty,$

and let $\overline{\mathsf{H}_{\wp}(\mathbb{X},\mathbb{Y})}$ be its closure in
$\mathscr{W}^{1,2}$. Further we have that
$\overline{\mathsf{H}_{\wp}(\mathbb{X},\mathbb{Y})}\subset\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$.

Now as in [6, Theorem 1.1.] is can be proved the following modification of the
previous proposition

######  Proposition 2.5.

[13, Lemma 2.15] Let $\mathbb{X}$ and $\mathbb{Y}$ be bounded doubly connected
planar domains. Assume that the boundary components of $\mathbb{X}$ do not
degenerate into points. There exists
$h\in\overline{\mathsf{H}_{\wp}(\mathbb{X},\mathbb{Y})}$ such that

${\mathcal{E}^{\wp}}[h]=\mathsf{E}^{\wp}(\mathbb{X},\mathbb{Y})=\inf\\{\mathcal{E}^{\wp}[h]:h\in\mathsf{H}_{\wp}(\mathbb{X},\mathbb{Y})\\}.$

###  2.1. Stationary mappings

We call a mapping $h\in\overline{\mathsf{H}_{\rho}(\mathbb{X},\mathbb{Y})}$
_stationary_ if

(2.1) $\frac{d}{dt}\bigg{|}_{t=0}{\mathcal{E}^{\wp}}[h\circ\phi_{t}^{-1}]=0$

for every family of diffeomorphisms
$t\to\phi_{t}\colon\mathbb{X}\to\mathbb{X}$ which depend smoothly on the
parameter $t\in\mathbb{R}$ and satisfy $\phi_{0}=\operatorname{id}$. The
latter mean that the mapping
$\mathbb{X}\times[0,\epsilon_{0}]\ni(t,z)\to\phi_{t}(z)\in\mathbb{X}$ is a
smooth mapping for some $\epsilon_{0}>0$. We now have.

######  Lemma 2.6.

[13] Let $\mathbb{X}=A(r,R)$ be a circular annulus, $0<r<R<\infty$, and
$\mathbb{Y}$ a doubly connected domain. If
$h\in\mathscr{W}^{\wp}(\mathbb{X},\mathbb{Y})$ is a stationary mapping, then

(2.2)
$\wp^{2}(h(z))h_{z}\overline{h_{\bar{z}}}\equiv\frac{c}{z^{2}}\qquad\text{in
}\mathbb{X}$

where $c\in\mathbb{R}$ is a constant.

Notice that the corresponding lemma in [13] is for deformations, but the proof
works as well for the class
$\overline{\mathsf{H}^{\wp}(\mathbb{X},\mathbb{Y})}\subset\mathfrak{D}^{\wp}(\mathbb{X},\mathbb{Y})$.

##  3\. Proof of the main result

We need the following important result concerning the local Lipschitz
character of certain mappings proved by Iwaniec, Kovalev and Onninen in [8].

######  Proposition 3.1.

Let $h\in\mathscr{W}^{1,2}(\mathbb{X})$ be a mapping with nonnegative
Jacobian. Suppose that the Hopf product $h_{z}\,\overline{h_{\bar{z}}}$ is
bounded and Hölder continuous. Then $h$ is locally Lipschitz but not
necessarily $\mathscr{C}^{1}$-smooth.

We need the also the following three results.

######  Lemma 3.2.

Every sense-preserving solution of Hopf equation so that $\mathrm{Hopf}(f)$ is
bounded that maps $\mathbb{A}(1,R)$ into $\mathbb{Y}$, mapping the inner/outer
boundary to inner/outer boundary is $(K,K^{\prime})$ quasiconformal, where

$K=1\ \ \text{and}\ \
K^{\prime}=\sup_{z\in\mathbb{A}(1,R)}\frac{|c|}{|z|^{2}\wp^{2}(f)}.$

###### Proof.

We have

$\begin{split}|Df(z)|^{2}&=|f_{z}|^{2}+|f_{\bar{z}}|^{2}\\\
&\leqslant(|f_{z}|^{2}-|f_{\bar{z}}|^{2})+2|f_{\bar{z}}|^{2}\\\
&\leqslant(|f_{z}|^{2}-|f_{\bar{z}}|^{2})+2|f_{\bar{z}}f_{z}|.\end{split}$

Since

$|f_{\bar{z}}f_{z}|=\frac{\mathrm{Hopf}(f)}{|\wp(f(z))|^{2}}=\frac{|c|}{|z|^{2}\wp^{2}(f)}\leqslant
K^{\prime}.$

This implies the claim. ∎

######  Proposition 3.3.

[4, Theorem 12.3] Let $K>1,K^{\prime}\geqslant 0$ and assume that
$f:\mathbb{X}\to\mathbf{C}$ is a $(K,K^{\prime})-$quasiconformal mapping so
that $|f(z)|\leqslant M$, $z\in\mathbb{X}$ and assume that
$\mathbb{X}^{\prime}\Subset\mathbb{X}$ and let
$d=\mathrm{dist}(\mathbb{X}^{\prime},\partial\mathbb{X})$. Then there is a
constant $C=C(K)$ so that

$|f(z)-f(z^{\prime})|\leqslant
C(K)(M+d\sqrt{K^{\prime}})|z-z^{\prime}|^{\beta},\
z,z^{\prime}\in\mathbb{X}^{\prime}$

where

$\beta=K-\sqrt{K^{2}-1}.$

If $K=1$, then the above theorem can be formulated for $K_{1}$ instead of $K$,
where $K_{1}>1$ is an arbitrary constant. For example for $K_{1}=5/4$, for
which we get $\beta=K_{1}-\sqrt{K_{1}^{2}-1}=1/2$.

We also need the following lemma

######  Lemma 3.4.

[13] Assume that $\mathbb{X}$ and $\mathbb{Y}$ are doubly connected domains,
and assume that
$a\colon\mathbb{X}\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\mathbb{A}(1,R)$
and
$b\colon\mathbb{A}(1,\rho)\overset{{}_{\textnormal{\tiny{onto}}}}{\longrightarrow}\mathbb{Y}$
are univalent conformal mappings and define
$\wp_{1}(w)={\wp(b(w))}{|b^{\prime}(w)|}$, $w\in\mathbb{A}(1,\rho)$. Then

  1. (a)

$\mathcal{E}^{\wp}[b\circ f\circ a]=\mathcal{E}^{\wp_{1}}[f]$ provided that
one of the two sides exist.

  2. (b)

$b\circ f\circ a\in D^{\wp}(\mathbb{X},\mathbb{Y})$ if and only if $f\in
D^{\wp_{1}}(\mathbb{A}(1,R),\mathbb{A}(1,\rho)))$.

  3. (c)

For Gauss curvature we have
$\mathcal{K}_{\wp}(b(w))=\mathcal{K}_{\wp_{1}}(w)$.

  4. (d)

$\wp$ is an allowable metric if and only if $\wp_{1}$ is an allowable metric.

  5. (e)

$W_{l}^{1,2}(\mathbb{A}(1,R),\mathbb{A}(1,\rho))=W^{1,2}(\mathbb{A}(1,R),\mathbb{A}(1,\rho)).$

  6. (f)

$b\circ f\circ a$ is $\wp$-harmonic if and only if $f$ is $\wp_{1}-$harmonic.

###### Proof of Theorem 1.1.

The existence part has been discussed in Proposition 2.5. We will apply a
self-improving argument. In view of Lemma 3.4, we can assume that
$\mathbb{X}=\mathbb{A}(1,R)$ for a constant $R>1$. Assume that
$f:\mathbb{A}(1,R)\to\mathbb{Y},$ so that

$\wp^{2}(f(z))f_{z}\bar{f}_{z}=\frac{c}{z^{2}}.$

Let $\Psi:\mathbb{Y}\to\mathbb{A}(1,\rho)$ be a conformal diffeomorphism and
let $\Phi=\Psi^{-1}$.

Then $F(z)=\Psi\circ f:\mathbb{A}(1,R)\to\mathbb{A}(1,\rho)$ is a minimizer of
$\wp_{1}-$ energy, between $\mathbb{A}(1,R)$ and $\mathbb{A}(1,\rho)$, where
$\wp_{1}(\zeta)={\wp(\Phi(\zeta))}{|\Phi^{\prime}(\zeta)|}$ and

(3.1) $\mathrm{Hopf}(F)(z)=\wp_{1}^{2}(F(z))F_{z}\bar{F}_{z}=\frac{c}{z^{2}}.$

Further for

$\tilde{F}(z)=\left\\{\begin{array}[]{ll}F(z),&\hbox{$1<|z|\leqslant R$;}\\\
\rho^{2}/\overline{F(R^{2}/\bar{z})},&\hbox{$R\leqslant|z|<R^{2}$}\\\
1/\overline{F(1/\bar{z})},&\hbox{$1/R\leqslant|z|\leqslant
1$},\end{array}\right.$

let

$\tilde{\wp}(w)=\left\\{\begin{array}[]{ll}\wp_{1}(w),&\hbox{$1\leqslant|w|\leqslant\rho$;}\\\
\wp_{1}(\rho^{2}/\bar{w}),&\hbox{$\rho<|w|<\rho^{2}$,}\\\
\wp_{1}(1/\bar{w}),&\hbox{$1/\rho<|w|\leqslant 1$.}\par\end{array}\right.$

Observe that $\tilde{F}$ is continuous. Namely, if $F(Re^{it})=\rho e^{is},$
then

$\rho^{2}/\overline{F(R^{2}/(Re^{-it}))}=\rho e^{is}.$

Thus $\tilde{F}:\mathbb{A}(1/R,R^{2})\to\mathbb{A}(1/\rho,\rho^{2})$ is
continuous and belongs to the same class as $F$, which mean it is in
$\mathscr{W}^{1,2}$. Then by direct calculation we get

$\mathrm{Hopf}(\tilde{F})=\tilde{\wp}^{2}(\tilde{F}(z))\tilde{F}_{z}\overline{\tilde{F}}_{z}=\left\\{\begin{array}[]{ll}\frac{c}{z^{2}},&\hbox{$1<|z|\leqslant
R$;}\\\
\frac{\rho^{4}}{|F(R^{2}/\bar{z})|^{4}}\frac{c}{z^{2}},&\hbox{$R\leqslant|z|<R^{2}$.}\\\
\frac{1}{|F(1/\bar{z})|^{4}}\frac{c}{z^{2}},&\hbox{$1/R\leqslant|z|<1$.}\end{array}\right.$

Let us demonstrate for example

$\mathrm{Hopf}(\tilde{F})=\frac{1}{|F(1/\bar{z})|^{4}}\frac{c}{z^{2}}$

for $1/R\leqslant|z|<1$.

We have

$\tilde{F}_{z}(z)=\frac{d}{dz}\left(\frac{1}{\overline{F}(\frac{1}{\bar{z}})}\right)=\frac{1}{\overline{F^{2}}(\frac{1}{\bar{z}})}\overline{F_{z}\left(\frac{1}{\bar{z}}\right)}\frac{1}{z^{2}}$

and

$\overline{\tilde{F}}_{z}(z)=\frac{d}{dz}\left(\frac{1}{{F}(\frac{1}{\bar{z}})}\right)=\frac{1}{{F^{2}}(\frac{1}{\bar{z}})}\overline{\overline{F}_{z}\left(\frac{1}{\bar{z}}\right)}\frac{1}{z^{2}}.$

Then from (3.1) we get

$\tilde{\wp}^{2}(\tilde{F}(z))\tilde{F}_{z}(z)\overline{\tilde{F}}_{z}(z)=\frac{1}{|F(\frac{1}{\bar{z}})|^{4}}\frac{cz^{2}}{z^{4}}=\frac{1}{|F(\frac{1}{\bar{z}})|^{4}}\frac{c}{z^{2}}.$

From Kellogg theorem for conformal mappings we know that $\Phi^{\prime}$ is
$\mathscr{C}^{\alpha}$ up to the boundary of $\mathbb{Y}$. This implies in
particular that the function

$\tilde{F}_{z}\overline{\tilde{F}}_{z}=\frac{\mathrm{Hopf}(\tilde{F})}{\tilde{\wp}^{2}(\tilde{F}(z))}$

is bounded in $\mathbb{A}(1/R,R^{2})$. From Lemma 3.2, in view of Lemma 2.3,
we obtain that $\tilde{F}$ is $(1,\tilde{K})$-quasiconformal. Now Proposition
3.3 implies that $F$ is $\frac{1}{2}-$Hölder continuous in $\mathbb{A}(1,R).$
More precisely

$|\tilde{F}(z)-\tilde{F}(z^{\prime})|\leqslant
C(1)\left(\rho^{2}+\frac{R^{2}-1}{R^{2}}\sqrt{\tilde{K}}\right)|z-z^{\prime}|^{1/2}.$

So

$|F(z)-F(z^{\prime})|\leqslant
C(1)\left(\rho^{2}+\frac{R^{2}-1}{R^{2}}\sqrt{\tilde{K}}\right)|z-z^{\prime}|^{1/2}.$

This implies that $\tilde{F}_{z}\overline{\tilde{F}}_{z}$ is $\alpha\cdot
1/2-$Hölder continuous as a composition of two mappings which are $\alpha-$
and $1/2-$Hölder continuous respectively in $\mathbb{A}(1/R,R^{2})$. Now
Theorem 3.1 implies that $\tilde{F}$ is locally Lipschitz. In particular $F$
is Lipschitz in $\overline{\mathbb{A}(1,R)}\Subset\mathbb{A}(1/R,R^{2})$.
Therefore $f:\mathbb{A}(1,R)\to\mathbb{Y}$ is Lipschitz continuous. Now Lemma
3.4 and Kellogg theorem for conformal mappings imply the claim. This finishes
the proof of the main result. ∎

## References

  * [1] L. V. Ahlfors, Lectures on quasiconformal mappings, Second edition. With supplemental chapters by C. J. Earle, I. Kra, M. Shishikura and J. H. Hubbard. University Lecture Series, 38. American Mathematical Society, Providence, RI, 2006.
  * [2] L. C. Evans _A new proof of local $C^{1,\alpha}$ regularity for solutions of certain degenerate elliptic P.D.E._ J. Differ. Equations 45, 356–373 (1982).
  * [3] K. Astala, T. Iwaniec, and G. Martin, Elliptic partial differential equations and quasiconformal mappings in the plane, Princeton University Press, Princeton, NJ, 2009.
  * [4] D. Gilbarg and N. Trudinger: _Elliptic Partial Differential Equations of Second Order_. 2 Edition, Springer 1977, 1983.
  * [5] E. Heinz: _On certain nonlinear elliptic differential equations and univalent mappings._ J. Anal. Math. 5, 197–272 (1957).
  * [6] T. Iwaniec, Tadeusz, J. Onninen _Mappings of least Dirichlet energy and their Hopf differentials._ Arch. Ration. Mech. Anal. 209, No. 2, 401–453 (2013).
  * [7] T. Iwaniec, N.-T. Koh, L.V. Kovalev, J. Onninen: _Existence of energy-minimal diffeomorphisms between doubly connected domains._ Invent. Math. 186, 667–707 (2011).
  * [8] T. Iwaniec, Tadeusz, J. Onninen, L. Kovalev: _Lipschitz regularity for inner-variational equations._ Duke Math. J. 162, No. 4, 643–672 (2013).
  * [9] J. Jost, Harmonic maps between surfaces. Lecture Notes in Mathematics, 1062\. Springer-Verlag, Berlin, 1984. x+133 pp.
  * [10] J. Jost, _Minimal surfaces and Teichmüller theory_. Yau, Shing-Tung (ed.), Tsing Hua lectures on geometry and analysis, Taiwan, 1990-91. Cambridge, MA: International Press. 149-211 (1997).
  * [11] D. Kalaj: _Lipschitz property of minimisers between double connected surfaces_. J. Geom. Anal. 30, 4150–4165(2020)
  * [12] D. Kalaj: _Kellogg’s theorem for diffeomophic minimisers of Dirichlet energy between doubly connected Riemann surfaces_ , arXiv:2011.04629.
  * [13] D. Kalaj: _Energy-minimal diffeomorphisms between doubly connected Riemann surfaces._ Calc. Var. Partial Differ. Equ.51 , 465–494 (2014).
  * [14] D. Kalaj, B. Lamel: _Minimisers and Kellogg’s theorem._ Math. Ann. 377, No. 3-4, 1643–1672 (2020).