# A unified linear intrinsic alignment model for elliptical and disc galaxies
and the resulting ellipticity spectra

Basundhara Ghosh1,⋆ Kai Nußbaumer2 Eileen Sophie Giesel2 Björn Malte
Schäfer2,♯ 1 Department of Physics, Indian Institute of Science, Bengaluru
560012, Karnataka, India
2Zentrum für Astronomie der Universität Heidelberg, Astronomisches Rechen-
Institut, Philosophenweg 12, 69120 Heidelberg, Germany

###### Abstract

Alignments of disc galaxies were thought to result from tidal torquing, where
tidal field of the cosmic large-scale structure exert torquing moments onto
dark matter haloes, determining their angular momentum and ultimately the
orientation of galactic discs. In this model, resulting intrinsic ellipticity
correlations are typically present on small scales, but neither observations
nor simulations have found empirical evidence; instead, simulations point at
the possibility that alignments of disc galaxies follow a similar alignment
model as elliptical galaxies, but with a weaker alignment amplitude. In our
article we make the case for the theory of linear alignments resulting from
tidal distortions of the galactic disc, investigate the physical properties of
this model and derive the resulting angular ellipticity spectra, as they would
appear as a contribution to weak gravitational lensing in surveys such as
Euclid’s. We discuss in detail on the statistical and physical properties of
tidally induced alignments in disc galaxies, as they are relevant for
mitigation of alignment contamination in weak lensing data, comment on the
consistency between the alignment amplitude in spiral and elliptical galaxies
and finally, estimate their observability with a Euclid-like survey.

###### keywords:

galaxies: intrinsic alignments, gravitational lensing, tidal fields, cosmic
large-scale structure

sf

††thanks: ⋆<EMAIL_ADDRESS>♯ bjoern.malte.schaefer@uni-
heidelberg.de

##  1 Introduction

The subject of galaxy alignments has been of prime importance and the
correspondence of linear alignment models to elliptical galaxies and quadratic
alignment models to spiral galaxies is well known (Schaefer, 2009; Joachimi et
al., 2015; Kiessling et al., 2015; Kirk et al., 2015; Troxel and Ishak, 2014;
Lamman et al., 2023). In general, any dark matter halo that contains galaxies
undergoes a distortion due to gravitational tidal fields when embedded in the
large-scale structure. The stellar component is oriented according to the
perturbation of its dynamical equilibrium, and this orientation would depend
on the alignment model being followed. If the tidal shear is small enough, the
observed shape or the degree of ellipticity can be expected to have a linear
relationship with the magnitude of the tidal field, described by a single
alignment parameter (Hirata and Seljak, 2004). In the particular case of
spiral galaxies, the galactic disc has an ellipticity that depends on its
angle of inclination relative to the line of sight of the observer, which is
determined by its angular momentum. Unlike in the case of elliptical galaxies
where ellipticity is a measure of the shape and depends on the strength of the
tidal shear, the ellipticity now is determined by the direction of the angular
momentum. Although the relation between the angular momentum and the
surrounding large-scale structure is not well-understood, one of the widely
accepted mechanisms is that of tidal torquing, where even before a galaxy is
formed, a torquing moment is exerted on the protohalo by the surrounding
gravitational field (Crittenden et al., 2001). This imparts a default angular
momentum to the host halo, and as the galactic disc forms later, it
incorporates this angular momentum. In this kind of alignment, the tidal field
in its lowest order enters in a quadratic form, hence the name quadratic
alignment model. While linear tidal torquing theory is quite successful, the
non-linear regime is not easy to describe analytically. The linear alignment
model, on the other hand, is well-described even if the large-scale structure
enters a nonlinear regime (nonlinear alignment model, see Hirata et al., 2007;
Bridle and King, 2007). There have been studies to analyse the observability
of intrinsic alignments of spiral galaxies. Crittenden et al. (2001); Schaefer
and Merkel (2012); Chisari et al. (2015); Rodriguez-Gomez et al. (2017) derive
ellipticity spectra for angular momentum-induced alignments. They should
generate rather strong signals for Euclid-like (Amendola et al., 2018) surveys
(Schäfer and Merkel, 2017). Combining the alignment models for the two major
galaxy types would lead to a rich phenomenology and the question of separating
the two mechanisms (Tugendhat and Schäfer, 2018; Reischke and Schäfer, 2019).

Tidal shearing and tidal torquing have been recently tested in Zjupa et al.
(2022); Delgado et al. (2023) using morphologically selected elliptical and
spiral galaxies from the IllustrisTNG simulation (Nelson et al., 2018) in the
range of redshifts $0\leq z<1$. In this work it was found that elliptical
galaxies show strong alignment signals further increasing with mass and
redshift. Contrary to that, spiral galaxies exhibit a strong alignment only at
high redshifts near $z=1$, and their alignment follows rather the tidal
shearing instead of the tidal torquing model. Weak lensing surveys take a
similar attitude, reporting a detection of a ellipticity cross-correlation
between lensing and the shape of elliptical galaxies as a signature of tidal
shearing. They place bounds on the cross-correlation of lensing and the shapes
of spiral galaxies, implying that their alignment must be much weaker
(Hildebrandt et al., 2017; Johnston et al., 2019; Samuroff et al., 2023),
casting doubts on the validity of the quadratic alignment model.

In this work, we propose a unified alignment model that can cater to both
elliptical and spiral galaxies: If a such a model naturally accounts for the
much weaker alignments of spirals, it would not be in conflict with
observations and could account for the results from simulations. In
particular, we investigate the possibility of a unified linear alignment model
for disc galaxies, that follows an approach analogous to our computations for
the ellipticity of elliptical galaxies in Ghosh et al. (2021). The underlying
principle is that disc galaxies subjected to a gravitational tidal field
undergo a distortion effect that can cause them to assume an elliptical shape,
resulting in a change of the circular orbit of a star at the edge of a galaxy
to an elliptical orbit. The star, which is bound by the NFW halo potential,
experiences tidal forces, and the ellipticity of the galaxy scales with the
tidal strength. An important ingredient for computing the ellipticity of
galaxies is the alignment parameter which depends on the scale at which
alignment due to tidal shear takes place. In case of elliptical galaxies, this
scale can be assumed to be that of the Sérsic radius, following a de
Vaucouleur’s profile with Sérsic index $n=4$, whereas for spiral galaxies, it
should be typically the edge of the galactic disc. We perform a perturbative
calculation to compute the effect of a quadrupolar tidal gravitational field
on a circular orbit, which naturally generates an ellipticity proportional to
the magnitude of the tidal field. We then simulate the orbit of a star in an
NFW-halo, which is tidally perturbed again by a quadrupolar field, and this in
turn is in agreement with our perturbative results. We eventually see that the
strength of alignment for disc galaxies is smaller than that of elliptical
galaxies. In following this idea, we want to provide a physical idea behind
linear alignments for disc galaxies as observed in simulations (Hilbert et
al., 2017; Jagvaral et al., 2022) and overpredicted amplitude of tidal
torquing models which is in disagreement with observations (Heymans et al.,
2013; Kilbinger et al., 2009), as well as bridge to effective alignment models
(Blazek et al., 2019, 2011, 2015). With a unified model for alignments of disc
and elliptical galaxies, there would be no need to differentiate between the
two types, leading to simplified effective models.

Typically, we present results for a spatially flat standard $\Lambda$CDM-
cosmology with conservative parameter choices: $\Omega_{m}=0.3$ (with a
baryonic contribution of $\Omega_{b}=0.04$), $\Omega_{\Lambda}=1-\Omega_{m}$,
$n_{s}=0.96$, $\sigma_{8}=0.8$ and a value of $h=0.7$. We assume an implicit
summation over repeated indices. After introducing the linear alignment model
for disc galaxies, investigating its physical properties and juxtaposing it
with an analogous alignment model for elliptical galaxies in Sect. 2, we
compute in Sect. 3 the resulting angular ellipticity spectra for a Euclid-like
survey with a redshift distribution $n(z)\mathrm{d}z\propto
z^{2}\exp(-z/z_{0})$ with a median redshift of 0.9. A summary and discussion
in Sect. 4 concludes our investigation.

##  2 Linear alignments for disc galaxies

###  2.1 Scaling arguments

Tidal fields introduce a distortion of an elliptical galaxy by differential
gravitational acceleration. Assuming that the galaxy is in virial equilibrium
with a velocity dispersion $\sigma^{2}$, the ellipticity is given by the
second derivatives of the gravitational potential in units of the velocity
dispersion,

$\epsilon_{\mathrm{el}}=r_{\mathrm{vir}}^{2}\partial^{2}\frac{\Phi}{\sigma_{\mathrm{vir}}^{2}}.$
(1)

As the second moments of the brightness distributions are computed out to
distances of $r_{\mathrm{vir}}$, the constant of proportionality, i.e. the
alignment parameter $D_{\mathrm{IA}}$, has units of area and should be related
to the squared size of the galaxy. One can compute dimensionless correction
factors to this relation for brightness profiles of the Sérsic-family,
resulting in small correction for exponential profiles, but rather large
influences in the case of de Vaucouleurs-profiles (Ghosh et al., 2021; Giesel
et al., 2022a). There should be a similar linear effect of tidal fields on the
shape of a disc galaxy: The galactic disc extends typically to the NFW-scale
radius $r_{s}$, which is a fraction $1/c(M_{\mathrm{vir}},z)$ of the virial
radius $r_{\mathrm{vir}}$, giving rise to the phenomenon of flat rotation
curves. The quantity $c(M_{\mathrm{vir}},z)$ is the concentration parameter
which relates the virial mass $M_{\mathrm{vir}}$ of a halo to its redshift
$z$. A quadrupolar tidal field distorts the edge of the galactic disc,
compressing it in one direction and elongating it in the perpendicular
direction, making it elliptical. The same scaling argument as in eqn. (1)
suggests a relation of the form

$\epsilon_{\mathrm{sp}}=r_{s}^{2}\partial^{2}\frac{\Phi}{\upsilon_{\mathrm{circ}}^{2}}$
(2)

where $r_{s}$ denotes the edge of the galactic disc.
$\upsilon_{\mathrm{circ}}$ as the circular velocity replaces $\sigma$ as the
velocity dispersion, describing elasticity or susceptibility of the stellar
component to tidal forces. To complete the scaling argument, we would like to
state that gravitational lensing shear introduces an ellipticity

$\gamma\propto 2\partial^{2}\frac{\Phi}{c^{2}}.$ (3)

The fact that tidal gravitational fields give rise to both lensing and
intrinsic alignments causes serious problems in the analysis of lensing data.
Whereas the stars as test particles of the tidal field determine the magnitude
of the observable effect by either $\sigma^{2}$ or
$\upsilon_{\mathrm{circ}}^{2}$ in intrinsic alignments, the velocity scale for
photons as test particles is of course $c^{2}$.

###  2.2 Perturbative calculation

From scaling arguments alone one would expect an effect of tidal gravitational
fields on the shape of an otherwise circular galactic disc. Whether the effect
is real and present at linear order, can be shown in a straightforward
perturbative calculation: The circular motion at the edge of the galactic disk
is a superposition of two harmonic oscillations in perpendicular directions
with the same angular frequency $\omega=\upsilon_{\mathrm{circ}}/r_{s}$ and a
phase shift of $\pi/2$. The effective potential for motion in the $x,y$-plane
consists of a harmonic potential responsible for the circular motion,
perturbed by a quadrupole with fixed amplitude:

$\Phi=\frac{\omega^{2}}{2}(x^{2}+y^{2})+\frac{\Omega^{2}}{2}(x^{2}-y^{2})$ (4)

where the amplitude of the quadrupolar tidal field is controlled by
$\Omega\ll\omega$. Effectively, this relation encapsulates a single component
of the general tidal perturbation $\partial_{i}\partial_{j}\Phi$ perturbing
the potential $\Phi$

$\Phi\rightarrow\Phi+\frac{1}{2}r_{i}r_{j}\>\partial_{i}\partial_{j}\Phi$ (5)

where we decompose the tensors $r_{i}r_{j}$ and $\partial_{i}\partial_{j}\Phi$
as

$r_{i}r_{j}=(x^{2}+y^{2})\sigma^{(0)}_{ij}+(x^{2}-y^{2})\sigma^{(1)}_{ij}+2xy\sigma^{(3)}_{ij}$
(6)

as well as

$\partial_{i}\partial_{j}\Phi=\Delta\Phi\sigma^{(0)}_{ij}+(\partial^{2}_{x}-\partial^{2}_{y})\Phi\sigma^{(1)}_{ij}+2\partial_{x}\partial_{j}\Phi\sigma^{(3)}_{ij}$
(7)

with the real-valued Pauli-matrices $\sigma^{(n)}$, such that the contraction
$r_{i}r_{j}\partial_{i}\partial_{j}\Phi$ becomes

$r_{i}r_{j}\partial_{i}\partial_{j}\Phi=(x^{2}+y^{2})\Delta\Phi+(x^{2}-y^{2})(\partial^{2}_{x}-\partial^{2}_{y})\Phi+2xy\partial_{x}\partial_{y}\Phi$
(8)

We focus on the second term in eqn. (4) with the tidal field amplitude
$\Omega^{2}=(\partial^{2}_{x}-\partial^{2}_{y})\Phi$. Interestingly, both the
unperturbed circular motion and the fixed tidal perturbation are mapped onto
harmonic potentials, albeit anisotropic in the second case. Considering a
single oscillation,

$x(t)=r_{s}\sin(\omega t)$ (9)

along the $x$-direction one obtains for the perturbed amplitude $\Delta(t)$
the equation of motion

$\ddot{\Delta}(t)=-\frac{\mathrm{d}\Phi}{\mathrm{d}x}\quad\mathrm{with}\quad\frac{\mathrm{d}\Phi}{\mathrm{d}x}=(\omega^{2}+\Omega^{2})x$
(10)

which is evaluated in the integration along the unperturbed trajectory $x(t)$,
effectively in a Born-type approximation:

$\ddot{\Delta}(t)=-(\omega^{2}+\Omega^{2})r_{s}\sin(\omega t)$ (11)

Integrating twice yields a perturbative solution to the trajectory. In
particular, evaluating $\Delta(t)$ at values 0 and $\pi/2$ for the phase
$\omega t$ corresponds to distortions $\propto 1\pm\Omega^{2}/\omega^{2}$ of
the semi-major and semi-minor axes, resulting in the ellipticity

$\epsilon=\frac{(1+\Omega^{2}/\omega^{2})-(1-\Omega^{2}/\omega^{2})}{(1+\Omega^{2}/\omega^{2})+(1-\Omega^{2}/\omega^{2})}=\frac{\Omega^{2}}{\omega^{2}}=\frac{r_{s}^{2}}{\upsilon_{\mathrm{circ}}^{2}}\partial^{2}\Phi$
(12)

as expected from scaling arguments.

An approximate energy consideration leads to a qualitatively identical result:
As the potential is perturbed by the tidal field according to the replacement

$\Phi=\frac{\omega^{2}}{2}x^{2}\quad\rightarrow\quad\Phi=\frac{\omega^{2}\pm\Omega^{2}}{2}x^{2}$
(13)

the amplitudes in $x$\- and $y$-direction that are reachable with the
potential energy $V=\omega^{2}r_{s}^{2}/2$ are changed to

$\Delta_{\pm}=r_{s}\sqrt{\frac{\omega^{2}}{\omega^{2}+\Omega^{2}}}\simeq
r_{s}\left(1\mp\frac{\Omega^{2}}{2\omega^{2}}\right)$ (14)

resulting again in an ellipticity $\epsilon$

$\epsilon=\frac{\Omega^{2}}{2\omega^{2}}=\frac{1}{2}\frac{r_{s}^{2}}{\upsilon_{\mathrm{circ}}^{2}}\partial^{2}\Phi$
(15)

The factor 2 appearing in comparison to the above argument can be traced back
to the virial theorem.

Yet another interpretation of the process would be to associate a time scale
to the potential of the halo, either the free-fall time scale
$t_{\mathrm{el}}=2\pi/\omega\simeq r_{\mathrm{vir}}/\sigma$ in the case of
elliptical galaxies, and the orbital time scale
$t_{\mathrm{sp}}=2\pi/\omega\simeq r_{s}/\upsilon_{\mathrm{circ}}$ for disc
galaxies, and set this in relation with the time scale that can be attributed
to the (quadrupolar) tidal field, $t_{\mathrm{tid}}=2\pi/\Omega$, as it is a
superposition of two harmonic potentials, with $\Omega^{2}$ reflecting
$\partial^{2}\Phi$. Then, the ellipticities are given by

$\epsilon_{\mathrm{el}}=\frac{t_{\mathrm{el}}^{2}}{t_{\mathrm{tid}}^{2}}\quad\mathrm{and}\quad\epsilon_{\mathrm{sp}}=\frac{t_{\mathrm{sp}}^{2}}{t_{\mathrm{tid}}^{2}}$
(16)

implying that the relative alignment strength in the same tidal field should
be
$\epsilon_{\mathrm{sp}}=t_{\mathrm{sp}}^{2}/t_{\mathrm{el}}^{2}\times\epsilon_{\mathrm{el}}$:
The ratio between free-fall time scale and orbital time scale should determine
the relative magnitude of the alignments.

###  2.3 Direct numerical simulation

The perturbative calculation neglects that a perturbed orbit experiences
changing gradients of the NFW-potential of the halo, and that the orbits can
not be closed as Bertrand’s theorem is violated. In particular the first
effect is potentially interesting: Density profiles $\propto 1/r^{2}$
naturally give rise to flat rotation curves, as equating the centrifugal
acceleration and the potential gradient according to
$\upsilon_{\mathrm{circ}}^{2}/r=\mathrm{d}\Phi/\mathrm{d}r$ suggests constant
$\upsilon_{\mathrm{circ}}$. While this is approximately true for NFW-profiles
around $r_{s}$, the calculation shows decreasing velocities in radial
direction. An intuitive interpretation of the violation of Bertrand’s theorem
in the case of a perturbed potential is that the effective frequencies
$\sqrt{\omega^{2}\pm\Omega^{2}}$ in the two directions are not identical
giving rise to a Lissajous figure. But for weak tidal fields,
$\Omega^{2}\ll\omega^{2}$, the effect is arbitrarily small.

A direct simulation of the orbit of stars as test particles at the edge of a
galactic disc with radius $r_{s}$ in a NFW-potential that is perturbed by a
quadrupolar tidal field confirms that there is a linear reaction of the shape
of the galactic disc and yields numerical results that are well in agreement
with the scaling and perturbative arguments. Using dimensionless variables
$x=r/r_{s}$ and $\tau=\omega t$ gives an equation of motion

$\frac{\mathrm{d}^{2}x}{\mathrm{d}\tau^{2}}=-\frac{\mathrm{d}}{\mathrm{d}x}\frac{\Phi}{\upsilon_{\mathrm{circ}}^{2}}\quad\text{with}\quad\upsilon_{\mathrm{circ}}^{2}=(\omega
r_{s})^{2}.$ (17)

Specifically, an NFW-potential would have potential gradients

$\frac{\mathrm{d}}{\mathrm{d}x}\Phi=\frac{GM_{\mathrm{vir}}}{Ar_{s}}\frac{x/(1+x)-\ln(1+x)}{x^{2}}$
(18)

which determines the circular velocity at $x=1$:

$\omega^{2}=\frac{GM_{\mathrm{vir}}}{r_{\mathrm{vir}}^{3}}\frac{B}{A}$ (19)

and therefore the amplitude of the potential, abbreviating

$B=1/2-\ln 2$

and

$A=\ln(1+c(M_{\mathrm{vir}},z))-\frac{c(M_{\mathrm{vir}},z)}{(1+c(M_{\mathrm{vir}},z))}$

Summarising, the equation of motion is given explicitly by

$\frac{\mathrm{d}^{2}x}{\mathrm{d}\tau^{2}}=-\frac{1}{B}\frac{x/(1+x)-\ln(1+x)}{x^{2}},$
(20)

which we solve numerically for a circular orbit at $x=1$ and
$\upsilon_{\mathrm{circ}}=1$, perturbed by tidal shear fields
$\Omega^{2}/2\>(x^{2}-y^{2})$ and $\Omega xy$.

Fig. 1 shows the effect of tidal fields onto a circular orbit at $r_{s}$,
specifically, the distortion into an ellipse caused by an anisotropic tidal
shear field. An estimate of the ellipticity $\epsilon$ from the semi-major and
semi-minor axes is exactly the tidal shear field, in our dimensionless choice
of coordinates: $\epsilon\simeq\Omega^{2}/\omega^{2}$ with $r_{s}=1$ and
$\upsilon_{\mathrm{circ}}=1$ such that $\omega=1/(2\pi)$. In fact, Fig. 1
suggests an ellipticity of $\simeq 0.1$ for a tidal shear field
$\Omega^{2}=0.04$, consistent with the expectation. Isotropic tidal shear
fields would cause an effect on the size of the galactic disc, which we
comment on in the Appendix.

Figure 1: Perturbed orbits (green) in comparison to unperturbed circular
orbits (blue) at $r=r_{s}$ (or $x=1$) in the potential of an NFW-halo with a
tidal field $\Omega^{2}(x^{2}-y^{2})/2$ leading to $\epsilon_{+}$ (top panel),
and a tidal field $\Omega^{2}xy$ leading to $\epsilon_{\times}$ (bottom
panel). The magnitude of the perturbing field is $\Omega^{2}=0.04$ in
dimensionless units.

###  2.4 Alignment parameters: numerical values and scaling

####  2.4.1 Perturbed dynamics in disc and elliptical galaxies

In essence, one would expect an alignment of disc galaxies following a linear
alignment model and the relation between the ellipticities of a disc galaxy
and an elliptical galaxy in the same tidal gravitational field to be

$\epsilon_{\mathrm{sp}}=\frac{r_{s}^{2}}{r_{\mathrm{vir}}^{2}}\frac{\sigma_{\mathrm{vir}}^{2}}{\upsilon_{\mathrm{circ}}^{2}}\times\epsilon_{\mathrm{el}}.$
(21)

Naively, the much smaller size $r_{s}$ of the galactic disc in comparison to
the virial radius $r_{\mathrm{vir}}$ would be responsible for a much weaker
alignment of disc galaxies relative to ellipticals. The alignment parameter
for discs should be smaller by a factor of $1/c(M_{\mathrm{vir}},z)^{2}$
relative to the alignment of ellipticals, with $c$ being the NFW-concentration
parameter. In fact,
$r_{s}^{2}/r_{\mathrm{vir}}^{2}=1/c(M_{\mathrm{vir}},z)^{2}\simeq 0.02\ldots
0.06$, at redshifts around unity, where Euclid’s galaxy distribution is
peaking, and for typical masses $10^{11\ldots 12}M_{\odot}/h$. The scaling of
$c(M_{\mathrm{vir}},z)$ with redshift and halo mass follows the numerical
model by Child et al. (2018), who condensed results from numerical simulations
into phenomenological relations summarised by Fig. 3. On the other hand, the
squared velocity ratio in the second term typically assumes values close to
one, specifically
$\sigma_{\mathrm{vir}}^{2}/\upsilon_{\mathrm{circ}}^{2}\simeq 0.9\ldots 0.95$
at a redshift of $z=1$ and for galaxy masses of $M=10^{11\ldots
12}M_{\odot}/h$ without a strong variability. Fig. 3 shows as well how the
ratio between circular velocity $\upsilon_{\mathrm{circ}}$ at the scale radius
$r_{s}$ and the velocity dispersion $\sigma_{\mathrm{vir}}^{2}$

$\frac{\upsilon_{\mathrm{circ}}^{2}}{\sigma_{\mathrm{vir}}^{2}}=\frac{c(M_{\mathrm{vir}},z)(\ln
2-1/2)}{\ln(c(M_{\mathrm{vir}},z)+1)-c(M_{\mathrm{vir}},z)/(c(M_{\mathrm{vir}},z)+1)},$
(22)

or more compactly,
$\upsilon_{\mathrm{circ}}^{2}/\sigma_{\mathrm{vir}}^{2}=c(M_{\mathrm{vir}},z)B/A$.
The ratio changes as a function of halo mass $M$ and redshift $z$ through the
concentration parameter $c(M_{\mathrm{vir}},z)$ according to the model
introduced above.

It is quite important to realise that there seems to be a fundamental
difference in the linear alignment of elliptical and spiral galaxies: While
the reaction of an elliptical galaxy as a virialised sphere to an external
field is homologous, the induced ellipticity of the orbit of a star is a
function of radius, as illustrated by Fig. 2. There is in fact a strong
dependence of ellipticity on radius, depicted in the range
$x=r/r_{s}=1/3\ldots 2$ caused by the fact that stars on more distant orbits
are less tightly bound and therefore more susceptible to tidal distortion. We
chose specifically this interval of radii to cover the range between the
Sérsic-radius on small scales up to the edge of the galactic disc, where the
circular velocity reaches the highest values, while emphasising that at even
larger radii the tidal effect would be stronger.

Figure 2: Scaling of the ellipticity $\epsilon$ with orbital radius given in
dimensionless units where $r=r_{s}$ corresponds to $x=1$, for a fixed tidal
shear field $\Omega^{2}=0.01,0.02,0.03,0.04$ (from bottom to top). The inset
shows the dimensionless rotation curve $\upsilon_{\mathrm{circ}}(x)$ for the
range of radii $x=0.3\ldots 2$ considered here.

From this point of view it becomes necessary to define a mean induced
ellipticity, weighted with the exponential profile $\propto\exp(-r/R)$ with
the Sérsic-radius $R$:

$\langle\epsilon\rangle=\frac{\int 2\pi
r\mathrm{d}r\>\epsilon(r)\exp(-r/R)}{\int 2\pi r\mathrm{d}r\>\exp(-r/R)},$
(23)

which leads with an assumed exponential dependence
$\epsilon(r)\propto\exp(+r/Q)$ suggested by Fig. 2 with a scale $Q$ to

$\langle\epsilon\rangle=\frac{1}{\left(1-R/Q\right)^{2}}$ (24)

which can in principle lead to large numbers if $Q\simeq R$, as the integral
diverges for this choice. Given the numerical results of Fig. 2, the averaging
will yield larger values for $\langle\epsilon\rangle$ compared to the value
for $\epsilon$ at $x=1$ by a few: $Q$ assumes values around unity, so that
$R/Q\simeq 1/3$, and consequently, $\langle\epsilon\rangle\simeq 9/4=2.25$.
Looking back at eqn. (21), this number corrects the expectation for
$\epsilon_{\mathrm{sp}}$ in relation to $\epsilon_{\mathrm{el}}$ from
$1/c(M_{\mathrm{vir}},z)^{2}\simeq 0.04$ to about 0.1, in agreement with the
results by Hilbert et al. (2017) found in a fluid-mechanical simulation of
galaxy formation and evolution in terms of the relative amplitudes of the
intrinsic alignment spectra separated by early- and late-type galaxies. More
recently, simulations of Jagvaral et al. (2022) find alignments of disc
galaxies that are not as weak as those reported by Hilbert et al. (2017),
which our model can accommodate by extending the radial averaging to larger
distances, with a stronger response of the orbits to perturbing tidal fields.
In summary, our analytical model does in fact explain ellipticities of disc
galaxies that scale linearly with tidal field, and predicts the correct
magnitude of the alignment. In the following, we work with a relative
alignment parameter of $4$, intermediate between Hilbert et al. (2017) and
Jagvaral et al. (2022).

####  2.4.2 Tidal effect on moments of the brightness distribution

Ellipticity as the actual observable pertains to the size of luminous
component of the galaxy and the brightness profile $\rho_{0}(r)$, for which we
assumed in previous works a Sérsic profile with characteristic radius $R$. The
perturbation $\Delta Q_{ab}$ of the brightness moments $Q_{ab}$ of an
elliptical galaxy can be computed with

$\Delta
Q_{ab}=\frac{1}{2\sigma_{\mathrm{vir}}^{2}}\frac{\int\mathrm{d}^{2}r\>\rho_{0}(r)\>r_{a}r_{b}\>r_{i}r_{j}}{\int\mathrm{d}^{2}r\>\rho_{0}(r)}\>\partial^{2}_{ij}\Phi=S_{abij}\>\partial^{2}_{ij}\Phi.$
(25)

From the second moments $Q_{ij}$ of the brightness distribution one defines
the complex ellipticity $\epsilon$,

$\epsilon\simeq\frac{\Delta Q_{xx}-\Delta
Q_{yy}}{Q_{xx}+Q_{yy}}+2\mathrm{i}\frac{\Delta Q_{xy}}{Q_{xx}+Q_{yy}}$ (26)

which directly reflects magnitude and orientation of the tidal shear field
$\partial^{2}_{ij}\Phi$, with the susceptibility $S_{abij}$ as the constant of
proportionality: It includes the virial velocity $\sigma_{\mathrm{vir}}^{2}$,
and the details of the brightness distribution $\rho_{0}$. For an elliptical,
$S_{abij}$ is fundamentally an area of size $R^{2}$ for a Sérsic profile with
scale $R$, and a numerical factor resulting from the integration that depends
strongly on Sérsic-index $n$ (Ghosh et al., 2021; Giesel et al., 2022a, b):
Typical numerical prefactors to $R^{2}$ resulting from the integration of the
Sérsic-profile are $S\simeq 10$ for exponential profiles with $n=1$ up to
$S\simeq 10^{3}$ for de Vaucouleurs-profiles with $n=4$. The reaction of a
galactic disc to a tidal shear field can in principle be encapsulated by the
same relation, with $\upsilon_{\mathrm{circ}}^{2}$ replacing
$\sigma_{\mathrm{vir}}^{2}$, and absorbing the averaging in eqn. (24) into
$S_{abij}$.

####  2.4.3 Consistency relations

For spiral galaxies we expect an identical argument to hold, with the choice
of an exponential profile typical for galactic discs, and a radius-averaged
ellipticity, as deduced in the previous chapter. But there are two major
differences: Tidal alignments of orbits in the galactic disc are strongly
radius-dependent, necessitating the averaging process outlined in eqn. (24).
Secondly, the squared Sérsic-radius, which is a sensible choice for the
alignment parameter of elliptical galaxies as it reproduces the alignment
amplitude and corresponds to the physical picture of a tidally distorted
virialised object, is remarkably similar for spiral and elliptical galaxies of
similar mass. Therefore, it cannot be the alignment parameter for a disc. If
one stipulates that $SR^{2}=r_{s}^{2}$ for disc galaxies and
$SR^{2}=r_{\mathrm{vir}}^{2}$ for ellipticial galaxies, the missing link
between the dynamical properties of the dark matter halo (encapsulated in
$r_{s}$ and $r_{\mathrm{vir}}$) and the light distribution of the luminous
component (described by $R$ and $n$) can potentially be established:
Numerically, this would imply for a Milky Way-sized halo of $\sim
10^{12}M_{\odot}/h$ with a virial radius of $\sim 100~{}\mathrm{kpc}/h$
(Dehnen et al., 2006), an NFW-scale radius $r_{s}\sim 10~{}\mathrm{kpc}/h$
(Lin and Li, 2019) and a typical Sérsic-radius of a few kpc (Goodwin et al.,
1997, 1998), which would not be too dissimilar for an elliptical of the same
mass.

####  2.4.4 Observational consequences

Through our scaling argument, one would assume the reaction of disc galaxies
to be very similar compared to that of elliptical galaxies. In addition,
because of the model’s linearity, there would be $(i)$ a $GI$-type cross-
correlation between the intrinsic shapes of disc galaxies and gravitational
lensing with the characteristic negative amplitude, as well as $(ii)$ an
$II$-type correlation between the intrinsic shapes of disc and elliptical
galaxies, with alignment parameters that follow consistently from the linear
tidal interaction model. Additionally, $(iii)$ one can speculate about an
effect of $\Delta\Phi$ as the trace of the tidal field on the size of the
galactic disc, giving rise to a size-density correlation analogous to that in
ellipticals: Galactic discs in low-density environments should be larger; an
effect which we elaborate on in the Appendix. Finally, $(iv)$ unlike angular
momentum-based alignment models, the linear alignment model would not give
rise to $B$-modes in the ellipticity field and $(v)$ there should be a nonzero
bispectrum of the ellipticity field at tree-level first order perturbation
theory, analogous to gravitational lensing. All these points show that linear
alignments for spiral galaxies are markedly different compared to quadratic
alignments, as they would result from an angular momentum-based model.

Figure 3: Scaling of the concentration parameter $c(M_{\mathrm{vir}},z)$ with
virial halo mass $M_{\mathrm{vir}}$ and redshift $z$ for two common halo
profile families (top panel) and dependence of the squared velocity ratio
$\sigma^{2}/\upsilon_{\mathrm{circ}}^{2}$ (bottom panel), both derived with
the fitting formula from 2018.

####  2.4.5 Tidal fields in the large-scale structure

With these considerations, the magnitude of intrinsic alignments tidally
generated by the large-scale structure can be estimated: The magnitude of the
perturbing tidal field is $\Omega^{2}$, which needs to be small compared to
the strength of the potential binding the stars to the galaxy: In fact, the
centrifugal balance condition implies

$\omega^{2}=\frac{GM_{\mathrm{vir}}}{r_{s}^{3}}\frac{B}{A}$ (27)

which is effectively the third Kepler-law modified by a numerical factor
$B/A\neq 1$. Then, with
$M_{\mathrm{vir}}=4\pi/3\>r_{\mathrm{vir}}^{3}\Delta\Omega_{m}\rho_{\mathrm{crit}}$
with the virial overdensity $\Delta\simeq 200$ and the relation
$r_{s}=r_{\mathrm{vir}}/c(M_{\mathrm{vir}},z)$ for NFW-haloes one arrives at:

$\omega^{2}=\frac{c(M_{\mathrm{vir}},z)^{3}}{2}\frac{B}{A}\Delta\Omega_{m}H_{0}^{2}$
(28)

such that the orbital time scale $1/\omega$ is shorter by a factor of
$c^{3}\Delta\Omega_{m}$ compared to the Hubble time scale $1/H_{0}$, where
$\Delta$ and $\Omega_{m}$ are universal, but $c(M_{\mathrm{vir}},z)$ remains a
weak function of $M_{\mathrm{vir}}$ and $z$. A typical perturbation of the
density field in the large-scale structure would have an amplitude of
$\delta\simeq 1$, resulting in a perturbation of the potential following the
Poisson-equation $\Delta\Phi=3/2\>\Omega_{m}H^{2}$, which gives immediately a
value of $\Omega^{2}=3/2\>\Omega_{m}H_{0}^{2}$ one can expect for the strength
of the tidal field. Then, the perturbations are necessarily weak,

$\frac{\omega^{2}}{\Omega^{2}}=\frac{c(M_{\mathrm{vir}},z)^{3}}{3}\frac{B}{A}\Delta\gg
1,$ (29)

driven by the fact that $\Delta$ and $c(M_{\mathrm{vir}},z)^{3}$ are much
larger than one. Reiterating this reasoning for elliptical galaxies would
start at virial equilibrium instead of centrifugal equilibirum,
$\omega^{2}=\sigma_{\mathrm{vir}}^{2}/r_{\mathrm{vir}}^{2}=GM_{\mathrm{vir}}/r_{\mathrm{vir}}^{3}$
such that

$\omega^{2}=\frac{1}{2}\Delta\Omega_{m}H_{0}^{2}$ (30)

and the tidal field exerted by the large-scale structure is again weak,
$\omega^{2}/\Omega^{2}=\Delta/3\gg 1$. On both cases, the galaxies are
sufficiently strongly bound such that typical tidal fields are weak
perturbations to their own gravitational fields. We suspect that the inverse
density scale introduced in the common parameterisation of the linear
alignment model for elliptical galaxies (Hirata and Seljak, 2004; Brown et
al., 2002; Blazek et al., 2019, 2011, 2015)

$\epsilon\simeq A_{I}\delta\quad\mathrm{with}\quad
A_{I}=-a_{1}\tilde{C}_{1}\frac{\rho_{\mathrm{crit}}\Omega_{m}}{D_{+}(z)}$ (31)

with $\tilde{C}_{1}$ fixed to $5\times
10^{-14}(\mathrm{Mpc}/h)^{3}/(M_{\odot}/h)$ effectively reproduces with the
background density $\rho_{\mathrm{crit}}\Omega_{m}$ the virial overdensity
$\Delta$, as it determines the ratio between internal and external
gravitational fields.

###  2.5 Tidal fields in 3 dimensions

Up to this point, the perturbing tidal shear field was restricted to the
galactic plane: Generalising this situation leads to a perturbing potential

$\Phi=\frac{\Omega^{2}}{2}\lambda^{(n)}_{ij}x_{i}x_{j}\quad\rightarrow\quad\partial_{i}\partial_{j}\Phi=\Omega^{2}\lambda^{(n)}_{ij}$
(32)

where a suitable basis $\lambda^{(n)}_{ij}$ could be provided by the symmetric
subset of Gell-Mann matrices (as a generalisation of the Pauli-matrices used
in the analogous expansion in eqns. 6 and 7), which form together with the
unit matrix a basis of all tidal shear fields $\partial_{i}\partial_{j}\Phi$.
There are six symmetric Gell-Mann matrices, and the effect of two of them,
specifically $\lambda^{1}$ and $\lambda^{(3)}$, was already investigated
above, as they contain the Pauli-matrices in their top left $2\times 2$
submatrix:

$\lambda^{(1)}=\left(\begin{array}[]{ccc}0&1&0\\\ 1&0&0\\\
0&0&0\end{array}\right),\quad\lambda^{(3)}=\left(\begin{array}[]{ccc}1&0&0\\\
0&-1&0\\\ 0&0&0\end{array}\right).$ (33)

The unit matrix, which is likewise an element of the Gell-Mann matrices, would
cause an isotropic change in size of the galactic disc (c.f. Appendix
Appendix: Isotropic increase in size). The two Gell-Mann matrices
$\lambda^{(4)}$ and $\lambda^{(6)}$,

$\lambda^{(4)}=\left(\begin{array}[]{ccc}0&0&1\\\ 0&0&0\\\
1&0&0\end{array}\right),\quad\lambda^{(6)}=\left(\begin{array}[]{ccc}0&0&0\\\
0&0&1\\\ 0&1&0\end{array}\right),$ (34)

are able to tilt the orbital plane of the disc, as exemplified by Fig. 4. This
effect can either increase or decrease the projected ellipticity, depending on
the unperturbed orientation of the galactic disc. It evolves dynamically, as
the effect of the tidal gravitational field perpendicular to the galactic disc
is cumulative. Finally, the Gell-Mann matrix $\lambda^{(8)}$,

$\lambda^{(8)}=\frac{1}{\sqrt{3}}\left(\begin{array}[]{ccc}1&0&0\\\ 0&1&0\\\
0&0&-2\end{array}\right),$ (35)

combines the effect of changing the size of the disc in an isotropic fashion
while at the same time causing oscillations above and below the unperturbed
orbital plane. Because the orbits are not closed in the perturbed potential
anymore, the orbital plane will undergo a precession motion which effectively
leads to an isotropisation of the system and a more spherical appearance with
low ellipticity as shown in Fig. 4.

Figure 4: Perturbed orbits (green) in comparison to unperturbed circular,
planar orbits (blue) in 3 dimensions, for $\lambda^{(4)}$ (top panel;
$\lambda^{(6)}$ would show the same alignment rotated by 90 degrees) and
$\lambda^{(8)}+\lambda^{(4)}$ (bottom panel, the linear combination was chosen
as the effect of $\lambda^{(8)}$ alone would be in the plane of the disc).

##  3 Ellipticity spectra

The primary observable of weak lensing surveys are (tomographic) ellipticity
spectra, which map out line of-sight-integrated fluctuations of the tidal
gravitational field in a wide redshift range reaching out beyond unity. As a
subdominant contribution to the lensing-induced ellipticity spectra there are
intrinsic alignments of elliptical galaxies, which reflect redshift-
distribution averaged correlations of tidal shear fields in the respective
redshift bin. With the results from the previous chapters, namely that disc
galaxies follow the same linear alignment model as ellipticals albeit with a
smaller alignment parameter, one can relate their intrinsic ellipticity
spectra by an overall scaling. Numerically, we follow the analysis outlined by
Ghosh et al. (2021), to which we refer the reader.

Angular ellipticity spectra for lensing $C_{ij}^{\gamma\gamma}(\ell)$, for
intrinsic alignments $C_{ij}^{\epsilon\epsilon}(\ell)$ and for their cross-
correlation $C_{ij}^{\gamma\epsilon}(\ell)$ are shown in Fig. 5, for a Euclid-
like survey with median redshift of 0.9. For clearer visualisation, we choose
to show a simplified setup with only three bins, with equal fractions of the
total galaxy population. While there are no cross-correlations across
different redshifts bins in intrinsic alignments, cross-correlations of
lensing do exist, but are excluded in the plots in order not to overcrowd the
figure. Parallel to the $GI$ and $II$-spectra for elliptical galaxies we show
the corresponding spectra of disc galaxies as well, both weighted by the
relative fraction of elliptical galaxies (making up roughly a third of the
galaxy population) and disc galaxies (accounting for the remaining two
thirds). The shape noise term assumes a standard value of
$\sigma_{\epsilon}=0.4$ for an individual galaxy, scaled down with
$n=4.727\times 10^{8}$ galaxies per steradian.

Figure 5: Angular ellipticity spectra for a simplified three-bin tomographic
survey with Euclid’s redshift distribution, comparing weak gravitational
lensing with intrinsic alignments of elliptical galaxies, and in particular
with alignments of disc galaxies according to a quadratic and to a linear
alignment model. We assume that spirals make up two thirds of the galaxy
population, and ellipticals the remaining third.

The ratios $C^{\epsilon\epsilon}_{ij}(\ell)/C^{\gamma\gamma}_{ij}(\ell)$ and
$C^{\gamma\epsilon}_{ij}(\ell)/C^{\gamma\gamma}_{ij}(\ell)$ between the
intrinsic alignment spectra and the weak lensing spectrum are shown in Fig. 6,
again for only three bins and for the choice $i=j$. Most importantly, one can
see by how much an angular-momentum based, quadratic alignment model would
overpredict the effect (Schäfer and Merkel, 2017), and how the linear
alignment model for spiral galaxies evades detection in earlier data sets
(Heymans et al., 2013). But most importantly, the dynamical model suggests
that the alignment amplitudes of elliptical and spiral galaxies should be
related, as the alignment parameter reflects effectively the area of the
luminous component of the systems in both cases, suggesting the question
whether marked correlation functions which combine intrinsic ellipticity
correlations with galaxy properties such as Sérsic-indices or Sérsic-radii can
isolate more information about the alignment mechanisms.

Figure 6: Ratios between the angular ellipticity spectra of intrinsic
alignments and weak gravitational lensing analogous to Fig. 5.

Starting with the numbers derived by Ghosh et al. (2021), who derive estimates
on the observability of ellipticity-lensing correlation in aligned ellipticals
in a Euclid-like survey: For ellipticals making up $q=1/3$ of the galaxy
population and an alignment parameter
$D_{IA}=-(10^{-5}~{}\mathrm{Mpc}/h)^{2}$, a Euclid-like redshift distribution
with a median redshift of $0.9$, 5 tomographic bins and integrating up all
multipoles from $\ell=10$ to $\ell=3000$ yields a statistical significance
$\Sigma/\sqrt{f_{\mathrm{sky}}}\simeq 10\ldots 15$ including effects on
nonlinear structure formation (and correspondingly higher tidal shear fields)
on small scales. Assuming that the measurement is limited by shape noise and
by cosmic variance from weak gravitational lensing, one can use the
proportionality of the spectra to scale down the significance for ellipticity
correlations of disc galaxies: Relative to the $GI$-spectrum for ellipticals
which is measurable at $\Sigma/\sqrt{f_{\mathrm{sky}}}\simeq 15$, the
$GI$-spectrum for spirals should generate a signal at a significance lower by
$\sqrt{(1-q)/q/\alpha}\simeq 0.7$. We assume the ratio of alignment parameters
$\alpha$ a value of $4$, intermediate between Hilbert et al. (2017) and
Jagvaral et al. (2022). Relative to the $II$-spectrum of elliptical galaxies
can be measured at a signficiance of about
$\Sigma/\sqrt{f_{\mathrm{sky}}}\simeq 10$, the $II$-spectrum of spirals alone
should generate a significance lower by $(1-q)/q/\alpha\simeq 0.5$, and the
$II$-type cross-correlation between spirals and ellipticals would be at
significances lower by $\sqrt{(1-q)/\alpha}\simeq 0.4$. This effectively puts
all spectra within reach of a Euclid-type survey.

##  4 Summary and discussion

The subject of this paper is the construction of a linear alignment model for
disc galaxies that conforms with the results of simulations of galaxy
formation and with observations of intrinsic alignments in weak lensing data.
In order to agree with observations, alignments of disc galaxies must be
weaker than in elliptical galaxies, and simulations do in fact suggest a weak
alignment at linear instead of quadratic order, with the interesting
consequence that $GI$-type alignments should exist in disc galaxies, very much
in contrast to the predictions by quadratic, angular-momentum based alignment
models.

  1. 1.

The first step was a perturbative calculation, with an associated energy
consideration, to derive the perturbative effect of a quadrupolar tidal
gravitational field on a circular orbit. We showed that the ellipticity
reflects the ratio between the squared orbital frequency of the unperturbed
orbit and the magnitude of the tidal gravitational field. This model produces
naturally an ellipticity which is proportional to the magnitude of the tidal
field. Being more specific, we then set up a direct simulation of the orbit of
a star inside the gravitational potential of an NFW-halo, tidally perturbed
again by a quadrupolar field. This direct simulation confirmed the
perturbative results. In addition, we could show that the induced ellipticity
depends on orbital radius, as stars at higher radii are less bound and are
more susceptible to tidal perturbation.

  2. 2.

With the alignment of both elliptical and spiral galaxies following a linear
alignment model and both being determined by the external tidal gravitational
field, the question of how the relative alignments are determined by galaxy
properties arises naturally. Elliptical galaxies react in a homologous way to
tidal fields, whereas the tidally-induced change in shape of spirals depends
on radius. The linear alignment model would now stipulate that the alignment
of spirals is weaker as the discs are smaller than the virial radius of the
host halo. At the NFW-scale radius, the ellipticity of a galactic disc is
$c(M_{\mathrm{vir}},z)^{2}$ times smaller than that of an elliptical, with
$c(M_{\mathrm{vir}},z)$ being the NFW-concentration parameter. Averaging the
ellipticity over all radii, weighted with an exponential brightness
distribution alleviates this number and sets the alignments strength of
spirals to be about an order of magnitude smaller than that of ellipticals.

  3. 3.

The physical picture behind linear alignments of disc galaxies suggests a
number of observational relations that could be tested with data: The
alignment should be stronger with higher Sérsic-index, and higher for systems
with larger discs. There should be a scaling of the alignment with the mass of
the system, driven by the change in NFW-concentration parameter
$c(M_{\mathrm{vir}},z)$, too: As low-mass systems tend to be more centrally
concentrated, their alignment should be weaker. Whether empirically introduced
laws to allow for redshift evolution can be replaced by the known evolution of
the NFW-concentration parameter $c(M_{\mathrm{vir}},z)$: Here, the alignment
would be highest for high $z$ and decrease towards the current epoch. To what
extent evolution of alignments with redshift and their scaling with mass are
driven by the dependence of the NFW-concentration parameter and what
additional astrophysical effects could play a role would be an interesting
question for surveys and for simulations.

  4. 4.

A further interesting consequence of the linear alignment model are changes in
shape: Both ellipticals and spirals should tend to be larger in low-density
regions. Flexion-like distortions from higher-order tidal shear fields on
galactic discs should exist but would be very small, given the result that the
effect induced by $\partial^{3}\Phi$ on ellipticals as investigated in Giesel
et al. (2022a) is already small. But physically, these distortions would be
able to warp a galactic disc, which would constitute a link between alignments
and disc dynamics.

  5. 5.

Ellipticity spectra resulting from the linear alignment model naturally
resemble those of elliptical galaxies as they are related just by a constant
prefactor. With their smaller amplitudes ellipticity correlations of spiral
galaxies evade constraints on alignments derived from current weak lensing
data and do not produce signals that can be measured at reasonable
significances. More important, though, are the much smaller amplitudes of the
linear alignment model in spiral galaxies compared to angular-momentum based
models: Typically, the spectra are smaller by more than an order of magnitude
at high multipoles. Additionally, the phenomenology of the spectra is quite
different: In contrast to the quadratic, angular-momentum based alignment
model which predicts no cross-correlation between the shapes of elliptical and
spiral galaxies, the linear alignment model would predict cross-correlations
between the shapes of the two galaxy types in a very natural way. Moreover,
the linear alignment model would produce pure $E$-mode correlations in the
ellipticity field, which would apply likewise for spiral and elliptical
galaxies. Again, the quadratic model would be markedly different, with a high
amplitude $B$-mode correlation being predicted. Lastly, one would expect a
nonzero cross-correlation between intrinsic shapes of both ellipticals and
spirals with gravitational lensing shear, with the typical anti-correlations
present in both galaxy types. Linking tidal alignment to galaxy properties one
should expect a dependence on Sérsic-index as more centrally concentrated
objects would show weaker alignments. Estimates suggest that all intrinsic
ellipticity spectra are within reach of Euclid.

Discs galaxies regularly have galactic bulges at their centres, which show
analogous alignments in simulations by Jagvaral et al. (2022), up to the fact
that they tend to have rather high Sérsic indices close to $n=4$, are rather
strongly bound with high velocity dispersions $\sigma$ and are rather compact
with a small Sérsic radii $R$; numbers of
$\sigma=140~{}\mathrm{km}/\mathrm{s}$ and $R=1/2~{}\mathrm{kpc}$ have been
reported for Milky Way (Valenti et al., 2018). Their amplitude seems to be
similar to that of elliptical galaxies, and is much stronger than that of the
host disc galaxy. Whether this is a consequence of the homologous distortion
of virial systems remains to be investigated. We intend to develop the topic
of tidally distorted galactic discs and the resulting ellipticity correlations
further by considering the full 3-dimensional tidal field, which would
effectively bridge between the concept of disc warps and intrinsic alignments,
or perhaps even link alignments with tidally-induced bar formation in spiral
galaxies.

## Acknowledgements

We would like to thank Robert Reischke and Joerg Jaeckel for their thoughtful
remarks.

##### Funding information

BG would also like to acknowledge financial support from the DST-INSPIRE
Faculty fellowship (Grant no. DST/INSPIRE/04/2020/001534) and the DAAD
Research Stays for University Academics and Scientists Fellowship (Ref No.
91879665), and Heidelberg University for hosting during the research visit.
ESG acknowledges support from the German National Academic Foundation.

## Appendix: Isotropic increase in size

Fig. 7 shows the effect of an isotropic tidal field $\Delta\Phi$ to the orbit
of a star at $r=r_{s}$. As expected, there is a change in size of the galactic
disc, analogous to the effect on elliptical galaxies (Ghosh et al., 2021) or,
equivalently, to the effect of lensing convergence (Bartelmann and Schneider,
2001). In both cases, galaxies would appear larger in low-density environments
and smaller in high-density environments. Visually, the effect of an isotropic
tidal field on the size of the disc seems to be smaller compared to the shear
depicted in Fig. 1, but one perceives only half of the effect as both semi-
axes are changed by the same factor $1+\Omega^{2}/\omega^{2}/2$ and unlike in
the case of ellipticity, $1+\Omega^{2}/\omega^{2}/2$ relative to
$1-\Omega^{2}/\omega^{2}/2$. The same argument as for disc-averaged
ellipticities would apply here as well, as less tightly bound stars at larger
radii would experience a stronger tidal enlargement of their orbital radii. In
consequence, one should expect an size alignment effect with an equal
amplitude as the ellipticity alignment effect, in analogy to the result for
size alignments in elliptical galaxies (Ghosh et al., 2021). As the size-
alignments for elliptical galaxies are difficult to measure even for surveys
such as Euclid’s, it is doubtful if the analogous measurement for disc
galaxies could yield a statistically significant detection, extrapolating the
results by Ghosh et al. (2021). Repeating the estimates in Sect. 3 in
conjunction with the results of Ghosh et al. (2021) yields the result that
correlations between size perturbations and weak lensing convergence are
likewise difficult to measure with data from a Euclid-like survey.

Figure 7: Perturbed orbits (green) in comparison to unperturbed circular
orbits (blue) at $r=r_{s}$ in the potential of an NFW-halo with an isotropic
tidal field $\Omega^{2}(x^{2}+y^{2})$ giving rise to size change of the disc
analogous to the lensing convergence $\kappa$. The magnitude of the perturbing
field is again $\Omega^{2}=0.04$ in dimensionless units.

## References

  * Schaefer (2009) B. M. Schaefer, Int. J. Mod. Phys. D 18, 173 (2009), 0808.0203.
  * Joachimi et al. (2015) B. Joachimi et al., Space Sci. Rev. 193, 1 (2015), 1504.05456.
  * Kiessling et al. (2015) A. Kiessling et al., Space Sci. Rev. 193, 67 (2015), [Erratum: Space Sci.Rev. 193, 137 (2015)], 1504.05546.
  * Kirk et al. (2015) D. Kirk et al., Space Sci. Rev. 193, 139 (2015), 1504.05465.
  * Troxel and Ishak (2014) M. A. Troxel and M. Ishak, Phys. Rept. 558, 1 (2014), 1407.6990.
  * Lamman et al. (2023) C. Lamman, E. Tsaprazi, J. Shi, N. N. Šarčević, S. Pyne, E. Legnani, and T. Ferreira (2023), 2309.08605.
  * Hirata and Seljak (2004) C. M. Hirata and U. Seljak, Phys. Rev. D 70, 063526 (2004), [Erratum: Phys.Rev.D 82, 049901 (2010)], astro-ph/0406275.
  * Crittenden et al. (2001) R. G. Crittenden, P. Natarajan, U.-L. Pen, and T. Theuns, Astrophys. J. 559, 552 (2001), astro-ph/0009052.
  * Hirata et al. (2007) C. M. Hirata, R. Mandelbaum, M. Ishak, U. Seljak, R. Nichol, K. A. Pimbblet, N. P. Ross, and D. Wake, Mon. Not. Roy. Astron. Soc. 381, 1197 (2007), astro-ph/0701671.
  * Bridle and King (2007) S. Bridle and L. King, New J. Phys. 9, 444 (2007), 0705.0166.
  * Schaefer and Merkel (2012) B. M. Schaefer and P. Merkel, Mon. Not. Roy. Astron. Soc. 421, 2751 (2012), 1101.4584.
  * Chisari et al. (2015) N. E. Chisari, S. Codis, C. Laigle, Y. Dubois, C. Pichon, J. Devriendt, A. Slyz, L. Miller, R. Gavazzi, and K. Benabed, Mon. Not. Roy. Astron. Soc. 454, 2736 (2015), 1507.07843.
  * Rodriguez-Gomez et al. (2017) V. Rodriguez-Gomez et al., Mon. Not. Roy. Astron. Soc. 467, 3083 (2017), 1609.09498.
  * Amendola et al. (2018) L. Amendola et al., Living Reviews in Relativity 21, 345 (2018), 1606.00180.
  * Schäfer and Merkel (2017) B. M. Schäfer and P. M. Merkel, Mon. Not. Roy. Astron. Soc. 470, 3453 (2017), 1506.07366.
  * Tugendhat and Schäfer (2018) T. M. Tugendhat and B. M. Schäfer, Mon. Not. Roy. Astron. Soc. 476, 3460 (2018), 1709.02630.
  * Reischke and Schäfer (2019) R. Reischke and B. M. Schäfer (2019), 1910.05994.
  * Zjupa et al. (2022) J. Zjupa, B. M. Schäfer, and O. Hahn, Mon. Not. Roy. Astron. Soc. 514, 2049 (2022), 2010.07951.
  * Delgado et al. (2023) A. M. Delgado et al., Mon. Not. Roy. Astron. Soc. 523, 5899 (2023), 2304.12346.
  * Nelson et al. (2018) D. Nelson et al. (2018), 1812.05609.
  * Hildebrandt et al. (2017) H. Hildebrandt, M. Viola, C. Heymans, S. Joudaki, K. Kuijken, C. Blake, T. Erben, B. Joachimi, D. Klaes, L. Miller, et al., MNRAS 465, 1454 (2017), 1606.05338.
  * Johnston et al. (2019) H. Johnston et al., Astron. Astrophys. 624, A30 (2019), 1811.09598.
  * Samuroff et al. (2023) S. Samuroff et al. (DES), Mon. Not. Roy. Astron. Soc. 524, 2195 (2023), 2212.11319.
  * Ghosh et al. (2021) B. Ghosh, R. Durrer, and B. M. Schaefer, Mon. Not. Roy. Astron. Soc. 505, 2594 (2021), 2005.04604.
  * Hilbert et al. (2017) S. Hilbert, D. Xu, P. Schneider, V. Springel, M. Vogelsberger, and L. Hernquist, Mon. Not. Roy. Astron. Soc. 468, 790 (2017), 1606.03216.
  * Jagvaral et al. (2022) Y. Jagvaral, S. Singh, and R. Mandelbaum, Mon. Not. Roy. Astron. Soc. 514, 1021 (2022), 2202.08849.
  * Heymans et al. (2013) C. Heymans et al., Mon. Not. Roy. Astron. Soc. 432, 2433 (2013), 1303.1808.
  * Kilbinger et al. (2009) M. Kilbinger et al., Astron. Astrophys. 497, 677 (2009), 0810.5129.
  * Blazek et al. (2019) J. Blazek, N. MacCrann, M. A. Troxel, and X. Fang, Phys. Rev. D 100, 103506 (2019), 1708.09247.
  * Blazek et al. (2011) J. Blazek, M. McQuinn, and U. Seljak, JCAP 05, 010 (2011), 1101.4017.
  * Blazek et al. (2015) J. Blazek, Z. Vlah, and U. Seljak, JCAP 08, 015 (2015), 1504.02510.
  * Giesel et al. (2022a) E. S. Giesel, B. Ghosh, and B. M. Schaefer, Mon. Not. Roy. Astron. Soc. 510, 2773 (2022a), 2107.09000.
  * Child et al. (2018) H. L. Child, S. Habib, K. Heitmann, N. Frontiere, H. Finkel, A. Pope, and V. Morozov, Astrophys. J. 859, 55 (2018), 1804.10199.
  * Giesel et al. (2022b) E. S. Giesel, B. Ghosh, and B. M. Schäfer, Mon. Not. Roy. Astron. Soc. 518, 5490 (2022b), 2208.07197.
  * Dehnen et al. (2006) W. Dehnen, D. McLaughlin, and J. Sachania, Mon. Not. Roy. Astron. Soc. 369, 1688 (2006), %****␣ia_spiral_ojastro.bbl␣Line␣300␣****astro-ph/0603825.
  * Lin and Li (2019) H.-N. Lin and X. Li, Mon. Not. Roy. Astron. Soc. 487, 5679 (2019), 1906.08419.
  * Goodwin et al. (1997) S. P. Goodwin, J. Gribbin, and M. A. Hendry (1997), astro-ph/9704216.
  * Goodwin et al. (1998) S. P. Goodwin, J. Gribbin, and M. A. Hendry, The Observatory 118, 201 (1998).
  * Brown et al. (2002) M. L. Brown, A. N. Taylor, N. C. Hambly, and S. Dye, Mon. Not. Roy. Astron. Soc. 333, 501 (2002), astro-ph/0009499.
  * Valenti et al. (2018) E. Valenti, M. Zoccali, A. Mucciarelli, O. A. Gonzalez, F. Surot, D. Minniti, M. Rejkuba, L. Pasquini, G. Fiorentino, G. Bono, et al., Astronomy & Astrophysics 616, A83 (2018), ISSN 1432-0746.
  * Bartelmann and Schneider (2001) M. Bartelmann and P. Schneider, Physics Reports 340, 291 (2001).


# Revisit of open clusters UPK 39, UPK 41 and PHOC 39 : a new binary open
cluster found

Xianhao Ye CAS Key Laboratory of Optical Astronomy, National Astronomical
Observatories, Chinese Academy of Sciences, Beijing 100101, China  School of
Astronomy and Space Science, University of Chinese Academy of Sciences,
Beijing 100049, China  Jingkun Zhao CAS Key Laboratory of Optical Astronomy,
National Astronomical Observatories, Chinese Academy of Sciences, Beijing
100101, China  Terry D. Oswalt Embry-Riddle Aeronautical University, 1
Aerospace Boulevard, Daytona Beach, FL 32114, USA<EMAIL_ADDRESS>Yong
Yang CAS Key Laboratory of Optical Astronomy, National Astronomical
Observatories, Chinese Academy of Sciences, Beijing 100101, China  School of
Astronomy and Space Science, University of Chinese Academy of Sciences,
Beijing 100049, China  Gang Zhao CAS Key Laboratory of Optical Astronomy,
National Astronomical Observatories, Chinese Academy of Sciences, Beijing
100101, China  School of Astronomy and Space Science, University of Chinese
Academy of Sciences, Beijing 100049, China

(Revised 25 July, 2022; Accepted 27 July, 2022)

###### Abstract

We investigate the three open clusters near Aquila Rift cloud, named as UPK 39
(c1 hereafter), UPK 41 (c2 hereafter) in Sim et al. (2019) and PHOC 39 (c3
hereafter) in Hunt & Reffert (2021), respectively. Using photometric
passpands, reddening, and extinction from Gaia DR3, we construct the color-
absolute-magnitude diagram (CAMD). Using isochrone fits their ages are
estimated as $6.3\pm 0.9$, $8.1\pm 1.4$ and $21.8\pm 2.2$ Myr, respectively.
Their proper motions and radial velocities, estimated using data from Gaia and
LAMOST are very similar. From their orbits, relative distances among them at
different times, kinematics, ages, and metallicities, we conclude that c1 and
c2 are primordial binary open cluster, which are likely to have been formed at
the same time, and c3 may capture c1, c2 in the future.

Open star cluster (1160); Stellar kinematics (1608)

††journal: AJ††software: Astropy (Astropy Collaboration et al., 2013), DBSCAN
(Ester et al., 1996), Galpy (Bovy, 2015), HDBSCAN (Campello et al., 2013;
McInnes et al., 2017), Matplotlib (Hunter, 2007), Numpy (van der Walt et al.,
2011), Pandas (McKinney, 2010), Scipy (Virtanen et al., 2020), Scikit-learn
(Pedregosa et al., 2012), Topcat (Taylor, 2005)

##  1 Introduction

Open clusters (OCs), gravitationally bound stars originally formed from giant
molecular clouds (GMCs; Lada & Lada 2003), are building blocks of the Milky
Way. Catalogs of OCs have been compiled for over a century (Dreyer, 1888).
High quality astrometric and photometric data from Gaia (Gaia Collaboration et
al., 2016, 2018, 2021), combined with new highly efficient tools such as
DBSCAN (Ester et al., 1996), HDBSCAN (Campello et al., 2013; McInnes et al.,
2017), UPMASK (Unsupervised Photometric Membership Assignment in Stellar
Clusters) (Krone-Martins & Moitinho, 2014), have dramatically increased the
number of OCs (Cantat-Gaudin et al., 2018, 2019; Castro-Ginard et al., 2018,
2019, 2020, 2022; Liu & Pang, 2019; Sim et al., 2019; Bica et al., 2019;
Ferreira et al., 2019, 2021; He et al., 2021, 2022; Hunt & Reffert, 2021; Qin
et al., 2021). A complete census of OCs within the solar neighborhood will
provide a sound basis to investigate a number of scientific questions. OCs
could support us investigate other scientific program as well, such as the
metallicity gradient of the Galaxy and radial migration (Netopil et al., 2022;
Zhang et al., 2021; Chen & Zhao, 2020), spiral arms (Castro-Ginard et al.,
2021b), and moving groups (Zhao et al., 2009; Liang et al., 2017; Yang et al.,
2021). With more and more identified OCs, delicate analysis for them is needed
to provide further understanding of the Galaxy.

Figure 1: Pluses (+) represent stars in the three OCs c1(blue), c2(orange) and
c3(green) of our study, located in 5D phase space. The corresponding color is
adopted through this paper. Other known OCs from literature in this part of
the sky are represented by colored triangles.

In the region around the Aquila Rift cloud, we identified three OCs with
similar proper motions (PMs). Two of them, c1, c2, were first found by Sim et
al. (2019) and the other, c3 Hunt & Reffert (2021). Sim et al. (2019)
identified the centers of c1 and c2, respectively, as
$\left(\mu_{\alpha}^{*},\mu_{\delta}\right)=\left(3.29\pm 0.39,-8.65\pm
0.41\right)$ $\mathrm{mas}\cdot\mathrm{y}^{-1}$ and $\left(2.51\pm
0.22,-8.13\pm 0.23\right)$ $\mathrm{mas}\cdot\mathrm{y}^{-1}$. Hunt & Reffert
(2021) located the center of c3 as,
$\left(\mu_{\alpha}^{*},\mu_{\delta}\right)=\left(1.89\pm 0.05,-8.75\pm
0.05\right)$ $\mathrm{mas}\cdot\mathrm{y}^{-1}$. Sim et al. (2019) estimated
the ages of c1 and c2 to be about $2.8$ Myr and $7.1$ Myr, respectively. As
these three OCs are close to each other, young, and have similar PMs, we are
interested in their detailed dynamic properties, such as the possibility they
are gravitationally bound or interacting.

A certain percentage ($\sim 12\%$) of OCs comprise binary or multiple systems
in solar neighborhood (de La Fuente Marcos & de La Fuente Marcos, 2009),
either primordial or captured. At the present rapid pace of OC discovery in
general, one can expect new binary and multiple OCs to be identified.

In this paper, we estimate the ages of c1, c2 and c3 from isochrone fitting.
We also estimate their three-dimensional space velocities. From these data and
metallicity estimates we examine the possibility that they are or were
physically associated.

The structure of this paper is as follows : Sec. 2 describes the data
extracted from Gaia EDR3; In Sec. 3 we detail the procedures of member
selections of the three clusters and isochrone fitting using Gaia DR3
passpands; Their ages, kinematic properties, relative distances, and
metallicities are presented in Sec. 4; Our conclusions are given in Sec. 5;
This paper is summarized in Sec. 6.

##  2 Data

The Gaia EDR3 astrometric parameters (Gaia Collaboration et al., 2021;
Lindegren et al., 2021b)
$\left(\alpha,\delta,\varpi,\mu_{\alpha}^{*},\mu_{\delta}\right)$ are used to
identify member candidates of the three clusters. We first select sources
within a $40^{\circ}\times 40^{\circ}$ zone on sky, and restrict the distance
via the parallaxes $\varpi$. The relative errors of $\varpi$ and PMs are
restricted to $10\%$ to ensure the qualities of data. Furthermore, we apply
the renormalized unit weight error (RUWE; Lindegren et al. 2018, 2021b) to
refine our selections. The following constraints were used to define our
primary sample:

  * •

260 $<$ $\alpha$ $<$ 300, -20 $<$ $\delta$ $<$ 20, 1/0.700 $<$ $\varpi$ $<$
1/0.250;

  * •

parallax over error $>$ 10;

  * •

$\sigma_{\mu_{\alpha}^{*}}/|\mu_{\alpha}^{*}|<0.1$;

  * •

$\sigma_{\mu_{\delta}}/|\mu_{\delta}|<0.1$;

  * •

RUWE $<$ 1.4.

Using the above criteria, 965,430 sources which contain all three Gaia
passbands $\left(G,G_{\mathrm{BP}},G_{\mathrm{RP}}\right)$ are retained as our
primary sample. However, only a small fraction of this sample have radial
velocities RVs in Gaia DR3 (Katz et al., 2022) and the number will drop
dramatically if we restrict the relative error of RV too strictly. It is more
reasonable to search for clusters members in 5D phase space
$\left(x,y,z,\kappa\cdot\mu_{\alpha}^{*}\cdot d,\kappa\cdot\mu_{\delta}\cdot
d\right)$ (3D Cartesian coordinates and 2D tangential velocity), where
$\kappa=4.74047$ and $d$ is distance. The photogeometric distance from Bailer-
Jones et al. (2021) is applied in the calculation and the inverse of parallax
$1/\varpi$ is used for the distance for those without the photometric distance
from Bailer-Jones et al. (2021). In addition, we correct $\varpi$ with the
code111https://gitlab.com/icc-ub/public/gaiadr3_zeropoint from Lindegren et
al. (2021a) to calculate the parallax zero point. Galpy (Bovy, 2015) is used
to calculate the Galactic Cartesian coordinates and velocities. We set the
radial distance and height of Sun in Galactocentric frame at
$R_{\odot}=8.3\quad\mathrm{kpc}$ (Gillessen et al., 2009),
$Z_{\odot}=0.027\quad\mathrm{kpc}$ (Chen et al., 2001), and its velocity as
$\left(U,V,W\right)_{\odot}=\left(11.1,12.24,7.25\right)\quad\mathrm{km\cdot
s^{-1}}$ (Schönrich et al., 2010) relative to the Local Standard of Rest
(LSR), where $V_{\mathrm{LSR}}=240.0\quad\mathrm{km\cdot s^{-1}}$ according to
Reid et al. (2014).

##  3 Methods

###  3.1 Member Selections in 5D phase space

There are a number of known OCs in the region covering our primary sample. We
collected catalogs from Cantat-Gaudin et al. (2018, 2019); Cantat-Gaudin &
Anders (2020a); Cantat-Gaudin et al. (2020b); Castro-Ginard et al. (2018,
2019, 2020); Kharchenko et al. (2013); Liu & Pang (2019); Sim et al. (2019);
He et al. (2021); Hunt & Reffert (2021) and show them in Fig. 1 as triangles
of different colors. Note we also use parallaxes gleaned from the literature
to restrict our search to OCs in the same distance range as our sample. We use
the same definition as in our previous papers (Ye et al., 2021a, b) to
calculate the number of neighbors of each star in 5D phase space. Stars within
a radius of $15$ pc in $\left(x,y,z\right)$ and a radius of $1.1$
$\mathrm{km\cdot s^{-1}}$ in tangential velocity are defined as the neighbors
of a given source. Applying a lower limit of $\mu+3\sigma$ (mean value plus
three times standard deviation of neighbors) retains $3,432$ stars, which
includes almost all the clusters in the literature. We then test both
HDBSCAN(Campello et al., 2013; McInnes et al., 2017) and DBSCAN(Ester et al.,
1996; Pedregosa et al., 2012) to search for member candidates of clusters
among those $3,432$ stars. Eventually, we focus on the clusters in the region
of $275^{\circ}<\alpha<285^{\circ}$, $-5^{\circ}<\delta<4^{\circ}$. Using
$\epsilon=0.5,\texttt{minPts}=30$ via DBSCAN in $\left(\alpha,\delta\right)$,
this procedure yields the three OCs marked as colored crosses in Fig. 1. The
blue and orange crosses correspond to c1 and c2. The green cluster corresponds
to c3.

###  3.2 Member refinement with Radial Velocity

Figure 2: Dec vs. RV for the clusters candidates with Gaia DR3 and LAMOST DR8
RVs. For Gaia RVs, only $\sigma_{\mathrm{RV}}$ less than the mean
$\sigma_{\mathrm{RV}}$ in each cluster were used. Data from LAMOST are marked
as triangles. The horizontal lines indicate the final estimated RVs for the
clusters which are the mean values.  Figure 3: [M/H]] vs. age. The gray
shadows show the coverages of PARSEC isochrones used in this paper. The
colored crosses present the top $1\%$ of isochrones that have smallest
$\bar{d^{2}}$ for the three clusters.

Both the literature PMs values and those determined for our OCs member
candidates are quite similar. However, we had no prior knowledge of their
radial velocities. In this paper, we found RVs for a few stars in Gaia DR3.
The mean errors in RVs of the member candidates are : $10.81$, $7.26$, and
$8.32$ $\mathrm{km\cdot s^{-1}}$ for c1, c2, and c3, respectively. We then set
these values as the upper limits for $\sigma_{\mathrm{RV}}$ in our subsequent
RV analysis. We found two stars from c1 and c2, respectively, have RVs in
LAMOST DR8 low resolution archive (Liu et al., 2015; Cui et al., 2012; Zhao et
al., 2012, 2006). In each cluster, some RVs deviate greatly from the median
value, especially at $G>14.5$ mag. Such stars are removed from our candidate
list. Finally, we adopt the RV for each cluster as the mean value of the
remaining members. Note that c1 contains one star for which only a LAMOST RV
is available; this star was included in the calculation of the mean RV for c1.
In Fig. 2, the remaining RVs from Gaia DR3 are shown in Dec vs. RV as dots in
different colors, and two RVs from LAMOST (in c1 and c2) are marked as
triangles. The adopted values are $-8.05$, $-9.53$, $-4.15$ $\mathrm{km\cdot
s^{-1}}$ for c1, c2 and c3, respectively, and are indicated with the dashed
lines in different colors in Fig. 2.

###  3.3 Color-absolute-magnitude Diagram and Isochrone Fitting

Figure 4: CAMD for c1, c2 and c3 (left to right), along with adopted ages and
the corresponding isochrones with [M/H] between $-0.5$ to $0.3$ dex (interval
$0.05$ dex) are shown.

The photometric parameters $\left(G,G_{\mathrm{BP}},G_{\mathrm{RP}}\right)$
and their errors, reddening, extinction used in this paper are from Gaia DR3,
which are obtained by cross-matching the member candidates with Gaia DR3 in
TOPCAT (Taylor, 2005). The errors in
$\left(G,G_{\mathrm{BP}},G_{\mathrm{RP}}\right)$ are from CDS
222https://vizier.cds.unistra.fr/viz-bin/VizieR-3?-source=I/355/gaiadr3. For
reddening E(BP-RP) and extinction AG, the uncertainties are taken as half the
difference between upper and lower confidence levels, same as the
uncertainties adopted in distance in this paper. There are $>73\%$, $>86\%$
and $>96\%$ stars having reddening and extinction data in Gaia DR3. With the
dereddened color $\left(G_{\mathrm{BP}}-G_{\mathrm{RP}}\right)_{0}$ and
absolute magnitude
$M_{\mathrm{G}}=G-5\lg\left(\frac{d[pc]}{10}\right)-A_{\mathrm{G}}$, where $d$
is the distance of each star, we estimated the ages of these clusters by
fitting isochrones to each cluster’s CAMD. Two stars in c2 that diverge from
the clearly coeval sequence were removed from the member candidate sample. The
procedure outlined above yields $61$, $103$, and $78$ member stars in each
cluster. In Sim et al. (2019), c1 and c2 were determined to have ages of about
$2.8$ Myr and $7.1$ Myr, respectively, according to isochrone fitting. The
main sequence of c3 is noticeably different from the other two clusters. We
use a series of PARSEC (version 1.2s) isochrones
333http://stev.oapd.inaf.it/cgi-bin/cmd (Bressan et al., 2012; Chen et al.,
2014, 2019) with Gaia EDR3 photometric data to fit with these three clusters.
The range of [M/H] of the isochrone grid is $-0.5\sim 0.3$ dex with a interval
of $0.05$ dex. The ages of isochrones range from $1$ Myr to $30$ Myr with a
interval of $0.1$ Myr. The closest isochrone for a cluster is adopted by
minimizing $\bar{d^{2}}$ in CAMD. $\bar{d^{2}}$ is defined as Eq. 1 in Liu &
Pang (2019). For a cluster with $n$ member candidates, $\bar{d}^{2}$ for a
given isochrone is :

$\displaystyle\bar{d^{2}}=\sum\limits_{i=1}^{n}\left|\mathbf{x}_{i}-\mathbf{x}_{i,nn}\right|^{2}/n\,,$
(1)

where $\mathbf{x}_{i}$ gives the position of member star $i$ in CAMD, and
$\mathbf{x}_{i,nn}$ represents the nearest neighbor of member $i$ in the
isochrone in the same parameter space. The nearest neighbor is identified by
the $k$-D tree approach in scipy (Virtanen et al., 2020). The influences of
the errors in $G_{\mathrm{BP}}$, $G_{\mathrm{RP}}$ and $G$ passpands have also
been considered. For each star, we use the errors in
$\left(G_{\mathrm{BP}}-G_{\mathrm{RP}}\right)_{0}$ and $M_{\mathrm{G}}$ to
randomly generate $100$ different data points following a Gaussian
distribution and calculate $\left|\mathbf{x}_{j}-\mathbf{x}_{j,nn}\right|$ for
$j$ in those data points. The mean value is adopted for star $i$ as
$\left|\mathbf{x}_{i}-\mathbf{x}_{i,nn}\right|$. Then $\bar{d^{2}}$ is
evaluated between each cluster and each isochrone. To obtain the ages, we use
the first $1\%$ of the isochrones which have the smallest $\bar{d^{2}}$. In
Fig. 3, we present the $1\%$ isochrone age vs. [M/H]. It is apparent in this
figure that c1 and c2 are similar both in age (c1 : $6.3\pm 0.9$ Myr, c2 :
$8.1\pm 1.4$ Myr) and [M/H] (c1 : $-0.3\pm 0.1$ dex, c2 : $-0.3\pm 0.1$ dex),
and c3 is much older and richer in [M/H] than the others.

##  4 Results

###  4.1 Age

The estimated ages of c1, c2, c3 from our best-fitting isochrones are $6.3\pm
0.9$, $8.1\pm 1.4$ and $21.8\pm 2.2$ Myr, respectively. In Fig. 4, cluster
members and isochrones with the corresponding best ages and [M/H] between
$-0.5$ to $0.3$ dex are presented.

###  4.2 Galactocentric Cylindrical velocity

Figure 5: Left panel : Orbits of c1, c2 and c3 traced back to their
birthplaces. c3 was close to c1 and c2 when they were born and c2 was near c1
in the time of the formation of c1. Right panel : Circles are where these
clusters were born and crosses are their locations after $400$ Myr. They did
and will orbit together in the thin disk. The arrow in each panel shows the
direction of motion.

Averaging over our clusters members, the PMs of c1, c2 and c3 are
$\left(\mu_{\alpha}^{*},\mu_{\delta}\right)=\left(3.24\pm 0.31,-8.59\pm
0.43\right)$ $\mathrm{mas}\cdot\mathrm{y}^{-1}$, $\left(2.65\pm 0.27,-8.35\pm
0.27\right)$ $\mathrm{mas}\cdot\mathrm{y}^{-1}$ and $\left(2.06\pm
0.20,-9.00\pm 0.26\right)$ $\mathrm{mas}\cdot\mathrm{y}^{-1}$, respectively.
The PMs and RVs imply similar space velocities for these three clusters in a
Galactocentric Cylindrical coordinates. We calculate that
$\left(V_{\mathrm{R}},V_{\phi},V_{\mathrm{Z}}\right)=\left(-5.18,238.33,-7.17\right)$
$\mathrm{km}\cdot\mathrm{s}^{-1}$, $\left(-3.54,236.09,-6.91\right)$
$\mathrm{km}\cdot\mathrm{s}^{-1}$ and $\left(-8.42,239.61,-3.61\right)$
$\mathrm{km}\cdot\mathrm{s}^{-1}$ for c1, c2 and c3, respectively, where
$V_{\mathrm{R}}>0$ represents moving away from Galactic center, $V_{\phi}$ is
the direction of Galactic rotation, and $V_{\mathrm{Z}}$ points to the
Galactic North Pole.

###  4.3 Orbit and Separation

Figure 6: Relative distances among c1, c2, c3 vs. time. The gray shaded areas
present the uncertainties of relative distances and the green shaded areas
denote their relative distances in multiples of $\bar{r_{\mathrm{t}}}$.

As the velocities of these three OCs are very similar, we are curious about
their Galactic orbits and three origins. We use MWPotential2014 in Galpy to
trace their birthplaces that correspond to their ages. In the left panel of
Fig. 5, we trace back c1, c2, c3 to their birthplaces, marked as filled
circles. The orange triangle shows c2’s position when c1 was born. The
positions of c3 in $6.3$ Myr and $8.1$ Myr ago are presented as other green
signs. The right panel of Fig. 5 presents the orbits from their birthplaces to
where they will be $400$ Myr henceforth. The circles represent the origins and
the crosses are their locations in the future. All the three clusters will
orbit in unison in the Galactic thin disk. The minimal distances among them
are on the order of dozens of parsecs. We use the standard deviations of
positions and velocities of member stars in each cluster to randomly generate
parameters following a Gaussian distribution. The adopted deviations in RVs
are taken as $1.1$ $\mathrm{km}\cdot\mathrm{s}^{-1}$, the same as we used in
Sec. 3.1 to extract member candidates from tangential velocities. We produce
three data points in each position parameters and five in each velocity
parameter. We then construct $3,375$ different orbits for each cluster and
calculate the relative distance among c1, c2, and c3 at different times.

In Fig. 6, the separations between two clusters are presented. The black
(solid, dotted, dashed) lines in the panels are the separations calculated by
the adopted positions and velocities for c1, c2, and c3. The gray shaded
regions indicate the uncertainties in the orbits. These lower and upper bounds
are $\mu-\sigma$ and $\mu+\sigma$ of the separations calculated with the
generated data points. If the relative distance between two clusters is less
than three times the outer radius, the mutual interaction becomes significant
(Innanen et al., 1972). de La Fuente Marcos & de La Fuente Marcos (2009) used
three times the mean tidal radius ($r_{\mathrm{t}}$) as the upper limit to
select paired OCs. The adopted mean value of tidal radius is $10$ pc in de La
Fuente Marcos & de La Fuente Marcos (2009). However, this value may be smaller
than the actual average value. In Kharchenko et al. (2013), the mean value of
$r_{\mathrm{t}}$ is $11.59$ pc as determined from an analysis of nearly
$3,000$ OCs. In recent research, parts of the known OCs have coronae or haloes
(Meingast et al., 2021). According to Tarricq et al. (2022), the mean tidal
radius is $\sim 30$ pc as inferred from an analysis of more than $300$ OCs in
the solar neighborhood. This same study suggested a few OCs have
$r_{\mathrm{t}}>65$ pc. Leaving out those OCs with too large $r_{\mathrm{t}}$,
the mean $r_{\mathrm{t}}$ is still $\sim 24$ pc in Tarricq et al. (2022). Some
of our study’s cluster members do not have reddening and extinction estimates,
which are important to estimating the stellar mass in CAMD using isochrones.
The member candidates of c1, c2, c3 are clearly incomplete in
$M_{\mathrm{G}}$, as seen in Fig. 4. Therefore, $r_{\mathrm{t}}$ calculated
based on only these member candidates may be underestimated. In this work, we
adopted an average tidal radius $\bar{r_{\mathrm{t}}}\sim 24$ pc. Therefore,
clusters separated by less than $3\times\bar{r_{\mathrm{t}}}$ pc may be
tidally interacting. In Fig. 6, the green shading with different
transparencies are the regions where relative distances are in
$1\times\bar{r_{\mathrm{t}}}$, $2\times\bar{r_{\mathrm{t}}}$ and
$3\times\bar{r_{\mathrm{t}}}$ pc.

The distributions of member candidates in Cartesian coordinates $xyz$ are
presented in Fig. 7 to show the dimensions of the clusters and
$\bar{r_{\mathrm{t}}}$ more clearly. In this figure, member stars in the $xy$
plane appear to be stretched along the line of sight, primarily due to
uncertainties in distance, reddening and extinction. The solid circles are
$r_{\mathrm{t}}$ calculated from the total cluster mass of c1, c2 and c3,
using Eq. 2 according to Pinfield et al. (1998),

$\displaystyle r_{\mathrm{t}}$ $\displaystyle=$
$\displaystyle\left[\frac{G\mathrm{M_{total}}}{2\left(A-B\right)^{2}}\right]^{1/3}\,.$
(2)

$\mathrm{M_{total}}$ in Eq. 2 is the accumulated stellar mass for a cluster.
The other constants are : Gravitational constant $G=4.3\times 10^{-6}$
$\mathrm{kpc}$ $(\mathrm{km}$ $\mathrm{s}^{-1})^{2}$
$\mathrm{M}_{\mathrm{\odot}}^{-1}$; the Oort constants $A$ = 15.3 $\pm$ 0.4
$\mathrm{kpc}^{-1}$ $\mathrm{km}$ $\mathrm{s}^{-1}$, $B$ = -11.9 $\pm$ 0.4
$\mathrm{kpc}^{-1}$ $\mathrm{km}$ $\mathrm{s}^{-1}$ from Bovy (2017). The
dashed circles represent the above three increments in $\bar{r_{\mathrm{t}}}$.
Some members of c1 are within the range of $\bar{r_{\mathrm{t}}}$ of c2, and
vice versa. However, c3 is farther apart from the other two clusters and is
unlikely to be tidally influenced by them at the present time.

Figure 7: Distributions of member stars of c1, c2 and c3 in Cartesian
coordinates $xyz$ centered at Sun. The solid circles are the tidal radii based
on the total cluster mass, and the dashed circles are
$\bar{r_{\mathrm{t}}}\sim 24$ pc.  Figure 8: Distributions of [M/H] vs. $G$
for c1, c2, and c3. The colored lines present the median values of [M/H] for
the clusters.

###  4.4 Metallicity

In Sec. 3.3, the metallicity [M/H] was estimated from isochrone fitting :
$-0.3\pm 0.1$ dex for c1, c2, and $-0.1\pm 0.1$ dex for c3. Here, we use [M/H]
(mh gspphot) from Gaia DR3 to evaluate the metallicities of these clusters. We
adopt the upper and lower bounds of [M/H] to represent the uncertainties of
[M/H] and select $21$ member candidates of each cluster with the smallest
[M/H] uncertainties, which lowers the mean uncertainties of [M/H] to $0.02$
dex for each cluster. In Fig. 8 we present those stars in [M/H] vs. $G$, and
the median values of [M/H] are $-0.26$ dex for c1, c2, and $-0.12$ dex for c3,
indicated by the dashed lines.

##  5 Conclusion

In this section, we discuss the probable relationships among these three OCs
based on their separations, kinematics, ages and metallicities. The properties
of c1, c2 and c3 are summarized in Tab. 1.

When c1 was formed, c2 was very close to it. The relative distance between
them is less than $2\times\bar{r_{\mathrm{t}}}$ pc at present. The extensions
of $\bar{r_{\mathrm{t}}}$ for c1 and c2 have some areas overlapped at present,
suggesting their tidally interacting. Their
$\left(V_{\mathrm{R}},V_{\phi},V_{\mathrm{Z}}\right)$, also differ by less
than $3$ $\mathrm{km}\cdot\mathrm{s}^{-1}$. In addition, their ages and
metallicities are almost the same. We conclude that c1 and c2 are a primordial
binary cluster, formed simultaneously.

On the other hand, c3 was formed earlier and has a higher metallicity. In
addition, the current separation between c3 and c1 (or c2) and their relative
space velocities suggest that there was very little chance that c3 was closer
to c1 (or c2) than $\sim 3\times\bar{r_{\mathrm{t}}}$ since c3 was formed.
Also, most massive member candidates in c3 is no more than about $4$
$\mathrm{M_{\odot}}$. Therefore, it is unlikely that c3 triggered the
formation of the other two clusters. We consider c3 is not part of a multiple
system with c1 and c2. However, Fig. 6 suggests that tidal exchanges of stars
among c1, c2 and c3 may occur in the future.

##  6 Summary

We have identified three OCs, named as c1, c2, c3 in this paper. Using data
from Gaia DR3, we calculate the dereddened colors
$\left(G_{\mathrm{BP}}-G_{\mathrm{RP}}\right)_{0}$ and absolute magnitudes
$M_{\mathrm{G}}$. Isochrone fits to the CAMD indicates that the ages of c1,
c2, c3 are about $6.3\pm 0.9$, $8.1\pm 1.4$ and $21.8\pm 2.2$ Myr,
respectively. All the cluster member candidates have similar PMs. From a few
RVs in Gaia and LAMOST, we have estimated the mean RVs for these clusters. The
difference between c1 and c2 is less than $3$
$\mathrm{km}\cdot\mathrm{s}^{-1}$. The Galactocentric Cylindrical velocities
$\left(V_{\mathrm{R}},V_{\phi},V_{\mathrm{Z}}\right)=\left(-5.18,238.33,-7.17\right)$
$\mathrm{km}\cdot\mathrm{s}^{-1}$, $\left(-3.54,236.09,-6.91\right)$
$\mathrm{km}\cdot\mathrm{s}^{-1}$ and $\left(-8.42,239.61,-3.61\right)$
$\mathrm{km}\cdot\mathrm{s}^{-1}$ for c1, c2 and c3, respectively. The above
data were used to compute the orbits of each OC and the relative distances
among them. Using metallicities for individual OC stars from Gaia DR3, the
median values of [M/H] are $-0.26$ dex, $-0.26$ dex, and $-0.12$ dex,
respectively. From relative distances, kinematics, ages, and metallicities, we
conclude c1 and c2 comprise a simultaneously formed primordial binary OC. We
also expect there may be tidal captures among c3 and c1, c2 in the future.

Table 1: Properties of c1, c2 and c3 from our member candidates.

Name | $gl$ | $\sigma_{gl}$ | $gb$ | $\sigma_{gb}$ | $d$ | $\sigma_{d}$ | $\mu_{\alpha}^{*}$ | $\sigma_{\mu_{\alpha}^{*}}$ | $\mu_{\delta}$ | $\sigma_{\mu_{\delta}}$ | RV | $V_{\mathrm{r}}$ | $V_{\phi}$ | $V_{\mathrm{Z}}$ | $\tau$ | $\sigma_{\tau}$ | [M/H]
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
| deg | kpc | $\mathrm{mas}\cdot\mathrm{y}^{-1}$ | $\mathrm{km}\cdot\mathrm{s}^{-1}$ | Myr | dex
c1 | 31.28 | 0.36 | 5.25 | 0.25 | 0.426 | 0.008 | 3.24 | 0.31 | -8.59 | 0.43 | -8.05 | -5.18 | 238.33 | -7.17 | 6.3 | 0.9 | -0.26
c2 | 31.86 | 0.49 | 2.98 | 0.37 | 0.468 | 0.017 | 2.65 | 0.27 | -8.35 | 0.27 | -9.53 | -3.54 | 236.09 | -6.91 | 8.1 | 1.4 | -0.26
c3 | 27.17 | 0.39 | 2.59 | 0.23 | 0.377 | 0.009 | 2.06 | 0.20 | -9.00 | 0.26 | -4.15 | -8.42 | 239.61 | -3.61 | 21.8 | 2.2 | -0.12

Note. — $\sigma_{gl}$, $\sigma_{gb}$, $\sigma_{d}$,
$\sigma_{\mu_{\alpha}^{*}}$ and $\sigma_{\mu_{\delta}}$ are the standard
deviations from cluster members. $\tau$ is the age of a cluster.

This study is supported by the National Natural Science Foundation of China
under grant No. 11988101, 11973048, 11927804, 11890694 and National Key R&D
Program of China No. 2019YFA0405502. We acknowledge the support from the 2m
Chinese Space Station Telescope project : CMS-CSST-2021-B05. This work
presents results from the European Space Agency (ESA) space mission Gaia. Gaia
data are being processed by the Gaia Data Processing and Analysis Consortium
(DPAC). Funding for the DPAC is provided by national institutions, in
particular the institutions participating in the Gaia MultiLateral Agreement
(MLA). The Gaia mission website is https://www.cosmos.esa.int/gaia. The Gaia
archive website is https://archives.esac.esa.int/gaia. Guoshoujing Telescope
(the Large Sky Area Multi-Object Fiber Spectroscopic Telescope LAMOST) is a
National Major Scientific Project built by the Chinese Academy of Sciences.
Funding for the project has been provided by the National Development and
Reform Commission. LAMOST is operated and managed by the National Astronomical
Observatories, Chinese Academy of Sciences. TDO gratefully acknowledges
support from the U.S. National Science Foundation grant AST-1910396.

## References

  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33. doi: 10.1051/0004-6361/201322068
  * Bailer-Jones et al. (2021) Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., et al. 2021, AJ, 161, 147. doi: 10.3847/1538-3881/abd806
  * Bica et al. (2019) Bica, E., Pavani, D. B., Bonatto, C. J., et al. 2019, AJ, 157, 12. doi: 10.3847/1538-3881/aaef8d
  * Bovy (2017) Bovy, J. 2017, MNRAS, 468, L63. doi: 10.1093/mnrasl/slx027
  * Bovy (2015) Bovy, J. 2015, ApJS, 216, 29. doi: 10.1088/0067-0049/216/2/29
  * Bressan et al. (2012) Bressan, A., Marigo, P., Girardi, L., et al. 2012, MNRAS, 427, 127. doi: 10.1111/j.1365-2966.2012.21948.x
  * Campello et al. (2013) Campello, R. J. G. B., Moulavi, D., & Sander, J. 2013, in Advances in Knowledge Discovery and Data Mining, eds. J. Pei, V. S. Tseng, L. Cao, H. Motoda, & G. Xu (Berlin, Heidelberg: Springer, Berlin Heidelberg), 160
  * Cantat-Gaudin & Anders (2020a) Cantat-Gaudin, T. & Anders, F. 2020a, A&A, 633, A99. doi: 10.1051/0004-6361/201936691
  * Cantat-Gaudin et al. (2020b) Cantat-Gaudin, T., Anders, F., Castro-Ginard, A., et al. 2020b, A&A, 640, A1. doi: 10.1051/0004-6361/202038192
  * Cantat-Gaudin et al. (2018) Cantat-Gaudin, T., Jordi, C., Vallenari, A., et al. 2018, A&A, 618, A93. doi: 10.1051/0004-6361/201833476
  * Cantat-Gaudin et al. (2019) Cantat-Gaudin, T., Krone-Martins, A., Sedaghat, N., et al. 2019, A&A, 624, A126. doi: 10.1051/0004-6361/201834453
  * Castro-Ginard et al. (2022) Castro-Ginard, A., Jordi, C., Luri, X., et al. 2022, A&A, 661, A118. doi:10.1051/0004-6361/202142568
  * Castro-Ginard et al. (2020) Castro-Ginard, A., Jordi, C., Luri, X., et al. 2020, A&A, 635, A45. doi: 10.1051/0004-6361/201937386
  * Castro-Ginard et al. (2019) Castro-Ginard, A., Jordi, C., Luri, X., et al. 2019, A&A, 627, A35. doi: 10.1051/0004-6361/201935531
  * Castro-Ginard et al. (2018) Castro-Ginard, A., Jordi, C., Luri, X., et al. 2018, A&A, 618, A59. doi: 10.1051/0004-6361/201833390
  * Castro-Ginard et al. (2021b) Castro-Ginard, A., McMillan, P. J., Luri, X., et al. 2021b, A&A, 652, A162. doi: 10.1051/0004-6361/202039751
  * Chen et al. (2019) Chen, Y., Girardi, L., Fu, X., et al. 2019, A&A, 632, A105. doi: 10.1051/0004-6361/201936612
  * Chen et al. (2014) Chen, Y., Girardi, L., Bressan, A., et al. 2014, MNRAS, 444, 2525. doi: 10.1093/mnras/stu1605
  * Chen et al. (2001) Chen, B., Stoughton, C., Smith, J. A., et al. 2001, ApJ, 553, 184. doi: 10.1086/320647
  * Chen & Zhao (2020) Chen, Y. Q. & Zhao, G. 2020, MNRAS, 495, 2673. doi: 10.1093/mnras/staa1079
  * Cui et al. (2012) Cui, X.-Q., Zhao, Y.-H., Chu, Y.-Q., et al. 2012, Research in Astronomy and Astrophysics, 12, 1197. doi: 10.1088/1674-4527/12/9/003
  * de La Fuente Marcos & de La Fuente Marcos (2009) de La Fuente Marcos, R. & de La Fuente Marcos, C. 2009, A&A, 500, L13. doi: 10.1051/0004-6361/200912297
  * Dreyer (1888) Dreyer, J. L. E. 1888, MmRAS, 49, 1
  * Ester et al. (1996) Ester, M., Kriegel, H. P., Sander, J., A Density-Based Algorithm for Discovering Clusters in Large Spatial Databases with Noise. in Proceedings of the Second Interna- tional Conference on Knowledge Discovery and Data Mining, KDD‘ 96 (AAAI Press), page 226C231, 1996.
  * Ferreira et al. (2021) Ferreira, F. A., Corradi, W. J. B., Maia, F. F. S., et al. 2021, MNRAS, 502, L90. doi: 10.1093/mnrasl/slab011
  * Ferreira et al. (2019) Ferreira, F. A., Santos, J. F. C., Corradi, W. J. B., et al. 2019, MNRAS, 483, 5508. doi: 10.1093/mnras/sty3511
  * Gaia Collaboration et al. (2021) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2021, A&A, 649, A1. doi: 10.1051/0004-6361/202039657
  * Gaia Collaboration et al. (2018) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018, A&A, 616, A1. doi: 10.1051/0004-6361/201833051
  * Gaia Collaboration et al. (2016) Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016, A&A, 595, A1. doi: 10.1051/0004-6361/201629272
  * Gillessen et al. (2009) Gillessen, S., Eisenhauer, F., Trippe, S., et al. 2009, ApJ, 692, 1075. doi: 10.1088/0004-637X/692/2/1075
  * He et al. (2022) He, Z., Li, C., Zhong, J., et al. 2022, ApJS, 260, 8. doi:10.3847/1538-4365/ac5cbb
  * He et al. (2021) He, Z.-H., Xu, Y., Hao, C.-J., et al. 2021, Research in Astronomy and Astrophysics, 21, 093. doi: 10.1088/1674-4527/21/4/93
  * Hunt & Reffert (2021) Hunt, E. L. & Reffert, S. 2021, A&A, 646, A104. doi: 10.1051/0004-6361/202039341
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90. doi: 10.1109/MCSE.2007.55
  * Innanen et al. (1972) Innanen, K. A., Wright, A. E., House, F. C., et al. 1972, MNRAS, 160, 249. doi: 10.1093/mnras/160.3.249
  * Katz et al. (2022) Katz, D., Sartoretti, P., Guerrier, A., et al. 2022, arXiv:2206.05902
  * Kharchenko et al. (2013) Kharchenko, N. V., Piskunov, A. E., Schilbach, E., et al. 2013, A&A, 558, A53. doi: 10.1051/0004-6361/201322302
  * Krone-Martins & Moitinho (2014) Krone-Martins, A. & Moitinho, A. 2014, A&A, 561, A57. doi: 10.1051/0004-6361/201321143
  * Lada & Lada (2003) Lada, C. J. & Lada, E. A. 2003, ARA&A, 41, 57. doi: 10.1146/annurev.astro.41.011802.094844
  * Liang et al. (2017) Liang, X. L., Zhao, J. K., Oswalt, T. D., et al. 2017, ApJ, 844, 152. doi: 10.3847/1538-4357/aa7cf7
  * Lindegren et al. (2021a) Lindegren, L., Bastian, U., Biermann, M., et al. 2021a, A&A, 649, A4. doi: 10.1051/0004-6361/202039653
  * Lindegren et al. (2021b) Lindegren, L., Klioner, S. A., Hernández, J., et al. 2021b, A&A, 649, A2. doi: 10.1051/0004-6361/202039709
  * Lindegren et al. (2018) Lindegren, L. 2018, technical note GAIA-C3-TN-LU-LL-124
  * Liu & Pang (2019) Liu, L. & Pang, X. 2019, ApJS, 245, 32. doi: 10.3847/1538-4365/ab530a
  * Liu et al. (2015) Liu, X.-W., Zhao, G., & Hou, J.-L. 2015, Research in Astronomy and Astrophysics, 15, 1089. doi: 10.1088/1674-4527/15/8/001
  * McInnes et al. (2017) McInnes, L., Healy, J., & Astels, S. 2017, JOSS, 2, 205. doi: 10.21105/joss.00205
  * McKinney (2010) McKinney, W. 2010, in Proc. of the 9th Python in Science Conf. (SciPy 2010), ed. S. van der Walt & J. Millman, 56, doi: 10.25080/Majora-92bf1922-00a
  * Meingast et al. (2021) Meingast, S., Alves, J., & Rottensteiner, A. 2021, A&A, 645, A84. doi: 10.1051/0004-6361/202038610
  * Netopil et al. (2022) Netopil, M., Oralhan, İ. A., Çakmak, H., et al. 2022, MNRAS, 509, 421. doi: 10.1093/mnras/stab2961
  * Pedregosa et al. (2012) Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2012, https://arxiv.org/abs/1201.0490
  * Pinfield et al. (1998) Pinfield, D. J., Jameson, R. F., & Hodgkin, S. T. 1998, MNRAS, 299, 955. doi: 10.1046/j.1365-8711.1998.01754.x
  * Qin et al. (2021) Qin, S.-M., Li, J., Chen, L., et al. 2021, Research in Astronomy and Astrophysics, 21, 045. doi: 10.1088/1674-4527/21/2/45
  * Reid et al. (2014) Reid, M. J., Menten, K. M., Brunthaler, A., et al. 2014, ApJ, 783, 130. doi: 10.1088/0004-637X/783/2/130
  * Röser et al. (2019) Röser, S., Schilbach, E., & Goldman, B. 2019, A&A, 621, L2. doi: 10.1051/0004-6361/201834608
  * Röser & Schilbach (2019) Röser, S. & Schilbach, E. 2019, A&A, 627, A4. doi: 10.1051/0004-6361/201935502
  * Schönrich et al. (2010) Schönrich, R., Binney, J., & Dehnen, W. 2010, MNRAS, 403, 1829. doi: 10.1111/j.1365-2966.2010.16253.x
  * Sim et al. (2019) Sim, G., Lee, S. H., Ann, H. B., et al. 2019, Journal of Korean Astronomical Society, 52, 145. doi: 10.5303/JKAS.2019.52.5.145
  * Tarricq et al. (2022) Tarricq, Y., Soubiran, C., Casamiquela, L., et al. 2022, A&A, 659, A59. doi: 10.1051/0004-6361/202142186
  * Taylor (2005) Taylor, M. B. 2005, adass XIV, 347, 29
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science and Engineering, 13, 22. doi: 10.1109/MCSE.2011.37
  * van Leeuwen (2009) van Leeuwen, F. 2009, A&A, 497, 209. doi: 10.1051/0004-6361/200811382
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261. doi: 10.1038/s41592-019-0686-2
  * Yang et al. (2021) Yang, Y., Zhao, J., Zhang, J., et al. 2021, ApJ, 922, 105. doi: 10.3847/1538-4357/ac289e
  * Ye et al. (2021a) Ye, X., Zhao, J., Liu, J., et al. 2021a, AJ, 161, 8. doi: 10.3847/1538-3881/abc61a
  * Ye et al. (2021b) Ye, X., Zhao, J., Zhang, J., et al. 2021b, AJ, 162, 171. doi: 10.3847/1538-3881/ac1f1f
  * Zhang et al. (2021) Zhang, H., Chen, Y., & Zhao, G. 2021, ApJ, 919, 52. doi: 10.3847/1538-4357/ac0e92
  * Zhao et al. (2009) Zhao, J., Zhao, G., & Chen, Y. 2009, ApJ, 692, L113. doi: 10.1088/0004-637X/692/2/L113
  * Zhao et al. (2012) Zhao, G., Zhao, Y.-H., Chu, Y.-Q., et al. 2012, Research in Astronomy and Astrophysics, 12, 723. doi: 10.1088/1674-4527/12/7/002
  * Zhao et al. (2006) Zhao, G., Chen, Y.-Q., Shi, J.-R., et al. 2006, Chinese J. Astron. Astrophys., 6, 265. doi: 10.1088/1009-9271/6/3/01