# Galaxy flows within 8,000 km s-1 from Numerical Action methods

Edward J.Shaya  Astronomy Department, University of Maryland, College Park, MD
20742, USA<EMAIL_ADDRESS>R. Brent Tully  Institute for Astronomy,
University of Hawaii, 2680 Woodlawn Drive, Honolulu, HI 96822, USA  Daniel
Pomarède  Institut de Recherche sur les Lois Fondamentales de l’Univers, CEA
Université Paris-Saclay, 91191 Gif-sur-Yvette, France  Alan Peel  Astronomy
Department, University of Maryland, College Park, MD 20742, USA:
<EMAIL_ADDRESS>
(Accepted to ApJ)

###### Abstract

The trajectories since $z=4$ of systems of galaxies (‘halos’) with $cz<8,000$
km s-1 are found through Numerical Action reconstructions. A set of 9,719
halos from a 2MASS group catalog and Cosmicflows-3 catalogs are given
attention. Present distances are adjusted to minimize departures from observed
redshifts. For those with the most precisely determined distances, compromises
are made between distance and redshift agreement. $H_{0}$ is varied from 69 to
77 km s-1 Mpc-1 with $\Omega_{m}$ set by the baryon acoustic oscillation
constraint from the Planck Satellite. A best fitting amplitude of the mass-to-
light relation is found. A uniform density associated with the interhalo
medium accounts for the matter not in halos. The solution paths provide the
histories of the formation of the nearby large structures and depict how the
voids emptied. Assuming no local over/underdensity, the best model has
$H_{0}=73$ km s-1 Mpc-1 with nearly the same density arising from interhalo
matter (IHM) as from halos. We examine local over/underdensities by varying
the IHM density and find a valley of best fit models along
$H_{0}=73.0(1+0.165\delta)$ km s-1 Mpc-1. Friedmann models with distinct
densities internal and external to the study region give a similar
relationship. The fraction of matter in the IHM seen in n-body simulations
roughly matches that in our $H_{0}=72$ scenario. Videos have been created to
visualize the complexities of formation of large-scale structures. Standard
n-body calculations starting from the first time-steps as tests of the NAM
solutions, and continue until cosmic scale factor $a=2$ provide glimpses into
the future.

##  1 Introduction

A grand goal in cosmology would be to know the full trajectories of all the
pieces of all the galaxies within a volume of the universe around us large
enough to be considered representative. Imaginary macrocosms can be simulated
(Springel et al., 2006, 2018; Vogelsberger et al., 2014; Schaye et al., 2015).
There are simulations based on constrained initial conditions that can
generate major clusters ($>10^{14}M_{\odot}$) and voids at locations that
approximate the observed nearby universe (Klypin et al., 2003; Sorce et al.,
2016). However stochastic processes prevent the accurate placements of
entities on scales of small groups and individual galaxies that accurately
mimic the real world.

Peebles (1989) introduced an alternative way forward. The comoving coordinate
orbits of particles follow paths that are extrema or saddle points of the
action, $S$, the integral of the Lagrangian, $\mathscr{L}$, over time, $t$:

$\delta S=\delta\int_{0}^{t_{0}}\mathscr{L}dt=0$ (1)

where for collisionless particles the Lagrangian is the comoving kinetic minus
the comoving potential energy in the system:

$\mathscr{L}=T-V$ (2)

The boundary constraints include the requirement that all initial peculiar
velocities are small. There have been various studies employing numerical
action methods to model galaxy orbits (Shaya et al., 1995; Peebles et al.,
2001; Shaya & Tully, 2013; Shaya et al., 2017). Nusser & Branchini (2000), in
particular, studied the accuracy of NAM and found that it reconstructs
velocities faithfully in high density regions where deviations from the
Zel’dovich solution are large.

The ability to model the non-linear dynamics of infall regimes is an important
advantage of the numerical action method. On scales of a few Mpc, most
galaxies and groups are, at least slightly, in a non-linear interaction with
nearby masses. The non-linear regime infall of test particles provide strong
constraints on the masses of major attractors over a wide range of radii. The
ubiquitous voids, being underdensities of $\delta=(\rho/\rho_{m}-1)\sim-1$,
obviously generate motions at their peripheries (i.e., most everywhere) not in
the linear regime, defined as $|\delta|\ll 1$.

The orbits are followed from an early time, here from $z=4$, age of the
universe $\sim$1.6 Gyr. Back then, the entities we see today were in many
pieces. The premise that is adopted is that the orbits trace the centers of
mass of those many pieces that have assembled into the current halos, and the
center of mass is unaffected by the details of the merger tree. Limitations of
orbit reconstructions must be acknowledged. For one thing, there is usually
insufficient information to untangle the orbits through multiple shell
crossings in collapsed regions.111The case of the Local Group might be the
exception (Shaya & Tully, 2013), especially as proper motions become more
reliably established. Hence, the units that will be considered in the
construction of our sample catalog are collapsed halos that often include a
multitude of galaxies. The dominion of a halo extends to roughly the virial
radius. Entities that are infalling, not yet at shell crossings, are kept
distinct.

The viability of the model depends on whether its sample catalog reasonably
captures the distribution of matter. Calculations become oppressive if the
number of entities is large, but a trial solution of $\sim$10,000 entities can
now be found by high end desktop computers in less than a day. There can be
considerable work afterward, though, untangling tight interactions and
choosing among multiple path solutions.

There should be as close to complete three-dimensional coverage as possible
around the origin, in this case the Local Group. However, there is missing
information in the zone of obscuration on the Galactic equator. Given the
computational penalty for carrying excessive entries, very small galaxies and
groups are excluded. Then the matter outside of the luminous halos must be
taken into account. Large n-body simulations show that the amount of matter in
sheets, filaments, halos too small to form and retain significant numbers of
stars, or free streaming is comparable to the amount of matter in luminous
bound systems. This assertion is supported by a study of the unbound particles
in the Millennium Run n-body calculation (Shattow & Croton, 2015). Studies of
voids show that they are not completely empty of mass. They are threaded by
networks of matter with average density at a few tenths of mean density
(Gottlöber et al., 2003).

There are two ways that entities (virialized halos or individual galaxies) can
enter the list considered for orbit reconstruction. The first is because they
cross an intrinsic luminosity threshold that infers they are important
constituents of the mass distribution. The second is because they have
measured distances that are considered accurate above a threshold. The more
entries with accurately known distances the better. These are test particles
responding to and informing us of the mass distribution.

The numerical action methodology generates orbits that are physically possible
within the constraints of the distribution of constituents. However, there may
be a number of plausible options, particularly in crowded regions. The
specifics of individual orbits are tenuous. However, mass estimates of the
major components and the overall density are reasonably robust (within the
specifics of the assumed cosmology). Also, the overall flow patterns are
robust.

The current study is an extension of Shaya et al. (2017) where the volume
considered was limited at 38 Mpc $\sim$2,850 km s-1 and 1,382 orbits were
followed. That volume encompassed the traditional Local Supercluster (de
Vaucouleurs, 1953) with the Virgo Cluster the dominant mass concentration and
the Local Void a dominant dynamical factor (Tully et al., 2019). In this new
study, orbits are reconstructed to a limit of 8,000 km s-1 in the CMB frame,
$\sim$107 Mpc, a volume 22 times greater, and 9,719 entities are given
attention.

There were two considerations in establishing the volume limit. One is the
practical reason that information on distances is thinning out around this
limit in the base Cosmicflows-3 catalog (Tully et al., 2016). The other is
related to the known structural elements in our vicinity. The 8,000 km s-1
limit (CMB frame) comfortably encloses both the so-called Great Attractor
region (Dressler et al., 1987a; Lynden-Bell et al., 1988), the core of
Laniakea Supercluster (Tully et al., 2014), and the Perseus-Pisces complex
(Haynes & Giovanelli, 1988). See Figure 1. These two high density regions are
given particular attention in this study. A third complex, the Great Wall (de
Lapparent et al., 1986), crosses the 8,000 km s-1 boundary and any analysis of
it is tentative. Aside from the Great Wall structure, the 8,000 km s-1
boundary generally runs through uncomplicated void-like regions.

Figure 1: Distribution of 2MASS galaxies within 100 Mpc in a projection onto
the supergalactic equator. Major clusters are given colors; red if in the
Laniakea Supercluster; blue if associated with the Perseus-Pisces filament;
and green if in the Great Wall. The zone of Galactic obscuration projects into
the $\pm 5^{\circ}$ wedge.

In §2 there will be a presentation of the input data to the model. In §3 there
will be a discussion of the numerical action methodology and the cosmology
assumptions. §4 presents statistical results from the numerical action study.
Cosmographic representations of results are given in §5. In §6 the evolution
of structure is followed into the future. There are brief conclusions in §7.

##  2 Data

The numerical action calculation minimally needs constituents that provide
good representation of the large scale clustered mass distribution through
apparent brightnesses and redshifts. Most entities we will consider have
redshifts with high S/N and sufficient accuracies. In addition, to invoke
goodness criteria, good distances are needed for a substantial fraction of the
entities, and these should be well distributed over the entire sample volume.
This study has awaited the availability of the copious assemblage of distance
determinations given in the Cosmicflows-3 compilation (Tully et al., 2016).

###  2.1 Mass distribution

The description of the mass distribution draws on two galaxy group catalogs.
One draws on the collection of redshifts (Huchra et al., 2012) for 2MASS, the
Two Micron All Sky Survey (Skrutskie et al., 2006) that is quasi-complete to
$K_{s}=11.75$ outside a small strip along the Galactic plane at
$|b|<5^{\circ}$ . The creation of a 2MASS-based group catalog began with a
study of selected well studied groups and clusters over the wide range of
masses from $10^{11}$ to $10^{15}M_{\odot}$ (Tully, 2015a). That study
provided scaling relationships between luminosity, group velocity dispersion,
$\sigma_{p}$, the characteristic radius of second turnaround, $R_{2t}$, and
group mass.

The scaling relationships provide the tools to build a group catalog based on
the 2MASS redshift compilation to $K_{s}=11.75$ (Tully, 2015b). The 2MASS
compilation is close to optimal because of its photometric integrity and its
sensitivity to older stellar populations that represent most of the baryonic
mass. Most of the mass in stars in the domain of the local universe that
concerns us (outside the gap on the Galactic plane) is probably well
represented in the 43,000 galaxies in the 2MASS redshift catalog.

Figure 2: Percent inclusion of mass associated with the 2MASS $K_{s}<11.75$
redshift catalog as a function of the mass cut on groups. Numbers at intervals
along the black curve give the cumulative numbers of groups as a function of
mass. The red curve gives equivalent information for the normalized mass
within the $\pm 5^{\circ}$ obscured zone.

The entries in the group catalog are collapsed halos that can range from rich
clusters with hundreds of associated 2MASS galaxies down to individual
isolated galaxies. To begin the construction of the catalog, the intrinsically
brightest galaxy in the 2MASS collection (assuming Hubble distances) was given
attention and gathered within its domain all other galaxies as expected from
the scaling relationships. Specifically, galaxies are linked if their
velocities agree within $2\sigma_{p}$ and their projected positions lie within
$R_{2t}$. The inclusion of new members dictates iterations based on the
increased luminosity. The masses statistically associated with halo
luminosities are taken to grow with a weak power of luminosity: $M\propto
L^{1.15}$. Once the composition of the halo associated with the brightest
galaxy is secured, the most luminous galaxy in the remainder of the 2MASS
catalog is given consideration, and so on, until all the galaxies have been
linked to a halo, including halos of a single member.

The cumulative mass fraction as a function of halo mass is shown in Figure 2.
It is seen that almost all the mass in collapsed halos is contained within
halos above $10^{12}~{}M_{\odot}$. There are $\sim$10,000 such halos within
8,000 km s-1.

There are details, such as arise when halos are in close proximity and may be
merged or split, and instances of what must be high velocity members. These
complexities are discussed by Tully (2015b). An issue of relevance is
adjustments that are made for lost luminosity as a function of distance. The
assumption is made that the summed luminosity in radial shells after
adjustment should be constant. A formula for a correction factor was
empirically determined to achieve this criterion (see Fig. 5 and Eq. 9 in
Tully (2015b)).

By 10,000 km s-1 the apparent magnitude cutoff of $K_{s}=11.75$ for the
catalog is in the vicinity of the absolute magnitude $M^{\star}_{K_{s}}$, at
the exponential break in the Schechter function (Schechter, 1976) (taking
$M^{\star}_{K_{s}}=-23.55+5{\rm log}$h, h = $H_{0}/100$ km s-1 Mpc-1, and
slope parameter $\alpha_{K_{s}}=-1.0$). Correction factor adjustments are
below a factor two up until this velocity but then rapidly increase with
redshift. Corrections for missing light are moderate within our limit of 8,000
km s-1.

The 2MASS based group catalog is less than ideal nearby because the space
density of the high surface brightness galaxies that enter the 2MASS
compilation is low and because redshifts break down as a way of establishing
relative galaxy distances. Nearby, within 3,500 km s-1, there is the
preferable group catalog of Kourkchi & Tully (2017). The groups in their
catalog are based on the same scaling relationships. It includes all 15,000
galaxies that were known within the velocity limit. For many small galaxies,
$K_{s}$ magnitudes had to be inferred from optical fluxes (Jarrett et al.,
2003; Karachentsev & Kutkin, 2005). Measured distances were used when
available, allowing the rejection of group interlopers and giving improved
clarity to group memberships down to $10^{10}M_{\odot}$.

There remains the problem of the zone of obscuration, a gap in the 2MASS
compilation of 9% of the sky. Locally, within 2,850 km s-1, this problem was
addressed by Shaya et al. (2017) with the addition of sources from the main
2MASS catalog (Jarrett et al., 2000) at $|b|<5^{\circ}$ combined into groups.
Consideration of blind HI surveys (Kerr & Henning, 1987; Staveley-Smith et
al., 2016) gave assurance that important obscured structures were given
representation and these same elements are included in this new study.

Beyond 2,850 km s-1 out to 8,000 km s-1, the $\pm 5^{\circ}$ Galactic
equatorial gap widens in physical space and there must be compensation for
what is missed. There are known clusters and individual galaxies at these low
latitudes from x-ray, infrared, and radio surveys and follow up redshift
observations (Kraan-Korteweg, 2005; Kocevski et al., 2007). These entities are
given reasonable mass estimates and included. However, these known pieces are
insufficient, as judged by the subsequently averaged density at low latitudes.
To compensate, their masses are amplified by a factor four to recover
densities associated with halos roughly equal to values in unobscured regions.
It is supposed that their inflated masses account for missing mass in their
proximity. Clearly, the orbits found for them and those adjacent should be
considered unreliable.

###  2.2 Accurate distances

Galaxy distance measures are drawn from the Cosmicflows-3 compilation (Tully
et al., 2016). Contributions come from six methodologies. Two were
foundational for the absolute scale: the Cepheid period-luminosity relation
(Leavitt & Pickering, 1912; Riess et al., 2019) and the tip of the red giant
branch method (Lee et al., 1993; Rizzi et al., 2007). The latter, TRGB, is
particularly noteworthy because the $\sim$500 distances with the $\sim$5%
accuracy it has provided results in an outstandingly precise view of structure
and peculiar motions within 10 Mpc.

Two other methodologies provide distance estimates in a wide $30-300$ Mpc
range that are accurate at the level of 10% but are limited in number: surface
brightness fluctuations (Tonry & Schneider, 1988; Jensen et al., 2015) and
Type Ia supernovae (Phillips, 1993). It can be anticipated that contributions
from these techniques will be increasingly relevant in the future.

It is the two remaining methodologies that contribute the most numerically
although their accuracies per galaxy are only $20-25\%$: the early-type
fundamental plane (Dressler et al., 1987b; Djorgovski & Davis, 1987; Magoulas
et al., 2012) and spiral galaxy luminosity-linewidth (Tully & Fisher, 1977;
Kourkchi et al., 2020) relations. These techniques contribute distances to
$\sim$97% of the 17,647 entries in Cosmicflows-3. Within 8,000 km s-1, the
concern of this study, the coverage around the sky outside the zone of
obscuration is reasonably uniform.222At $8,000-15,000$ km s-1 the coverage in
Cosmicflows-3 is greater in the celestial south because of the important
contribution from the fundamental plane component of the Six Degree Field
survey (Springob et al., 2014).

The accumulation of distance estimates within groups and clusters is important
beyond the linkage of methodologies. The separate distance moduli determined
by different methods for a given galaxy or for members of a group can be
averaged, with appropriate weights. The entities in our catalog are assigned
velocities averaged over all known members and distances averaged in the
modulus over all members with suitable measurements. Consequently, the
distances to some halos have high statistical accuracy. In the extreme, for
the Virgo Cluster there are 170 measures. The statistical error on the cluster
distance is 1%. Systematic errors dominate. Figure 3 provides the number of
halos with better than a given uncertainty in distance modulus. The 924 halos
with distance errors $\mu_{e}<0.30$ mag (left of vertical dotted line)
nominally satisfy our criterion for model goodness.

Figure 3: Cumulative number of distance modulus errors for halos in our input
catalog.

###  2.3 Influences beyond 8,000 km/s

It can be expected that external tidal influences are a concern for the
present study, particularly toward the edges of the sample and perhaps for the
overall dipole motion within the volume. Unfortunately, our current knowledge
of structure beyond 8,000 km s-1 is crude and depends on the direction. The
reasonably good information to $\sim$16,000 km s-1 in the celestial south
provided by the Six Degree Field Redshift Survey (Jones et al., 2009) is
insufficiently balanced by our more limited knowledge of the celestial north.
Fortunately, the bulk dipole motion of the volume plays no role in the paths
determined in this work as the paths are in the frame of the center of mass of
the sample. Furthermore, because this study has a much larger radius than was
used in Shaya et al. (2017), the higher order gravitational influences are
greatly diminished.

We take recourse once again to the 2MASS redshift survey to $K_{s}=11.75$
(Huchra et al., 2012) that gives homogeneous representation of the entire sky,
save for the Galactic plane cutout, and to the group catalog built from that
survey (Tully, 2015b). We represent overdense structure in a shell
$8,000-16,000$ km s-1 in the CMB frame by the 44 clusters in the group catalog
with masses greater than $10^{15}M_{\odot}$. Clearly this description is
approximate. Selection corrections to masses are very large at these redshifts
and potentially significant components might be entirely missed. We work on
the hypothesis that the richest clusters delineate regions of the most
substantial overdensities. We hold these elements fixed in co-moving
coordinates and following linear theory increase their mass with time to reach
masses observed today.

##  3 Numerical Action Method

The Numerical Action Method (NAM) for determining the orbits of the components
of the large scale structure over cosmological time invokes the variational
principle in comoving coordinates as described in Peebles (1989, 1990, 1994),
wherein the physical trajectories are the paths of the extrema and saddle
points of the action integral. All of the derivatives of the action integral
with respect to variations of these paths at discrete points in time result in
zero. Positions trace the center of mass over time of the particles that today
are bound to each halo.

The boundary conditions at the first time step in the action hold that the 3d
velocities are low compared to $H_{0}d$. This condition is met because the
velocities, as determined by the differences in positions at the first two
time steps, following linear theory, are proportional to the gravitational
accelerations at the first time step. The boundary conditions at the final
time steps are either the 3d positions or the 2d sky position plus redshifts.

One can start with a random walk of positions backwards from the last step for
all particles and take turns improving the orbits by making small positional
changes that reduce the sum of the square of all the derivatives (SOS)
involving that particle. As this sum is driven to as close to zero as
possible, a process called ‘path relaxation’, the path and the resultant
redshift or final position converge to ever more precise values. The
directions and amplitudes of positional changes that relax a path is informed
by having the closed form second derivatives of the action integral. Each
particle is treated separately in turn while ignoring its effects on the other
particles until the next one is to be treated. Many iterations over all
particles are required to reach consistency of the mutual interactions.

###  3.1 Path Calculations

There is a wide range of precision in our knowledge of the distances to halos
in our catalog with the majority having no velocity independent distance
determinations. On the other hand, we typically have redshifts with small
uncertainties. A model can be built by using redshifts to determine distances
and, for those entities that have good distance measurements, the model
distances can be used to determine a goodness metric. We find paths that have
best agreement between model redshift, $cz_{m}$, and observed redshift,
$cz_{o}$.

It is imperative that the model positions are not free to drift too far from
their measured distances because that would likely imply a gravitational field
that is off from reality. We therefore treat the halos with the highest
precision distances differently. When we have distances with $\mu_{e}\leq 0.2$
mag, we chose distances that minimize a goodness metric that weighs the
ability of the model to reproduce both observed distances and redshifts.

We use the New NAM (NNAM) algorithms described in Peebles (2010) for speed.
This procedure solves for redshifts given position and then, as desired, we
allow individual distances to ‘drift’ to improve the redshift agreement. When
searching for simple redshift agreement, we solve for the paths at two
different distances to determine the rate of change of redshift with distance,
and calculate the size of the step to take based on the current disagreement.
This procedure usually converges quickly since $cz$ is nearly linear with
distance over most of space. When $cz$ agreement is better than 4 km s-1, we
stop. Determination of a best fit (distance, velocity) pair is substantially
more time consuming. Starting at the measured distance, we take a 0.3 Mpc step
and grow the step size by 4% with each step that improves the fit. When a step
causes $\chi^{2}$ to rise, the direction reverses and the step size is halved.
When the step size is reduced to 0.003 Mpc, we stop. But, because there are
triple-valued regions that can trap the final position at the wrong location,
we also try two more sequences starting at both $\pm 3\sigma$ locations in
distance.

Two issues complicating solutions are that sometimes NAM finds a very complex
path that is not warranted and, at some locations, there may not be a
solution. For complex paths, we go through a procedure discussed in § 3.4 on
untying knots. For locations with no solution, our procedure is to continue
stepping along until a solution can be found and consider this the start of a
new sequence. Throughout, we keep track of the best fit path thus far.

As in the procedure to find a first solution at initial positions, the
distance adjustment procedure must be repeated on all halos until no halo is
making significant excursions. Typically the whole set of orbits reaches a
stable solution after a dozen or so iterations.

###  3.2 Goodness criteria

The goodness metric applied to halos is:

$\chi^{2}=\left(\frac{5\log_{10}({d_{m}}/{d_{o}})}{\mu_{e}}\right)^{2}+\left(\frac{cz_{m}-cz_{o}}{cz_{e}}\right)^{2}$
(3)

where $d_{m}$ and $d_{o}$ are model and observed distances. The error in
redshift is taken to be $cz_{e}=60$ km s-1 to account for the possibility that
the average velocity of dark matter in a cluster may differ from the average
of the cluster galaxies and because the unweighted average that we use
(because the data is incomplete) may differ from the mass weighted value
appropriate for center of mass motion.

Since there are multiple plausible paths, care must be taken when comparing
⟨$\chi$⟩ of different scenarios to ensure that the variation is not just due
to different paths taken by random chance. As part of our procedure, we
compare individual $\chi=\sqrt{\chi^{2}}$ values in the most similar previous
scenarios to see which $\chi$ values have made big changes for the worse.
Then, by an interactive procedures described in § 3.4, attempts are made to
see if the better path can be achieved in the present scenario.

###  3.3 Setting up the parameters

We use 100 time steps to go from $z=4.0$ to $z=0$, evenly spaced in the scale
factor $a(t)$. On the scale of the distances between entities in this model,
motions earlier than $z=4$ are well described by linear perturbation theory,
but at later times non-linear interactions are very common. We solve for
cosmological scenarios keeping $H_{0}$ in the range from 69 to 77 km s-1
Mpc-1. The Planck Satellite measurement of the acoustic angular scale on the
sky sets the coupling between the mean density of matter in the universe and
the Hubble Constant: $\Omega_{m}h^{3}=\xi_{P}=0.096345$ (Planck Collaboration
et al., 2020).333$\Omega_{m}$ is the density of matter, $\rho_{m}$, with
respect to the critical density for a universe closed by matter, $\rho_{c}$.
The calculations assume the standard flat $\Lambda$CDM universe (i.e.,
$\Omega_{\Lambda}+\Omega_{m}=1$, where
$\Omega_{\Lambda}\equiv\frac{c^{2}\Lambda}{3H_{0}^{2}}$ and $\Lambda$ is the
cosmological constant) with negligible non-gravitational forces or dynamical
friction. The measured distances and redshifts are never altered. That is,
while $H_{0}$ is uncertain, we do know distances (with uncertainties) within
the constraints of our zero point calibration.

The $K-$band luminosity to mass mapping for halos is set by Eq. 7 of Tully
(2015b):

$M=\kappa\left(\frac{L_{K}}{10^{10}L_{\odot}}\right)^{1.15}10^{10}~{}M_{\odot}.$
(4)

The coefficient $\kappa$ in this relation was found to be $43h$ in that study.
The linear dependence on $H_{0}$ arises because masses are generally derived
from equilibrium dynamics that are independent of time, but proportional to
distance, while luminosity is proportional to distance squared. Hence, mass
divided by light is inversely proportional to $H_{0}$. In this study,
distances are fixed by their velocity independent measurements and masses are
not in equilibrium configurations, so a more appropriate dependence of
$\kappa$ on $H_{0}$ needs to be found.

To account for mass not bound to the regions illuminated by major galaxies,
such as in voids, sheets, filaments, low mass halos, or streaming particles, a
uniform density interhalo medium,
$\Omega_{\mathrm{IHM}}\equiv\rho_{\mathrm{IHM}}/\rho_{c}$, is added to NAM.
With the transformation from physical to comoving coordinates, there is a
fictitious outward acceleration term on all particles from empty space
(Peebles, 1993, Eq. 20.27) given by the second term in the expression for
peculiar gravitational acceleration

$\mathbf{g_{i}}=\frac{G}{a^{2}}\sum_{j}m_{j}\mathbf{\frac{x_{j}-x_{i}}{|x_{j}-x_{i}|^{3}}}+\frac{1}{2}\Omega_{m}H_{0}^{2}a\mathbf{x_{i}}.$
(5)

To implement the IHM component, the second term is modified to

$\frac{1}{2}\Bigl{(}\Omega_{m}-\frac{a\delta_{0}+1}{\delta_{0}+1}\Omega_{\mathrm{IHM}}\Bigr{)}H_{0}^{2}a\mathbf{x_{i}}.$
(6)

This formulation assumes that any present overdensity of the study volume,
$\delta_{0}$, grew to its current amplitude with time dependence a(t). The IHM
needed in our models must be primarily dark matter since the baryon density is
known to be limited to about $0.021<\Omega_{b}h^{2}<0.025$ (Fields et al.,
2014).

###  3.4 How to untie a knotted orbit

Because we are dealing with many orbits, fairly rare events can be quite
numerous. A problem arises when neighboring paths undergo very tight
interactions. In reality, such interactions would typically result in a merger
after some sloshing back and forth of the orbits. But, since our calculations
do not include dynamical friction, mergers do not occur. Typically, a
completed run would have several dozen orbits tied up. To avoid these unwanted
interactions, code was added to help find and interactively untie knots. A
menu of tools was added to the NAM code to allow the user to select individual
paths and choose to 1) set a new distance, 2) find paths backward in time from
present position and redshift (‘backtrack’), 3) check on the quality of the
path ($\chi$ and SOS), 4) relax the path (reduce the SOS) to a prescribed
value, 5) step a chosen halo in distance until model redshift agrees with
observed, 6) write out path solutions now, in case one is about to make things
worse.

The first step in the procedure to untie an interaction is to ‘backtrack’ at
least one of the pair. A backtrack consists of stepping back in time from the
present with a velocity consistent with the observed redshift given a spatial
direction. The entire $4\pi$ radians of present velocity directions are
sampled at $\sim$5° intervals, and only paths with SOS below a set threshold
are displayed. One of these is chosen, usually the one with lowest SOS if it
is not a complicated orbit. This selection unties the two paths, and sometimes
is close to an actual solution, although usually the SOS is still far from its
target goal. Next, one returns to the distance adjustment procedure and allows
each to drift until the proper redshift is met. If the two are still
interacting even mildly, the distance adjustments need to be iterated until
the paths are stable.

Additional IDL code was written to locate pairs that are tightly interacting.
The program finds, for each time step, the pairs with the smallest
separations. Some pairs come very near to each other but do not suffer a
strong interaction, so a list of halos to exclude can be provided. After
untangling a full set of interacting pairs, the knot finder is run again to
identify any remaining interactions.

Yet another method to find tight interactions that we used with each scenario
is to run an n-body forward integration, with a large number of timesteps,
starting with the initial NAM position and velocity, and to compare the final
positions with the NAM solution. Tight interactions show up as large
disagreements ($>0.3$ Mpc) in the final positions because the NAM time steps
are too large to accurately follow close interactions.

##  4 Results

We began our exploration by setting the density of the sample to the mean
density of the universe. For each of $H_{0}=69$, 73, and 77 km s-1 Mpc-1,
three different amounts of galaxy density were tested,
$\Omega_{\mathrm{halos}}/\Omega_{m}=0.3$, 0.5, and 0.7. To attain this, first
$\Omega_{\mathrm{IHM}}$ was set to 0.7, 0.5, and 0.3 of $\Omega_{m}$. Then the
coefficient $\kappa$ on the luminosity to light mapping was altered until
$\delta=\Omega_{\mathrm{local}}/\Omega_{m}-1\leq\pm 0.01$. If, at the end of
an iteration over all halos, $\delta$ wandered by more than 1%, $\kappa$ was
adjusted. The resulting distribution of mass density within spheres centered
on the Local Group is shown in Figure 4 for the case of $H_{0}=73$,
$\Omega_{\mathrm{halos}}/\Omega_{m}=0.5$. From 40 to 70 Mpc there is a net
underdensity but the density settles to a constant by 70 Mpc. The other cases
gave similar distributions.

Figure 4: Components of density with distance - The sum of the mass density
with respect to critical density within radius R from the Local Group in
halos, $\Omega_{\mathrm{halos}}$, plus the IHM, $\Omega_{\mathrm{IHM}}$, add
up to $\Omega_{\mathrm{local}}$ and is close to $\Omega_{m}$ at 100 Mpc in
model with $H_{0}=73$, $\Omega_{m}=0.248$, and $\Omega_{\mathrm{IHM}}=0.124$.
The jump at 16 Mpc arises from the inclusion of the Virgo Cluster.

For each scenario, we calculated the goodness metric for the precision
distance halos ($\mu_{e}\leq 0.3$ mag) that are presented in Table 1: Column 5
is the average value of $\chi$; Column 6 is the median $\chi$. Columns 7 and 8
have the average and median $\chi_{d}$ of just the distance term in Eq. 3.
These values indicate, in the best cases, that the NAM distances for the
roughly 8,000 halos with poor or no distance measurements have fairly high
reliability and should be preferred to their Hubble distances. Column 9 is the
Virgo Cluster $\chi$ to be discussed in § 4.6 and Column 10 has the number of
halos with $\chi<7$ used in the statistics.

Figure 5 shows the distribution of $cz$ agreement for the 9,157 halos either
lacking a distance determination or with a distance uncertainty $\mu_{e}>0.2$
mag in our fiducial run with $H_{0}=73$, mass split evenly between halos and
the IHM, and $\delta=0$. We quit seeking better $cz$ agreement when the
agreement is $<4$ km s-1. Figure 6 shows the cumulative distribution of $cz$
disagreement for the 527 halos with higher precision distances ($\mu_{e}\leq
0.2$ mag) that have present distances adjusted to minima in $\chi$ (fiducial
case). The vast majority found locations where the $cz$ disagreement is $<70$
km s-1. A few halos fail to find $cz$ agreement near their observed distances,
but there are many reasons why that may happen. For instance, heavily grouping
galaxies can rob a galaxy of a nearby galaxy critical to attaining its
$cz_{o}$ near its true distance. Also, in regions around massive structures,
there can be multiple solutions with the same radial projected redshifts.
However, for the majority of the halos only one reasonable solution is found,
probably because most are in the nearly linear perturbation regime.

Figure 7 presents the distribution of $\chi$ for elements with high precision
distances for the fiducial case. The distribution is close to that of a half-
normal distribution (red curve) with the same number of measurements and
variance. Figure 8 presents the individual $\chi$ values as a function of
$cz_{o}$ for that case, with horizontal lines at the mean and median. Values
below 0.1 are shown at 0.1 to keep them in the plot.

Figure 5: Histogram of differences between model and observed redshifts of all halos that have $\mu_{e}>0.2$ mag or no $\mu_{e}$ ($H_{0}=73$ run). Figure 6: Cumulative distribution of differences between model and observed redshifts for halos placed at minima in $\chi$ ($\mu_{e}\leq 0.2$). Table 1: Run Parameters and Output Input | Output
---|---
$H_{0}$ | $\Omega_{m}$ |  $\vrule depth=6.0pt,height=12.0pt,width=0.0pt\Omega_{\mathrm{IHM}}$ | $\delta$ | $\langle\chi\rangle$ | $\mathrm{Median}(\chi)$ | $\langle\chi_{d}\rangle$ | $\mathrm{Median}(\chi_{d})$ | $\chi(Virgo)$ | N
$\delta\simeq 0$ runs
69 | 0.2933 | 0.205 | -0.000 | 1.230 | 0.990 | 1.117 | 0.844 | 2.620 | 1195
69 | 0.2933 | 0.147 | -0.001 | 1.217 | 0.989 | 1.101 | 0.847 | 3.639 | 1195
69 | 0.2933 | 0.088 | 0.000 | 1.315 | 1.071 | 1.200 | 0.965 | 4.304 | 1196
71 | 0.2692 | 0.150 | 0.003 | 1.150 | 0.901 | 1.042 | 0.797 | 1.930 | 1195
73 | 0.2477 | 0.173 | 0.000 | 1.159 | 0.943 | 1.051 | 0.814 | 1.161 | 1196
73aafiducial scenario | 0.2477 | 0.124 | -0.001 | 1.120 | 0.902 | 1.006 | 0.772 | 2.640 | 1193
73 | 0.2477 | 0.074 | -0.001 | 1.149 | 0.891 | 1.047 | 0.790 | 2.879 | 1195
75 | 0.2284 | 0.100 | 0.000 | 1.158 | 0.915 | 1.061 | 0.806 | 2.116 | 1196
77 | 0.2195 | 0.148 | 0.001 | 1.305 | 1.054 | 1.178 | 0.941 | 0.248 | 1194
77 | 0.2195 | 0.105 | 0.000 | 1.230 | 0.993 | 1.130 | 0.880 | 1.180 | 1196
77 | 0.2195 | 0.063 | -0.004 | 1.224 | 1.000 | 1.131 | 0.886 | 1.738 | 1195
$\delta<0$ runs
69 | 0.2933 | 0.150 | -0.091 | 1.187 | 0.959 | 1.084 | 0.836 | 2.767 | 1195
69 | 0.2933 | 0.120 | -0.181 | 1.156 | 0.936 | 1.058 | 0.802 | 2.479 | 1195
69 | 0.2933 | 0.090 | -0.283 | 1.137 | 0.914 | 1.042 | 0.780 | 1.908 | 1195
69 | 0.2933 | 0.060 | -0.382 | 1.143 | 0.903 | 1.051 | 0.801 | 1.453 | 1195
69 | 0.2933 | 0.045 | -0.429 | 1.152 | 0.920 | 1.055 | 0.800 | 1.445 | 1195
69 | 0.2933 | 0.030 | -0.477 | 1.178 | 0.955 | 1.079 | 0.849 | 1.528 | 1195
71 | 0.2692 | 0.122 | -0.099 | 1.126 | 0.890 | 1.031 | 0.782 | 2.354 | 1195
71 | 0.2692 | 0.108 | -0.151 | 1.120 | 0.887 | 1.024 | 0.789 | 2.207 | 1195
71 | 0.2692 | 0.090 | -0.215 | 1.121 | 0.904 | 1.021 | 0.780 | 2.043 | 1195
$\delta>0$ runs
75 | 0.2284 | 0.112 | 0.051 | 1.150 | 0.891 | 1.059 | 0.809 | 2.174 | 1196
75 | 0.2284 | 0.127 | 0.098 | 1.141 | 0.902 | 1.050 | 0.797 | 2.313 | 1196
75 | 0.2284 | 0.139 | 0.147 | 1.137 | 0.898 | 1.044 | 0.799 | 2.121 | 1196
75 | 0.2284 | 0.150 | 0.195 | 1.146 | 0.898 | 1.053 | 0.807 | 2.823 | 1196
77 | 0.2195 | 0.136 | 0.237 | 1.173 | 0.924 | 1.078 | 0.810 | 1.993 | 1195
77 | 0.2195 | 0.150 | 0.300 | 1.167 | 0.927 | 1.065 | 0.809 | 2.273 | 1196
77 | 0.2195 | 0.161 | 0.351 | 1.165 | 0.924 | 1.065 | 0.812 | 2.294 | 1196
77 | 0.2195 | 0.172 | 0.398 | 1.164 | 0.922 | 1.066 | 0.809 | 2.446 | 1196
77 | 0.2195 | 0.190 | 0.484 | 1.181 | 0.932 | 1.080 | 0.835 | 3.327 | 1195
Table 2: Input/Output Halo Catalog for
$H_{0}=73,\delta=0,\Omega_{\mathrm{halos}}=0.5\Omega_{m}$ Scenario

Input | Output
---|---
1PGC | Halo | Name | SGL | SGB | $N_{g}$ | $o_{d_{o}}$ | $d_{o}$ | $\mu_{e}$ |  $\vrule depth=6.0pt,height=12.0pt,width=0.0ptcz_{o}$ | $M_{K_{s}}$ | $d_{m-o}$ | $cz_{m-o}$ | Log(M) | $V_{sgx}$ | $V_{sgy}$ | $V_{sgz}$ | $\chi$
|  |  | deg | deg |  |  | Mpc | mag |  $\vrule depth=6.0pt,height=12.0pt,width=0.0pt\mathrm{km~{}s^{-1}}$ | Mag | Mpc | km s-1 | $M_{\odot}$ | km s-1 | km s-1 | km s-1 |
44715 | 100001 | Coma Cl | 89.6226 | 8.1461 | 136 | 106 | 96.4 | 0.05 | 6926 | -30.199 | -0.0 | 1 | 15.246 | -162 | 183 | -217 | 0.03
57612 | 200002 | ESO137-006 | 188.2039 | 6.9076 | 167 | 0 | 65.3 | 0.00 | 4767 | -29.964 | 0.6 | -3 | 15.147 | -105 | 142 | -95 | 0.00
16124 | 0 |  | 2.6934 | -23.1887 | 6 | 0 | 93.9 | 0.00 | 6854 | -29.789 | 1.0 | -1 | 15.067 | -141 | -149 | 25 | 0.00
12429 | 200001 | Perseus Cl | 347.7159 | -14.0594 | 180 | 44 | 67.3 | 0.06 | 5546 | -29.756 | 2.5 | -30 | 15.079 | 205 | -115 | -47 | 1.43
70042 | 200450 | NGC7426 | 330.2136 | 33.2224 | 6 | 1 | 141.3 | 0.40 | 5812 | -29.725 | -69.1 | 0 | 14.357 | 189 | -118 | -51 | 0.00
166288 | 0 |  | 174.1307 | -9.8597 | 38 | 0 | 74.9 | 0.00 | 5465 | -29.289 | 4.8 | 1 | 14.889 | 89 | 194 | 144 | 0.00
36487 | 100005 | Leo Cl | 92.0255 | -10.4950 | 61 | 28 | 88.3 | 0.06 | 6458 | -29.189 | -0.3 | 1 | 14.778 | -76 | 408 | 193 | 0.12
63247 | 200228 | ESO460-004 | 216.2660 | 44.1999 | 11 | 2 | 137.4 | 0.35 | 7282 | -28.974 | -40.5 | 2 | 14.333 | -109 | -100 | -96 | 0.00
4149 | 200018 | Abell 2877 | 251.7068 | -13.6204 | 32 | 25 | 102.8 | 0.09 | 6876 | -28.871 | -10.0 | 24 | 14.533 | 114 | 10 | -48 | 2.49
9201 | 200013 | Abell 347 | 343.2559 | -5.2785 | 54 | 13 | 78.0 | 0.12 | 5940 | -28.801 | 2.1 | -5 | 14.628 | -189 | 64 | -483 | 0.48
2801981 | 0 |  | 184.6380 | 6.3997 | 17 | 0 | 67.1 | 0.00 | 4896 | -28.639 | 9.5 | -2 | 14.660 | 543 | 76 | -78 | 0.00
11188 | 200012 | Abell 400 | 312.1102 | -27.2311 | 44 | 33 | 92.9 | 0.07 | 7366 | -28.626 | 3.7 | -19 | 14.561 | 164 | -25 | 48 | 1.26
6962 | 200003 | Abell 262 | 334.9977 | -2.0249 | 63 | 37 | 66.7 | 0.06 | 5102 | -28.621 | -5.3 | 91 | 14.437 | 371 | -189 | -216 | 3.37
… data continues … |  |  |  |  |  |  |

Note. — Column descriptions: (1) PGC ID of brightest halo member; (2) Halo ID;
(3) Name of brightest halo member; (4) Supergalactic longitude; (5)
Supergalactic latitude; (6) Number of galaxy members brighter than
$K_{s}=11.75$ mag; (7) Number of distance measurements; (8) Observed distance
or Hubble distance if no measures; (9) Distance modulus error; (10) Observed
radial velocity in Local Group frame; (11) Output absolute $K_{s}$-band
magnitude ($M_{K_{s}}=3.27$ is assumed for the Sun); (12) Model minus observed
distance; (13) Model minus observed radial velocity; (14) Model Log10 of Mass;
(15,16,17) Supergalactic vx,vy,vz in frame of sample (18) Goodness measure.

Note. — Table 2 is published in its entirety in the machine-readable format. A
portion is shown here for guidance regarding its form and content.

Table 2 provides input and output data on each of the individual halos of the
fiducial $H_{0}=73,\Omega_{\mathrm{halos}}/\Omega_{m}=0.5$ calculation with
local density at the mean density. Only a few rows are presented here, but the
entire machine readable table is available at the online journals, and it has
been installed into the Extragalactic Distances Database
(https://edd.ifa.hawaii.edu).

Figure 7: Distribution of $\chi$ values for the $\mu_{e}\leq 0.3$ halos in the
fiducial scenario. Bin size is 0.2. Red curve is expected half-normal
distribution with 1194 measurements and assuming mean of 1. The means (black:
data; red: half-normal distribution) are solid vertical lines and medians are
dotted vertical lines. Figure 8: Individual $\chi$ values - for halos with
distance errors $\mu_{e}\leq 0.3$ versus observed redshift for the fiducial
scenario. Values below $\chi=0.1$ are set to that value to keep them in the
plot. Horizontal lines are placed at the mean and median values. Figure 9:
⟨$\chi$⟩ as a function of $\mathbf{H_{0}}$ and fraction of mass density in
halos. Goodness metrics, mean $\chi$ values, are shown next to plus marks for
3 runs at $H_{0}$= 69, 73, and 77 km s-1 Mpc-1 ($\Omega_{m}h^{3}=\xi_{P}$) and
mass-to-light ratios varied to obtain $\Omega_{\mathrm{halos}}/\Omega_{m}$ =
0.3, 0.5, and 0.7. Additionally, there are single runs at $H_{0}=71$ and 75.
Adding the IHM density brings the total density to the global mean. The
⟨$\chi$⟩ values are fit with a 2d second order polynomial, shown as contours.
The minimum $\langle\chi\rangle$ of the fit occurs at $H_{0}$=72.95. The
minima at each $H_{0}$ is shown as a dashed line.

Figure 10: Parameters along the minima $\chi$ - The top curve traces along the
minima ⟨$\chi$⟩ at each $H_{0}$ value of the fit to $\delta=0$ solutions
(previous figure). Bottom, the $\Omega_{\mathrm{halos}}/\Omega_{m}$ along the
trace is shown by solid black line. The red dash line provides the approximate
value of coefficient $\kappa$ in the luminosity to mass mapping (use axis on
the right).

Figure 9 presents the goodness metric, ⟨$\chi$⟩ for halos with $\mu_{e}\leq
0.3$, as functions of both $\Omega_{\mathrm{halos}}/\Omega_{m}$ and $H_{0}$.
The values for 11 model runs at $\delta=0$ are shown at the plus signs. Halos
with $\chi>7$ were omitted from the mean (usually just an object or two). We
also omitted halos with observed distances $d_{o}<3$ Mpc because they are too
sensitive to small changes in the model velocity of the Local Group. A second
order polynomial in two dimensions is fit to the data and is shown as
contours. The minima at each $H_{0}$ were found and they occur along the
dashed line. The overall minimum of the fit is very close to the fiducial
case. Minimum ⟨$\chi$⟩ versus $H_{0}$ is shown in the top section of Fig. 10.
The typical 1$\sigma$ error bar is shown at its minimum at $H_{0}=73.0$. In
addition, we have transformed the location of the minima into $\kappa$, the
coefficient in the mass-to-light ratio formula, seen as the red dashed line
and using the right hand axis. A value of
$\kappa=M/L_{K}=21.9~{}M_{\odot}/L_{\odot}$ is preferred at $H_{0}=73$. The
greater than linear increase of mass with $H_{0}$ is attributed to the
decrease in age of the universe and the increase in $\Omega_{\Lambda}$.

###  4.1 $\delta\neq 0$

Tully et al. (2016) demonstrated (their Figure 21) that, with the
Cosmicflows-3 distances used here, global outflow from the region within 8,000
km s-1 is anticipated for values of $H_{0}$ less than $\sim$74 km s-1 Mpc-1
and infall for larger $H_{0}$ values. Demonstrably, within a radius of 30 Mpc
we live in a local overdensity (see Figure 4). On the full scale of this
study, there have been suggestions that we live in an underdense region
(Keenan et al., 2013; Shanks et al., 2019; Böhringer et al., 2020).

Figure 11: Dependence of ⟨$\chi$⟩ on over/underdensity at select values of
$\mathbf{H_{0}}$ - A series of runs with different $\Omega_{\mathrm{IHM}}$ at
each $H_{0}$ value were fitted to find where the minimum ⟨$\chi$⟩ occurs as a
function of local (100 Mpc) overdensity in the models. The red lines are 3rd
order polynomial fits. Vertical dashed lines drop from the minima of each fit
and the values of $H_{0}$ are given next to this line. The horizontal dashed
lines are placed at the minima of the fits. The circle is placed at the
minimum of the solutions with $\delta=0$ in the previous 2 figures. The green
ticks at the bottom are analytical predictions from integration of the
Friedmann equation (§ 4.2).

Figure 12: Mass density components along valley of lowest ⟨$\chi$⟩, This
diagram delineates the valley where ⟨$\chi$⟩ is a minimum as the IHM density
is varied for each $H_{0}$ (see previous figure). (Top) The crosses reflect
the resulting fraction of density in halos, and the green line through them is
the fit, Eq. 9. The squares give the fraction of density attributed to the
IHM, and the dash-dotted blue line uses the fit, Eq. 10. The triangles show
the sum of the two mass components. The fiducial $H_{0}=73$ run with
$\delta=0$ is included and appears to lie along the ⟨$\chi$⟩ valley. Bottom
plots are the same as the top but densities are divided by the critical
densities rather than the global mean densities. The thin dotted horizontal
line in the top figure sits at the level of $\Omega_{\mathrm{IHM}}/\Omega_{m}$
expected from n-body simulations (discussed in § 4.3).

The ⟨$\chi$⟩ arising from assuming a range of $\Omega_{\mathrm{IHM}}$ values
in NAM runs, relaxing the requirement of no overall over/under density in the
study volume, are shown at $H_{0}$ values of 69, 71, 75 and 77 in Fig. 11.
These strings of values are fit with third-order polynomials. The vertical
dashed lines are placed at the fractional overdensity value where ⟨$\chi$⟩
reaches its minima. A circle is placed in Figure 11 at the minimum value of
the $\delta=0$ runs in Figure 9 that occurs at $H_{0}=73$, and the 1$\sigma$
errorbar is added.

It is to be noted that solutions with acceptable values of ⟨$\chi$⟩ can be
found from $H_{0}=69$ to 75. Increasing $\Omega_{\mathrm{IHM}}$ to the optimal
value at $H_{0}=75$ results in a ⟨$\chi$⟩ that is only about $1\sigma$ more
than our best case. The values of $\rho_{\mathrm{IHM}}/\rho_{m}$,
$\rho_{\mathrm{halos}}/\rho_{m}$, and $\rho_{\mathrm{local}}/\rho_{m}$ at the
minima in ⟨$\chi$⟩ are shown in Figure 12. For $H_{0}=69$, a very low value of
$\Omega_{\mathrm{IHM}}$ is required, implying a best fit underdensity of
$\sim$32%. These optimal scenarios, summarized in Figure 12, are described by
the relationship between $H_{0}$ and $\delta$:

$\frac{H_{0}}{H_{\delta=0}}=1+0.165\delta$ (7)

where $H_{\delta=0}$ = 73.0. This result can be compared with the related
problem of the fractional change in a galaxy’s expansion rate $H_{i}$ expected
in linear perturbation theory (Lahav et al., 1991; Peebles, 1993, p. 116),

$\frac{H_{0}-H_{i}}{H_{0}}=-\frac{v_{pec}}{H_{0}d}=\frac{\Omega_{m}^{0.55}}{3}\delta$
(8)

The linear theory coefficient on $\delta$ for $\Omega_{m}=0.248$ is 0.155.

The values of $\Omega_{\mathrm{halos}}$ can be thought of as an input
parameter in that the mass-to-light relation $\kappa$ was varied in accordance
with the fit in Fig 10. It is nearly constant and fit by:

$\Omega_{\mathrm{halos}}=0.122+3.9\times 10^{-4}(H_{0}-73.0).$ (9)

An expression fitting the run of $\Omega_{\mathrm{IHM}}$ is:

$\Omega_{\mathrm{IHM}}=0.612\left(\frac{H_{0}}{H_{limit}}-1\right)$ (10)

where $H_{limit}=61.2$ is the $H_{0}$ value below which the density of the IHM
at the bottom of the ⟨$\chi$⟩ valley becomes negative.

###  4.2 Top hat solutions to expansion rate discrepancy

Since the overall minimum ⟨$\chi$⟩ for $\delta=0$ cases occurs very close to
the fiducial $H_{0}=73$, we can expect that in each flat universe with
($H_{0}$, $\Omega_{m}$) pairings that the optimal density for the local region
will be close to one that obtains a local expansion of 73 km s-1 Mpc-1. Thus,
the density required to embed inside a flat universe of
$(h,\Omega_{m},\Omega_{\Lambda})=(h,\xi_{P}/h^{3},1-\xi_{P}/h^{3})$ a
homogeneous subregion with $\Omega_{m}=\Omega_{m}^{73}$ with expansion rate of
$h=0.73$ provides additional estimates of the parameters of the valley of
minimum ⟨$\chi$⟩ and a check on the NAM analysis. By setting the age and the
cosmological constant, $\Lambda\propto\Omega_{\Lambda}h^{2}$, of the embedded
scenarios to be the same as in the exterior, we can solve for the desired
subregion density $\Omega_{m}^{73}$ in:

$t_{0}(h,\Omega_{m},\Omega_{\Lambda})=t_{0}(0.73,\Omega_{m}^{73},\Omega_{\Lambda}h^{2}/0.73^{2}).$
(11)

The subregion scenario has no closed form formula for $t_{0}$, but it can be
found by integrating the solution for the Friedmann equation:

$H_{0}t_{0}=\int_{0}^{1}\frac{da}{\left[\Omega_{m}a^{-1}+\Omega_{\Lambda}a^{2}+(1-\Omega_{m}-\Omega_{\Lambda})\right]^{\frac{1}{2}}}$
(12)

To convert $\Omega_{m}^{73}$ to an overdensity, we note that mass density
varies $\propto\Omega_{m}h^{2}$, so the overdensity is:

$\delta=\frac{0.73^{2}\Omega_{m}^{73}h}{\xi_{P}}-1$ (13)

For $H_{0}=(67,69,71,75,77)$ one finds the estimates for $\delta$ are (-0.456,
-0.321, -0.169, 0.188, 0.396) as shown in Figure 11 as green ticks at the
bottom. For negative $\delta$, the NAM calculated best $\delta$s are in very
good agreement with these values, but at the higher values there is
substantial disagreement. This discrepancy from the top hat model is to be
expected, because adding density to regions in the non-linear regime will
alter velocities much more than in low density regions. This limits how much
IHM can be added before some regions badly overreact.

###  4.3 Matter outside of halos in simulations

N-body simulations can predict the amount of matter in and not in luminous
independent halos. We can check if our results are consistent with
simulations, and whether simulations favor a portion of the range of our
results. There are numerous papers on the mass function of halos found in
simulations (Sheth & Tormen, 1999; Jenkins et al., 2001; Warren et al., 2006;
Tinker et al., 2008; Watson et al., 2013; Diemer, 2020). The study of observed
halo masses by Kourkchi & Tully (2017) suggests that the luminous halos in our
catalog can be associated with almost all of the mass in halos above $2\times
10^{10}M_{\odot}$.

Figure 13: Integrated mass in halos above mass M - Integration of the mass
function of halos out to the splashback radius from n-body calculations
downward from the most massive divided by $\rho_{m}$. The divide between
luminous and non-luminous halos happens around $10^{10}M_{\odot}$. The IHM
fraction would be 1 minus this fraction.

A recent study by Diemer (2020) considers the mass within the second
turnaround (splashback) radius of halos in simulations, a choice that matches
the radius used in the definition of our halos. We chose to integrate their
mass spectrum, as found in the Colossus toolkit (Diemer, 2018), that includes
particles within the radius that includes 90% of the second turnaround points.
Integrating downward from the most massive halos (Figure 13), one finds 59% of
the total mass is in halos down to $2\times 10^{10}M_{\odot}$ for $h=0.7$.
Hence 41% of mass lies in smaller halos or at large in sheets, filaments, and
voids, a value most consistent with our $H_{0}\sim$72 result, as shown in
Figure 12. Establishing an accurate errorbar for this number is problematic at
this time.

###  4.4 Local Group velocity

The NAM models provide direct knowledge of the motion of the Local Group in
the frame of the center of mass of the whole sample, including the external
sources. For the $H_{0}=73$ model, the present $SGX,SGY,SGZ$ motion is
($-194$, 292, $-231$) km s-1 or 420 km s-1 into $SGL,SGB$ of (124°, $-34$°),
galactic $l,b$ of (264°, 39°). Subtracting this vector from the CMB dipole
leaves a residual velocity of 240 km s-1 into $SGL,SGB$ (165°, $-22$°), or
galactic $l,b$ of (293°, $12$°). It is not clear how much more motion comes
from large scale structure beyond our external catalog, or unrepresented voids
in the opposite direction (Hoffman et al., 2017) or from a more complete
accounting of the masses in the region covered by the external catalog.

###  4.5 Influence of external catalog

If we remove the 44 external masses in our external catalog and hold all paths
from the $H_{0}=73$ solution fixed except for the path of the Local Group, the
direct effect of the external catalog on the Local Group can be seen. The
Local Group $SGX,SGY,SGZ$ motion is $V_{LG}$ = ($-161$, 290, $-234$) km s-1 or
406 km s-1 into $SGL,SGB$ of (119°, $-35$°), galactic $l,b$ of (259°, 40°).
Subtracting the Local Group motion with the external masses from the result of
the previous section gives a Local Group motion arising from external masses
into $SGX,SGY,SGZ$ of $V_{LG}^{ext}$ = ($33$, -2, -3) or 33 km s-1 into
$SGL,SGB$ = (356°, -5°) or galactic $l,b$ = (141°, -4°), a direction into
$SGX$ and the zone of avoidance. Its low amplitude means that, apparently, the
objects from around the sky between 8,000 and $16,000$ km s-1 mostly cancel
one another and result in a very small net motion.

###  4.6 Virgo Cluster

The nearby and well studied Virgo Cluster has a relatively poor $\chi=2.15$ in
the fiducial $H_{0}=73$ model (see Column 9 of Table 1). The observed average
redshift for the cluster is 1058 km s-1 in the Local Group frame (1150 km s-1
heliocentric) and its distance with the assumed zero point is 15.8 Mpc,
resulting in a peculiar velocity of $-94$ km s-1. The fiducial model gives a
peculiar velocity of $-111$ km s-1, but to get so close in cz, Virgo was moved
out to 16.3 Mpc, about 1$\sigma$ off. When placed at its observed distance,
the peculiar velocity is $-241$ km s-1. Models with a lower mass get better
$\chi^{2}$, but the mass in the fiducial model already strains the lower mass
limit. The model mass of the cluster is $3.2\times
10^{14}(\frac{d}{15.8\mathrm{Mpc}})M_{\odot}$ while the virial theorem mass
estimate is $6.0\pm 0.9\times 10^{14}(\frac{d}{15.8\mathrm{Mpc}})M_{\odot}$
(Kashibadze et al., 2020). The $H_{0}=77$ cases reveal that higher $H_{0}$ can
improves both cz and mass, but only slightly.

The peculiar velocity of Virgo is predominantly the component of the velocity
of the Local Group in the Virgo direction because the radial component of the
cluster’s velocity is quite low in all cases. Therefore, issues related to the
path of cluster are probably not an important source of error. One might
suspect problems with the specific choice of path for the Local Group
depending on tugs from nearer masses, but a look at our neighbors show all
moving towards Virgo at nearly the same velocity. A possibility is that the
IHM density is biased (Bardeen et al., 1986) to unusually low densities in
this region rich in halos. The overdensity to the Virgo Cluster of $\sim$2.0
in our fiducial model, could be as low as $1.5$ without an IHM and that would
reduce the peculiar velocity by 25%.

On the observational side, the $\chi$ could be improved with a lower observed
velocity (unlikely) and/or a greater observed distance (more plausible). The
assumed statistical distance uncertainty is 3% or 0.5 Mpc. The depth of the
cluster from fore to backside of the second turnaround surfaces (roughly the
virial dimensions) is 4 Mpc, 8 times the distance uncertainty of the center.
The Virgo Cluster has structural complexity (Mei et al., 2007; Ferrarese et
al., 2012). Many galaxies are falling in today from the Southern Extension
(Tully & Shaya, 1984; Kourkchi & Tully, 2017). M49, the brightest galaxy in
the cluster, is near the geometric edge with its own entourage. The center of
mass of the cluster, with an extensive and still unsettled halo, may not be
precisely at the mean position of its galaxies. If the simplistic description
of Virgo as a single entity is relaxed, solutions could be found with reduced
$\chi$ but great ambiguity.

###  4.7 Centaurus Cluster

Figure 14: Possible Paths for Cen30 and Cen40 Clusters \- The Centaurus
Cluster is considered to be two clusters superimposed along the line of site.
Trajectories in comoving supergalactic coordinates have filled circles at 0.5
Gyr intervals and large open circles show the position at the present day.

Another cluster that defies a simple description is the Centaurus Cluster.
This cluster has been suspected to be two superposed clusters which have been
dubbed Cen30 and Cen45 because one component peaks in velocity at $\sim$3,000
km s-1 and the other peaks at $\sim$4,500 km s-1 (Lucey et al., 1986; Stein et
al., 1997). At first, we retained it as a single cluster in our runs at
$d_{o}=37\pm 1.1$ Mpc and with $cz$=3,142 km s-1. The model output distances
were all around 43 Mpc, which is unacceptable. We therefore broke it up,
guided by its double-peaked velocity distribution, into the two suspected
components. Based on 43 mainly fundamental plane and spiral luminosity-
linewidth distances, Cen30 was placed at 37.8 Mpc with $cz$ = 2,825 km s-1
while with 23 distances Cen45 was placed at 45.5 Mpc with $cz$=4,306 km s-1,
where Cen30 is about twice as massive as Cen45. The separation between the two
components based on all Cosmicflows-3 measures is $8\pm 2$ Mpc. It is to be
noted, though, that surface brightness fluctuation measures found the two
components to be at similar distances within the uncertainties Mieske & Hilker
(2003); Mieske et al. (2005). A special search of solutions was made by
repeatedly ‘backtracking’ their paths from the present positions and redshifts
to find a possible strong interaction, but the separation of about 8 Mpc is
too large. Interacting scenarios with separations of about 2 to 3 Mpc, at
best, were found that provide enough positive peculiar velocity to Cen30.

Figure 14 shows such an example of a close interaction in which Cen45 flies by
Cen30 and then goes out to $\sim$2 Mpc beyond it. Cases with greater
separation were found, but they required Cen45 having an earlier interaction
with yet another halo. Although possible, we have refrained from using such
complicated scenarios. Due to the large uncertainty of such models related to
the ambiguity in the relative distances of Cen30 and Cen45, we have removed
both of these halos from all statistics in this study.

##  5 Flow Patterns within 100 Mpc

Our discussion takes frequent recourse to the accompanying group of four
videos444The four videos are narrated by Prof. R. Brent Tully. and four
interactive models; see animated Figures 15, 16, 17, 18 and interactive
Figures 19, 20, 23, 28555The videos and interactive visualizations can also be
accessed at the following address: http://irfu.cea.fr/nam8k. In the opening
overview video (Figure 15) and interactive model (Figure 19), the first
thought of the view of 10,000 orbits is likely to be of a bad hair day. A
noticeable gap runs horizontally. This poorly represented part of the plot is
associated with the zone of obscuration of the Milky Way and will be remarked
frequently depending on the orientation of views. Elements of the model do lie
in this zone sufficient to maintain roughly constant density, as described in
§2, but they are numerically few. The eye also sees a great concentration of
orbits near the center. Near to home, within 10 Mpc, there are an abundance of
(mostly light weight) galaxies with accurately known distances from the tip of
the red giant branch method (Anand et al., 2021).

Figure 15: Overview video - The high-resolution video in UHD-4K format is
available online for viewing and download at
https://vimeo.com/pomarede/nam8k-overview

Figure 16: Video visualization of the orbits in Laniakea. The high-resolution
video in UHD-4K format is available online at
https://vimeo.com/pomarede/nam8k-laniakea

Figure 17: Video visualization of the orbits in Perseus-Pisces. The high-
resolution video in UHD-4K format is available online at
https://vimeo.com/pomarede/nam8k-perseus-pisces

Figure 18: Video of the Future. The high-resolution video in UHD-4K format is
available online at https://vimeo.com/pomarede/nam8k-intothefuture Figure 19:
Overview visualization displaying the entire sample of 9,719 model objects and
44 externals. Start interaction at https://bit.ly/NAM8Koverview. The red,
green, blue 1,000 km/s-long arrows emanating from the origin indicate the
three cardinal axes of the Supergalactic coordinate system.

At intervals throughout the videos, one sees the superposition of density
contours on the orbits. These surfaces are derived from a quasi-linear Wiener
filter analysis of Cosmicflows-3. The methodology is described in an earlier
analysis of Cosmicflows-2 data (Hoffman et al., 2018). These contours provide
a reference frame to aid in orientation, appreciation of scale, and to show
target destinations for flow patterns. The quasi-linear maps are chosen over
other options because of the excellent resolution that they afford. The quasi-
linear analysis with Cosmicflows-3 data will be described in detail by Hoffman
et al. (in preparation).

In the big picture, it is immediately evident that there are three overarching
flow pattern zones within 8,000 km s-1, which will be referred to as Laniakea,
Perseus-Pisces, and Great Wall, as seen in the interactive Figure 20. Each
will be discussed in turn. The first two are reasonably contained within the
boundary of the study. Only the closest parts of the Great Wall are being seen
so the current description of flows in the region of that structure is
perfunctory. Another dominant overview feature is the relative absence of mass
elements in the sector at negative $SGY$ (the sky south of the Galactic plane)
and negative $SGX$. See in particular the interactive model restricted to $\pm
20$ Mpc about $SGX=0$ (Interactive Figure 20).

Figure 20: Visualization of the Supergalactic Equator plane, where the final
positions of objects is restricted to the slice $-20<SGZ<+20$ Mpc. Start
interaction at https://bit.ly/NAM8KFigequator.

###  5.1 Laniakea

The concept of the Laniakea Supercluster was introduced by Tully et al.
(2014). The Milky Way lies within its embrace. The major over density
constituents, strongly linked together, are the classical Local Supercluster
(de Vaucouleurs, 1953), the so-called Great Attractor (Dressler et al.,
1987a), and the Pavo$-$Indus filament (di Nella et al., 1996). In all, there
are 18 Abell clusters within the envelope of Laniakea if one counts the
supplemental entries in the list of Abell et al. (1989) and the Virgo Cluster.
As must be expected, the overdense structures are encased in voids that take
up most of the volume of Laniakea.

The pattern of orbits in the region of the Local Supercluster were discussed
extensively by Shaya et al. (2017) and the main points are summarized here.
The dominant patterns within this nearby region are flows out of the Local
Void and toward the Virgo Cluster. The entire region is participating in a
flow toward the Great Attractor. Specifically, (1) essentially all galaxies
north of the supergalactic equator (the sector of the Local Void), those with
positive values of $SGZ$ in supergalactic coordinates, are moving toward
negative $SGZ$, whereas galaxies south of the supergalactic equator have
disorganized $SGZ$ motions. Also, (2) essentially all galaxies across the
region have motions toward negative $SGX$ and positive $SGY$, the direction
toward the Great Attractor and the apex of the CMB dipole.

Figure 21: Orthogonal views of orbits in the Great Attractor region. At left
is a view from $-SGX$ and at right is a view from $-$SGY. The orbit of the
Local Group is in red, terminating at the origin marked by orientation arrows:
red toward $+SGX$, green toward +SGY and blue toward +SGZ. Throughout this
presentation, the orientation arrows are 20 Mpc in length. Isodensity contours
of the quasi-linear model are shown faintly. Figure 22: Orbits in the Great
Attractor region - viewed from the north supergalactic pole (+SGZ).

We retain the nomenclature “Great Attractor” but what does this phrase mean?
We take it to refer to the rough terminus of flow patterns within the Laniakea
envelope. It is not a well defined location.

An overview of the Great Attractor region is shown in the two panels of Figure
21. The scenes are approximate extracts from the Part-1 (Overview, see Figure
15) video at 5:50 and the Part-2 (Laniakea, see Figure 16) video at 3:10. The
overwhelmingly dominant pattern is the downward flow of orbits toward negative
$SGZ$ for galaxies above the supergalactic equator at $SGZ=0$. This pattern is
the same as was seen in the Local Supercluster study, now exhibited on a much
larger scale. The space at positive $SGZ$ is the domain of the Local Void
blending into the Hercules Void (Tully et al., 2019). There is also a trend of
motions upward from galaxies at negative $SGZ$ toward $SGZ=0$ but with less
coherence across the full domain of these plots. Closer attention is required
to isolate patterns.

A prominent feature of local peculiar motions is a flow toward $-SGX$ in the
Local Group vicinity, and it is seen in the right panel of Figure 21 and in
Figure 22 (see the Part-2 video at 2:45, Figure 16) that this pattern persists
from $SGX=+10$ Mpc to $SGX=-45$ Mpc before breaking down. Detailed attention
suggests a convergence of flows roughly centered between the three clusters
A3565, A3574, and S753 at roughly $SGX=-55$ Mpc, as seen in the interactive
Figure 23. These entities are part of what was called the 4-cluster strand by
Courtois et al. (2013). A3565 and A3574 are each $\sim$8 Mpc from the more
central S753, separations comparable to each of the individual turnaround
radii about these clusters. It is plausible that these three clusters are
destined to merge although this scenario is not favored in our discussion in
§6. While there is a convergence of flow lines in their proximity, it is a bit
of a mystery that this location is not particularly prominent in the quasi-
linear density maps.

Figure 23: Zoom visualization of the Great Attractor region. Start 4D
interaction at https://bit.ly/NAM8KGA. In this visualization of the n-body
solution, the time evolution is run on the basis of one second for one Gyr.

Beyond the 4-clusters area flows are confused. The convergence of flows toward
the supergalactic equator is retained: toward negative $SGZ$ from above the
plane and toward positive $SGZ$ from below the plane. However there is
disorganization in the lateral $SGX$ motions. We are approaching the 8,000 km
s-1 edge of the orbit reconstruction volume here. External tides are only
crudely represented, and not far beyond in this sector is the Shapley
Concentration, the densest collection of rich clusters within $0.1c$ (Tully et
al., 1992). Evidently, the fringes of the influence of Laniakea Supercluster
are being probed. Deeper surveys are required for a reliable picture.

Figure 24: Coherent flow pattern through the region of the Antlia Cluster to
the region of the Centaurus Cluster. This Antlia strand is a part of the more
extensive Centaurus-Puppis-PP filament.

Returning to the core region of the Great Attractor, the most prominent
feature off the supergalactic equator is a structure running from the
Centaurus Cluster to the Antlia Cluster at negative $SGZ$, slanting to
positive $SGX$ but at almost constant $SGY$ (see Figure 24). This structure
was called the Hydra Wall by Fairall (1998) but the Hydra Cluster is to the
background. It is called the Antlia strand by Courtois et al. (2013); one of
five filamentary structures that converge on the Centaurus Cluster (two of the
others have already been given attention; the 4-clusters strand and the planar
Local Supercluster that runs between the Virgo and Centaurus clusters). The
flow pattern running up from Antlia Cluster through Centaurus Cluster is
coherent and pronounced. As seen in Figure 24 (the Part-3 Perseus-Pisces video
at 5:00, Figure 17), if followed beyond Antlia toward positive $SGX$ there is
a structure that extends all the way over to the Perseus-Pisces complex. This
feature will be given attention at a further point in the discussion.

Figure 25: The Pavo-Indus filament - a prominent structure that rises above
the supergalactic equator to embrace the Local Void.

The really major appendage off the Great Attractor core of Laniakea is the
Pavo-Indus filament running from the Norma Cluster to high supergalactic
latitudes. Fairall (1998) called this structure the Centaurus Wall. Some
ambiguity arises regarding its connection with the Great Attractor region
because of the intervening zone of obscuration (Kraan-Korteweg & Lahav, 2000;
Radburn-Smith et al., 2006). The suggestion that Pavo-Indus is connected to
the Perseus-Pisces complex by di Nella et al. (1996) has more recently been
confirmed by Pomarède et al. (2017) through a feature called the Arch. Pavo-
Indus is a cupped sheet-like structure, a curved partial wall enclosing the
Local Void (Tully et al., 2019). Flow patterns illustrated in Figure 25
(Part-2 Laniakea video at 5:00, Figure 16) are very coherent, rushing toward
the surface of the sheet from both sides and downward toward Norma Cluster and
the Great Attractor core.

###  5.2 Perseus-Pisces

Our vantage point within Laniakea Supercluster is very near the boundary with
the Perseus-Pisces complex. The full extent of this entity remains to be
defined because of the limitations in distance of current data.666We resist
calling the Perseus-Pisces complex a supercluster because we may not be
capturing the full extent of this structure. Following from the discussion of
Laniakea (Tully et al., 2014) we reserve the nomenclature of “supercluster” to
refer to the full gravitational basin of flow lines jointly decoupled from the
cosmic expansion. The full extent of the Perseus-Pisces structure remains to
be defined. However the flow toward the spine of the Perseus-Pisces filament
is dramatic, as seen at 6:32 of the Part-1 Overview video (Figure 15, see also
Figure 26). There is a coherent fore-side flow away from us toward positive
$SGX$ out of the Local Void. There is evident backside flow in our direction
toward negative $SGX$ although with less fidelity given the proximity to the
sample boundary. The Perseus Cluster is a pronounced focus of flows as can be
seen at either 7:00 in the Part-1 Overview video (Figure 15) or at 1:03 in the
Part-3 Perseus-Pisces video (Figure 17, see also Figure 27).

Figure 26: Flows toward the spine of the Perseus-Pisces filament - viewed from
the supergalactic north pole (+SGZ). The Local Group lies at the origin marked
by the red-green-blue orientation arrows. Figure 27: An orthogonal view of the
Perseus-Pisces region looking in from $\mathbf{-SGX}$. The Perseus Cluster is
a focus of flow patterns. The A400 Cluster is a secondary focus of flows. The
NGC 1600 and Cen-Vir-PP features are parts of bridges to Laniakea
Supercluster.

The Abell 400 Cluster is seen as a secondary focus of flows. Not too much
should be said about this feature because it is in proximity to the sample
boundary. Pomarède et al. (2017) have referred to structure in this region as
the Funnel. A case can be made that it is part of the extremely extensive
South Pole Wall (Pomarède et al., 2020).

The features labeled NGC 1600 and Cen-Vir-PP are parts of filaments running
toward us and connecting to our Laniakea Supercluster. These filamentary links
will be given attention later in the discussion.

###  5.3 Great Wall

The Great Wall (de Lapparent et al., 1986) is roughly confined to an $SGX-SGZ$
plane at constant $SGY$ that intersects our 8,000 km s-1 volume at its nearest
part in the vicinity of the Coma and A1367 clusters. Very important parts,
such as its extension to the major Hercules Cluster, are external to our study
region. There are coherent flows from the foreground toward positive $SGY$ and
the Great Wall (8:16 in the Part-1 Overview video Figure 15) arising out of
the Hercules Void, the separator between Laniakea and the Great Wall (Tully et
al., 2019). Whatever motions arise from the background are lost to this study.

###  5.4 Supercluster Connections

The most impressive bridge between the Laniakea and Perseus-Pisces
gravitational basins was called the Centaurus-Puppis-PP filament by Pomarède
et al. (2017). See at 5:00 in the Part-3 Perseus-Pisces video (Figure 17). At
the Laniakea end, near the Great Attractor core, this filament manifests as
the Antlia strand that was discussed in §5.1. At the Perseus-Pisces end, the
filament is picked up at the NGC 1600 knot mentioned in §5.2 and seen in
Figure 27. In between, this filament passes close enough for the very coherent
flow pattern to have been mapped in detail in the numerical action orbit study
of Shaya et al. (2017). As has become familiar with such bridges between
potential wells, there is a point along the filament where there is a flow
reversal (Pomarède et al., 2017). In this case, NGC 1600 is within the embrace
of Perseus-Pisces, with the reversal occurring somewhat toward negative $SGX$
of the group containing this galaxy.

A second very evident bridge filament was called Centaurus-Virgo-PP by
Pomarède et al. (2017) because the band goes from Centaurus Cluster through
the Virgo Cluster then onward to the Perseus Cluster. Like the previously
discussed Cen-Pup-PP filament, it is not well known because parts are hidden
by the zone of obscuration. As the filament enters within the dominant
influence of the Perseus-Pisces gravitational domain the flow pattern orients
toward the Perseus Cluster as seen in Figure 27 and the Part-3 Perseus-Pisces
video segment at 1:03 (see animated Figure 17).

The Part-3 Perseus-Pisces video at 6:45 also gives attention to a less well
defined structure (in parts, because of the considerable distance and its
concealment by obscuration) that Pomarède et al. (2017) called the Centaurus-
Arch-PP filament. This feature is the high $SGZ$ cap on the Local Void (Tully
et al., 2014, 2019) linking the Pavo-Indus and Perseus-Pisces filaments.

Courtois et al. (2013) also identified a filament emanating from the Centaurus
Cluster that splits into two branches. One passes through the Fornax Cluster
and Eridanus cloud on its way to the Perseus-Pisces complex and is labelled
For-PP at 5:00 of the Part-3 Perseus-Pisces video (animated Figure 17). The
other passes through NGC 6868 and the Telescopium-Grus cloud (Tully & Fisher,
1987) on its way to the Perseus-Pisces region.

The filaments discussed in this section can all be seen as fragments of an
envelope bounding the Local Void (Tully et al., 2019). The flow patterns out
of the Local Void are nicely seen at 6:45 and 7:32 of the Part-3 Perseus-
Pisces video (see animated Figure 17).

In the final scenes of the Part-1 Overview video (see animated Figure 15), the
envelope of Laniakea Supercluster is superimposed as defined by Tully et al.
(2014). Aside from the Arrowhead “mini-supercluster” (Pomarède et al., 2015)
that is too inconsequential to be well proscribed by this orbit analysis,
there are two neighboring structures that to a greater or lesser degree fall
within the 8,000 km s-1 confine of this study: Perseus-Pisces and Great Wall.
It is seen that the putative boundary of Laniakea actually extends beyond our
study region. That outer area is mostly the domain of the Sculptor Void (Tully
et al., 2019). At the interfaces between Laniakea and alternatively Perseus-
Pisces and the Great Wall, the boundaries of separable gravitational wells lie
in voids, albeit threaded with shearing filaments.

##  6 Next Hubble time

The orbits in our model can be computed into the future. There is an
asymmetry. Information on mergers in the past is lost, with the consequences
swallowed into the elements of our model constrained at the present. By
contrast, we can predict the occurrences of mergers in the future. Any
prediction is only as good as the representation of our orbits, so caution is
warranted in specific cases. Moreover, the elements are not the point-like
billiard balls implicit to our model. The elements are extended and composed
of many particles, such that they will tend to coalesce through dynamical
friction as their envelops overlap. In our Part-4 Into the Future video (see
animated Figure 18) and two interactive displays going into the future (see
interactive Figures 23 and 28), elements pass through each other or dance in
ways not to be believed.

The overarching picture for the future is boring. Except in a small fraction
of places, the velocities of attraction toward high density locations are
overwhelmed by the cosmic expansion, made worse by $\Lambda$ as the expansion
rate $H(t)$ never drops below $H_{0}\Omega_{0,\Lambda}^{1/2}$. This basic
truth is obscured by presentations in co-moving coordinates. (See the Part-4
video in Figure 18 at 0:53 to see the model evolve in physical coordinates.)
Most merging occurred in the distant past when entities were close to each
other in physical space. The onset of the influence of dark energy only
accentuated the drift to boredom.

Nevertheless, there will be some merging in the future. We continued the
computed paths in our $H_{0}=73$ km s-1 model until $a=2$, when the age of the
universe will be 23.66 Gyr and $H=64$. It is interesting to consider a few
cases where the events will be particularly substantial and the evidence is
reasonably strong. We will give attention to what is projected to happen along
the Perseus-Pisces chain of clusters (see interactive Figure 28), the core of
the Great Attractor (see interactive Figure 23), and our immediate region from
the Local Group to the Virgo Cluster.

Our numerical action orbit reconstruction emphatically confirms that the
Perseus Cluster (Abell 426, Halo 200001 (Tully, 2015b)) is the dominant
attractor within the Perseus-Pisces filament, as seen in the interactive
Figure 28. The convergence of orbits toward this cluster, already apparent in
orbits from the past, becomes more explicit with orbits into the future that
pass through or loop around Perseus (see 1:50 and following in the Part-4 Into
the Future video in the animated Figure 18). There is a secondary
gravitational basin along the filament involving the Pisces Cluster (NGC 410
Cluster: Halo 200005). Evidently, there will be a collision between Pisces and
the comparably large NGC 507 Cluster (Halo 200006) when the universe has an
age of 22 Gyr (3:10 in the Part-4 Into-the-Future video Figure 18). The rich
Abell 262 (Halo 200003) and NGC 315 (Halo 200070) clusters are in the
neighborhood.

Figure 28: Zoom visualization of the Perseus-Pisces region. Start 4D
interaction at https://bit.ly/NAM8KPP. In this visualization of the n-body
solution, the time evolution is run on the basis of one second for one Gyr.

In the densest part of the Great Attractor region, there are two foci of
convergence. The nearer is around the Centaurus Cluster (Abell 3526). In
recognition of the two kinematic components to the traditional Centaurus
Cluster, our modeling considers that there are two entities converging upon
each other roughly today (Lucey et al., 1986; Stein et al., 1997). The larger
portion with dominant galaxy NGC 4696 (1PGC 43296, Halo 100003)7771PGC is the
Principal Galaxies Catalog name (Paturel et al., 1996) of the dominant galaxy
in a group. Halo names are from Tully (2015b). has observed velocities
$\sim$3,100 km s-1. The smaller portion about dominant galaxy NGC 4709 (1PGC
43423, Halo 120003) has observed velocities $\sim$4,600 km s-1. The orbits of
these two elements are merging by construction. The pass through projected for
the future will surely be damped by dynamical friction. Meanwhile, other
elements will have merged or will be on the verge of merging with Centaurus by
twice the current age of the universe (Part-4 Into-the-Future video at 4:21 in
the animated Figure 18).

The other focus of convergence (seen in the same video sequence) is the
vicinity referred to by Courtois et al. (2013) as the ‘4 clusters strand’.
Along the strand lie Abell 3565 (1PGC 48040, Halo 100027), Abell 3574 (1PGC
49025, Halo 100011), and Abell S753 (1PGC 50100, Halo 100012). Abell 3537 lies
apart, but the cluster around NGC 5135 (1PGC 46974, Halo 100074) is a good
fourth element. Over a doubling of the age of the universe, elements merge
with these major clusters and the four draw closer in co-moving coordinates
but they move apart in physical coordinates. A more refined model at the
considerable distances of these clusters is needed to know their fate.

Closer to home, the Virgo Cluster is easier to study (Shaya et al., 2017).
Many distances are known with high accuracy and the associated galaxies
consequently have well determined peculiar velocities. Tip of the red giant
branch contributions with 5% accuracies to about 500 galaxies enrich the
analysis. Tidal influences are well described on relevant scales.

Infall and accretion into Virgo are seen following 5:45 in the Part-4 Into-
the-Future video (Figure 18)888The slingshot ejection of 3 elements including
M61 and NGC 4636 at $\sim$17 Gyr are artifacts of the point-like nature of the
model constituents. It is expected that these elements will be consumed by the
cluster.. We see a great deal of activity because the coverage is dense; the
situation must be similar about other rich clusters at greater distance. The
focus of the current study is on large scales, so we defer to the more local
numerical action model by Shaya et al. (2017) for Virgo infall details. To
summarize from that study, the mass within the collapsed cluster was found to
be $6.3\pm 0.8\times 10^{14}~{}M_{\odot}$ and the surface separating infall
from expansion at a=1 (encompassing the galaxies that reach the cluster by
$a=2$) is roughly triaxial, with radial dimensions 6.6 Mpc along the axis of
current major infall (through the Virgo Southern Extension), 7.9 Mpc toward
the relative voids at the Supergalactic poles, and an intermediate 7.2 Mpc on
the axis aligned with our sight-line. The mass in groups within turnaround
including the Virgo Cluster is $8.3\times 10^{14}~{}M_{\odot}$.

The video pauses at a=1 to illustrate the view from our position today. The
concentration of pink elements that represent major members of the Virgo
Southern Extension are found to be a transient occurrence today of elements
from diverse origins in our model.

The video leads us back home at 7:08. We live in a boring place, as is common
in the universe today. The merger between the Milky Way and Andromeda
anticipated in 6 Gyr (van der Marel et al., 2012) is concealed here by the
union of the two in the Local Group entity in our model. The most exciting
local event predicted is the merger between M83 and Centaurus A, combining
into the most massive halo within 8 Mpc. We are trending toward this
collapsing region in co-moving coordinates but will never be captured.

##  7 Summary

Most peculiar velocity studies invoke linear perturbation theory.
Cosmicflows-3 distance and velocity information (Tully et al., 2016) and 2MRS,
the 2MASS redshift survey to $K_{s}=11.75$ (Huchra et al., 2012), have been
used in linear regime studies (Graziani et al., 2019; Lilow & Nusser, 2021;
Hoffman et al., 2021) and there have been studies with variant samples
(Carrick et al., 2015; Jasche & Lavaux, 2019; Kitaura et al., 2021).

However in regions of higher density, like the Great Attractor core of
Laniakea Supercluster or the Perseus-Pisces complex, and in large voids, the
linear approximation is inappropriate. Orbit reconstructions from the
Numerical Action Method (Peebles, 1989; Peebles et al., 2001; Shaya et al.,
1995, 2017) provide a granular description of peculiar velocities up to the
point of shell crossing. Upon collapse, orbits become too confused to be
recovered, so the entities whose orbits we follow are the virial / second
turnaround domain of halos (be these called ‘groups’ or ‘individual
galaxies’).

The 9,719 entities drawn from the 2MRS $K_{s}=11.75$ catalog provides a good
description of the distribution of old stars, expected to correlate with mass,
over the unobscured 91% of the volume within the 8,000 km s-1 limit of this
study. The grouping into turnaround domains was discussed by Tully (2015a, b),
where justification is given for the linkage between observed light and
underlying mass described by Eq. 4 and found adequate in the more local
numerical action study by Shaya et al. (2017). Somewhat lower halo $M/L_{K}$
values are preferred from the present analysis, as illustrated graphically as
a function of $H_{0}$ in Fig. 10.

Mass in the 9% of the sky behind the Milky Way is poorly established. We
assign mass to known entities sufficient to maintain average density in the
zone of obscuration but the orbits of these and neighboring objects should not
be trusted. Likewise, our description of external tidal influences with 44
major clusters at $8,000-16,000$ km s-1 is crude. Were we to repeat this
exercise, we might use the description of matter in this velocity regime given
by Lilow & Nusser (2021).

The other essential ingredients for this analysis are the distance measures
and inferred peculiar velocities provided by Cosmicflows-3 (Tully et al.,
2016). Although individual errors are often large, averaging within groups
provides a corpus with small statistical uncertainties.

Orbits of our $\sim$10,000 constituents are constrained by the numerical
action methodology (Peebles, 1989) to give end-point velocities and distances
consistent with observations. However, we know both from observations of the
Ly$\alpha$ forest and from numerical simulations that there must be
significant mass outside the collapsed entities whose orbits we follow. Our
simplistic method to account for this component is to provide a uniform
density described by the parameter $\Omega_{\mathrm{IHM}}$. The total local
density is the sum of the contributions in halos and in the interhalo medium:
$\Omega_{\mathrm{local}}=\Omega_{\mathrm{halos}}+\Omega_{\mathrm{IHM}}$.

We find optimal models, assuming the Planck Satellite acoustic angular scale
constraint: $\Omega_{m}h^{3}=\xi_{P}=0.096345$ (Planck Collaboration et al.,
2020). We find a good fiducial model with
$\Omega_{\mathrm{local}}\simeq\Omega_{m}$, with $H_{0}\simeq 73.0$ km s-1
Mpc-1 and $\Omega_{m}\simeq 0.248$, with half the mass in the observed bound
mass of galaxies, groups, and clusters. The other half would reside in a non-
luminous medium lurking beyond these entities.

Good models with other values of $H_{0}$ can be found if the constraint of
$\delta=0$ is relaxed. We find, both by numerical models of the 8,000 km s-1
region with alternative values of $\Omega_{\mathrm{IHM}}$, and by analytical
models, how $\delta$ varies with global $H_{0}$ along the best fit valley in
⟨$\chi$⟩. Models with acceptable $\chi$ values can be found with $H_{0}$
almost as low as the range of Planck and other CMB determinations, but these
cases require that we live in an unusual underdensity ($\sim$30% to reach as
low as $H_{0}=69$) over a volume of at least 100 Mpc in radius.

From observations, on the scale of 30 Mpc, we live in the overdense region of
the Local Supercluster and its appendage to the Great Attractor complex. On
larger scales of $\sim$100 Mpc, some $K$-band and x-ray based studies suggest
(Keenan et al., 2013; Böhringer et al., 2020) that we may live in a
substantial hole. These studies are difficult to interpret, though, because
the concentration of stars and, particularly, thermal x-ray emission arising
from massive clusters, are likely biased tracers of mass (Kaiser, 1984).

The underdensities of the most immediate Local and Hercules voids (Tully et
al., 2019) are reasonably offset by the bounding Perseus-Pisces and Great Wall
structures. However, a very large fraction of the volume south of the Milky
Way is the domain of the Sculptor Void (Kauffmann & Fairall, 1991; Tully et
al., 2019), leaving room for doubt as to whether our study region is
underdense. Accordingly, there could be a preference for a lower value of
$H_{0}$ that accommodates an outflow. Considerations derived from simulations,
discussed in §4.3, favor the value of $H_{0}\sim$72 km s-1 because of the
expectation of the fraction of mass outside of the massive halos of about 41%.
However there must be reservations in a comparison of our model of discrete
halos on a smooth IHM pedestal with the complexity of the real universe. Other
external information regarding the Hubble Constant obtained from the far field
from Type Ia supernovae, consistent with our zero-point calibration (Neill et
al., 2014; Riess et al., 2021), place limits on $H_{0}$ in the $\sim$73 range.
Using a sample of 1300 SNIa, Kenworthy et al. (2019) find no evidence of a
change in Hubble Constant values across the redshift range $0.02<z<0.15$.

Turning to the orbits found in our best models, the most evident features are
the flow patterns that formed the three main concentrations of galaxies: the
core of Laniakea, the Perseus-Pisces spine, and the nearest sector of the
Great Wall. Shears in the flows occur in the voids: in the Local Void between
Laniakea and Perseus-Pisces, and in the Hercules Void between Laniakea and the
Great Wall. The largest of the nearby voids, the Sculptor Void, extends well
beyond our study region, creating a division between Laniakea and the South
Pole Wall (Pomarède et al., 2020).

Drilling down in detail, we document the dominant flow downward at all
positive $SGZ$ over the Laniakea region, including the substantial coherent
flow along the Pavo-Indus component of Laniakea. The focus of the flow within
the Laniakea core, the Great Attractor (Dressler et al., 1987a), is not well
specified. At best, we place it in the proximity of what we call the
4-clusters region (Courtois et al., 2013). The Centaurus and Norma clusters
are important constituents somewhat apart.

The flow histories are best represented by the accompanying videos and
Sketchfab interactive displays. It is interesting to take note of the
filamentary connections between the densest regions, particularly those
between Laniakea and Perseus-Pisces, such as the Centaurus-Puppis-PP filament.

Of course, while sectors of space can be parsed in co-moving coordinates
(Tully et al., 2014), the expansion of space in physical coordinates is
overwhelmingly dominant on all but very local scales in a flat $\Lambda$CDM
universe. The large scale structures continue to expand in one or two
dimensions becoming ever thinner. For the most part, halos have stopped
growing. However, our orbits, to the degree that they can be individually
believed, as extended in time until the scale of the universe doubles, and in
interesting regions, such as around the Perseus, Centaurus, and Virgo
clusters, present scenarios of a few merger events that will occur in the far
future.

##  8 acknowledgements

Above all, we recognise the critical contributions of P.J.E. (Jim) Peebles to
the development of the numerical action methodology. Jim modestly demurred
from being a co-author of this paper although it could be viewed as a
culmination of his inspiration. Any deficiencies in the exploitation rest with
the authors. We also thank our colleague Yehuda Hoffman who provided a file of
the 3D density field from an as yet unpublished quasi-linear modeling of
Cosmicflows-3 distances. This material is seen as the iso-density features of
the videos and interactive displays. Funding for the Cosmicflows project has
been provided by the US National Science Foundation grant AST09-08846, the
National Aeronautics and Space Administration grant NNX12AE70G, and multiple
awards to support observations with Hubble Space Telescope through the Space
Telescope Science Institute. The resources of the NASA-IPAC Extragalactic
Database and HyperLEDA hosted in Lyon, France have been indispensable.

## References

  * Abell et al. (1989) Abell, G. O., Corwin, Harold G., J., & Olowin, R. P. 1989, ApJS, 70, 1, doi: 10.1086/191333
  * Anand et al. (2021) Anand, G. S., Rizzi, L., Tully, R. B., et al. 2021, arXiv e-prints, arXiv:2104.02649.  https://arxiv.org/abs/2104.02649
  * Bardeen et al. (1986) Bardeen, J. M., Bond, J. R., Kaiser, N., & Szalay, A. S. 1986, ApJ, 304, 15, doi: 10.1086/164143
  * Böhringer et al. (2020) Böhringer, H., Chon, G., & Collins, C. A. 2020, A&A, 633, A19, doi: 10.1051/0004-6361/201936400
  * Carrick et al. (2015) Carrick, J., Turnbull, S. J., Lavaux, G., & Hudson, M. J. 2015, MNRAS, 450, 317, doi: 10.1093/mnras/stv547
  * Courtois et al. (2013) Courtois, H. M., Pomarède, D., Tully, R. B., Hoffman, Y., & Courtois, D. 2013, AJ, 146, 69, doi: 10.1088/0004-6256/146/3/69
  * de Lapparent et al. (1986) de Lapparent, V., Geller, M. J., & Huchra, J. P. 1986, ApJ, 302, L1, doi: 10.1086/184625
  * de Vaucouleurs (1953) de Vaucouleurs, G. 1953, AJ, 58, 30, doi: 10.1086/106805
  * di Nella et al. (1996) di Nella, H., Couch, W. J., Paturel, G., & Parker, Q. A. 1996, MNRAS, 283, 367, doi: 10.1093/mnras/283.2.367
  * Diemer (2018) Diemer, B. 2018, ApJS, 239, 35, doi: 10.3847/1538-4365/aaee8c
  * Diemer (2020) —. 2020, ApJ, 903, 87, doi: 10.3847/1538-4357/abbf52
  * Djorgovski & Davis (1987) Djorgovski, S., & Davis, M. 1987, ApJ, 313, 59, doi: 10.1086/164948
  * Dressler et al. (1987a) Dressler, A., Faber, S. M., Burstein, D., et al. 1987a, ApJ, 313, L37, doi: 10.1086/184827
  * Dressler et al. (1987b) Dressler, A., Lynden-Bell, D., Burstein, D., et al. 1987b, ApJ, 313, 42, doi: 10.1086/164947
  * Fairall (1998) Fairall, A. P. 1998, Large-scale structures in the universe
  * Ferrarese et al. (2012) Ferrarese, L., Côté, P., Cuillandre, J.-C., et al. 2012, ApJS, 200, 4, doi: 10.1088/0067-0049/200/1/4
  * Fields et al. (2014) Fields, B. D., Molaro, P., & Sarkar, S. 2014, arXiv e-prints, arXiv:1412.1408.  https://arxiv.org/abs/1412.1408
  * Gottlöber et al. (2003) Gottlöber, S., Łokas, E. L., Klypin, A., & Hoffman, Y. 2003, Monthly Notices of the Royal Astronomical Society, 344, 715, doi: 10.1046/j.1365-8711.2003.06850.x
  * Graziani et al. (2019) Graziani, R., Courtois, H. M., Lavaux, G., et al. 2019, MNRAS, 488, 5438, doi: 10.1093/mnras/stz078
  * Haynes & Giovanelli (1988) Haynes, M. P., & Giovanelli, R. 1988, Large-scale structure in the local universe: the Pisces-Perseus supercluster., ed. V. C. Rubin & G. V. Coyne (Princeton University Press), 31–70
  * Hoffman et al. (2021) Hoffman, Y., Nusser, A., Valade, A., Libeskind, N. I., & Tully, R. B. 2021, MNRAS, 505, 3380, doi: 10.1093/mnras/stab1457
  * Hoffman et al. (2017) Hoffman, Y., Pomarède, D., Tully, R. B., & Courtois, H. M. 2017, Nature Astronomy, 1, 0036, doi: 10.1038/s41550-016-0036
  * Hoffman et al. (2018) Hoffman, Y., Carlesi, E., Pomarède, D., et al. 2018, Nature Astronomy, 2, 680, doi: 10.1038/s41550-018-0502-4
  * Huchra et al. (2012) Huchra, J. P., Macri, L. M., Masters, K. L., et al. 2012, ApJS, 199, 26, doi: 10.1088/0067-0049/199/2/26
  * Jarrett et al. (2000) Jarrett, T. H., Chester, T., Cutri, R., et al. 2000, AJ, 119, 2498, doi: 10.1086/301330
  * Jarrett et al. (2003) Jarrett, T. H., Chester, T., Cutri, R., Schneider, S. E., & Huchra, J. P. 2003, AJ, 125, 525, doi: 10.1086/345794
  * Jasche & Lavaux (2019) Jasche, J., & Lavaux, G. 2019, A&A, 625, A64, doi: 10.1051/0004-6361/201833710
  * Jenkins et al. (2001) Jenkins, A., Frenk, C. S., White, S. D. M., et al. 2001, MNRAS, 321, 372, doi: 10.1046/j.1365-8711.2001.04029.x
  * Jensen et al. (2015) Jensen, J. B., Blakeslee, J. P., Gibson, Z., et al. 2015, ApJ, 808, 91, doi: 10.1088/0004-637X/808/1/91
  * Jones et al. (2009) Jones, D. H., Read, M. A., Saunders, W., et al. 2009, MNRAS, 399, 683, doi: 10.1111/j.1365-2966.2009.15338.x
  * Kaiser (1984) Kaiser, N. 1984, ApJ, 284, L9, doi: 10.1086/184341
  * Karachentsev & Kutkin (2005) Karachentsev, I. D., & Kutkin, A. M. 2005, Astronomy Letters, 31, 299, doi: 10.1134/1.1922528
  * Kashibadze et al. (2020) Kashibadze, O. G., Karachentsev, I. D., & Karachentseva, V. E. 2020, A&A, 635, A135, doi: 10.1051/0004-6361/201936172
  * Kauffmann & Fairall (1991) Kauffmann, G., & Fairall, A. P. 1991, MNRAS, 248, 313, doi: 10.1093/mnras/248.2.313
  * Keenan et al. (2013) Keenan, R. C., Barger, A. J., & Cowie, L. L. 2013, ApJ, 775, 62, doi: 10.1088/0004-637X/775/1/62
  * Kenworthy et al. (2019) Kenworthy, W. D., Scolnic, D., & Riess, A. 2019, ApJ, 875, 145, doi: 10.3847/1538-4357/ab0ebf
  * Kerr & Henning (1987) Kerr, F. J., & Henning, P. A. 1987, ApJ, 320, L99, doi: 10.1086/184983
  * Kitaura et al. (2021) Kitaura, F.-S., Ata, M., Rodríguez-Torres, S. A., et al. 2021, MNRAS, 502, 3456, doi: 10.1093/mnras/staa3774
  * Klypin et al. (2003) Klypin, A., Hoffman, Y., Kravtsov, A. V., & Gottlöber, S. 2003, ApJ, 596, 19, doi: 10.1086/377574
  * Kocevski et al. (2007) Kocevski, D. D., Ebeling, H., Mullis, C. R., & Tully, R. B. 2007, ApJ, 662, 224, doi: 10.1086/513303
  * Kourkchi & Tully (2017) Kourkchi, E., & Tully, R. B. 2017, ApJ, 843, 16, doi: 10.3847/1538-4357/aa76db
  * Kourkchi et al. (2020) Kourkchi, E., Tully, R. B., Anand, G. S., et al. 2020, ApJ, 896, 3, doi: 10.3847/1538-4357/ab901c
  * Kraan-Korteweg (2005) Kraan-Korteweg, R. C. 2005, Reviews in Modern Astronomy, 18, 48, doi: 10.1002/3527608966.ch3
  * Kraan-Korteweg & Lahav (2000) Kraan-Korteweg, R. C., & Lahav, O. 2000, A&A Rev., 10, 211, doi: 10.1007/s001590000011
  * Lahav et al. (1991) Lahav, O., Lilje, P. B., Primack, J. R., & Rees, M. J. 1991, MNRAS, 251, 128, doi: 10.1093/mnras/251.1.128
  * Leavitt & Pickering (1912) Leavitt, H. S., & Pickering, E. C. 1912, Harvard College Observatory Circular, 173, 1
  * Lee et al. (1993) Lee, M. G., Freedman, W. L., & Madore, B. F. 1993, ApJ, 417, 553, doi: 10.1086/173334
  * Lilow & Nusser (2021) Lilow, R., & Nusser, A. 2021, arXiv e-prints, arXiv:2102.07291.  https://arxiv.org/abs/2102.07291
  * Lucey et al. (1986) Lucey, J. R., Currie, M. J., & Dickens, R. J. 1986, MNRAS, 221, 453
  * Lynden-Bell et al. (1988) Lynden-Bell, D., Faber, S. M., Burstein, D., et al. 1988, ApJ, 326, 19, doi: 10.1086/166066
  * Magoulas et al. (2012) Magoulas, C., Springob, C. M., Colless, M., et al. 2012, MNRAS, 427, 245, doi: 10.1111/j.1365-2966.2012.21421.x
  * Mei et al. (2007) Mei, S., Blakeslee, J. P., Côté, P., et al. 2007, ApJ, 655, 144, doi: 10.1086/509598
  * Mieske & Hilker (2003) Mieske, S., & Hilker, M. 2003, A&A, 410, 445, doi: 10.1051/0004-6361:20031296
  * Mieske et al. (2005) Mieske, S., Hilker, M., & Infante, L. 2005, A&A, 438, 103, doi: 10.1051/0004-6361:20041583
  * Neill et al. (2014) Neill, J. D., Seibert, M., Tully, R. B., et al. 2014, ApJ, 792, 129, doi: 10.1088/0004-637X/792/2/129
  * Nusser & Branchini (2000) Nusser, A., & Branchini, E. 2000, MNRAS, 313, 587, doi: 10.1046/j.1365-8711.2000.03261.x
  * Paturel et al. (1996) Paturel, G., Garnier, R., Petit, C., & Marthinet, M. C. 1996, A&A, 311, 12
  * Peebles (1989) Peebles, P. J. E. 1989, ApJ, 344, L53, doi: 10.1086/185529
  * Peebles (1990) —. 1990, ApJ, 362, 1, doi: 10.1086/169237
  * Peebles (1993) —. 1993, Principles of Physical Cosmology
  * Peebles (1994) —. 1994, ApJ, 429, 43, doi: 10.1086/174301
  * Peebles (2010) —. 2010, arXiv e-prints, arXiv:1009.0496.  https://arxiv.org/abs/1009.0496
  * Peebles et al. (2001) Peebles, P. J. E., Phelps, S. D., Shaya, E. J., & Tully, R. B. 2001, ApJ, 554, 104, doi: 10.1086/321326
  * Phillips (1993) Phillips, M. M. 1993, ApJ, 413, L105, doi: 10.1086/186970
  * Planck Collaboration et al. (2020) Planck Collaboration, Aghanim, N., Akrami, Y., et al. 2020, A&A, 641, A6, doi: 10.1051/0004-6361/201833910
  * Pomarède et al. (2017) Pomarède, D., Hoffman, Y., Courtois, H. M., & Tully, R. B. 2017, ApJ, 845, 55, doi: 10.3847/1538-4357/aa7f78
  * Pomarède et al. (2020) Pomarède, D., Tully, R. B., Graziani, R., et al. 2020, ApJ, 897, 133, doi: 10.3847/1538-4357/ab9952
  * Pomarède et al. (2015) Pomarède, D., Tully, R. B., Hoffman, Y., & Courtois, H. M. 2015, ApJ, 812, 17, doi: 10.1088/0004-637X/812/1/17
  * Radburn-Smith et al. (2006) Radburn-Smith, D. J., Lucey, J. R., Woudt, P. A., Kraan-Korteweg, R. C., & Watson, F. G. 2006, MNRAS, 369, 1131, doi: 10.1111/j.1365-2966.2006.10347.x
  * Riess et al. (2021) Riess, A. G., Casertano, S., Yuan, W., et al. 2021, ApJ, 908, L6, doi: 10.3847/2041-8213/abdbaf
  * Riess et al. (2019) Riess, A. G., Casertano, S., Yuan, W., Macri, L. M., & Scolnic, D. 2019, ApJ, 876, 85, doi: 10.3847/1538-4357/ab1422
  * Rizzi et al. (2007) Rizzi, L., Tully, R. B., Makarov, D., et al. 2007, ApJ, 661, 815, doi: 10.1086/516566
  * Schaye et al. (2015) Schaye, J., Crain, R. A., Bower, R. G., et al. 2015, MNRAS, 446, 521, doi: 10.1093/mnras/stu2058
  * Schechter (1976) Schechter, P. 1976, ApJ, 203, 297
  * Shanks et al. (2019) Shanks, T., Hogarth, L. M., Metcalfe, N., & Whitbourn, J. 2019, MNRAS, 490, 4715, doi: 10.1093/mnras/stz2863
  * Shattow & Croton (2015) Shattow, G. M., & Croton, D. J. 2015, MNRAS, 452, 1779, doi: 10.1093/mnras/stv1464
  * Shaya et al. (1995) Shaya, E. J., Peebles, P. J. E., & Tully, R. B. 1995, ApJ, 454, 15, doi: 10.1086/176460
  * Shaya & Tully (2013) Shaya, E. J., & Tully, R. B. 2013, MNRAS, 436, 2096, doi: 10.1093/mnras/stt1714
  * Shaya et al. (2017) Shaya, E. J., Tully, R. B., Hoffman, Y., & Pomarède, D. 2017, ApJ, 850, 207, doi: 10.3847/1538-4357/aa9525
  * Sheth & Tormen (1999) Sheth, R. K., & Tormen, G. 1999, MNRAS, 308, 119, doi: 10.1046/j.1365-8711.1999.02692.x
  * Skrutskie et al. (2006) Skrutskie, M. F., Cutri, R. M., Stiening, R., et al. 2006, AJ, 131, 1163, doi: 10.1086/498708
  * Sorce et al. (2016) Sorce, J. G., Gottlöber, S., Yepes, G., et al. 2016, MNRAS, 455, 2078, doi: 10.1093/mnras/stv2407
  * Springel et al. (2006) Springel, V., Frenk, C. S., & White, S. D. M. 2006, Nature, 440, 1137, doi: 10.1038/nature04805
  * Springel et al. (2018) Springel, V., Pakmor, R., Pillepich, A., et al. 2018, MNRAS, 475, 676, doi: 10.1093/mnras/stx3304
  * Springob et al. (2014) Springob, C. M., Magoulas, C., Colless, M., et al. 2014, MNRAS, 445, 2677, doi: 10.1093/mnras/stu1743
  * Staveley-Smith et al. (2016) Staveley-Smith, L., Kraan-Korteweg, R. C., Schröder, A. C., et al. 2016, AJ, 151, 52, doi: 10.3847/0004-6256/151/3/52
  * Stein et al. (1997) Stein, P., Jerjen, H., & Federspiel, M. 1997, A&A, 327, 952
  * Tinker et al. (2008) Tinker, J., Kravtsov, A. V., Klypin, A., et al. 2008, ApJ, 688, 709, doi: 10.1086/591439
  * Tonry & Schneider (1988) Tonry, J., & Schneider, D. P. 1988, AJ, 96, 807, doi: 10.1086/114847
  * Tully (2015a) Tully, R. B. 2015a, AJ, 149, 54, doi: 10.1088/0004-6256/149/2/54
  * Tully (2015b) —. 2015b, AJ, 149, 171, doi: 10.1088/0004-6256/149/5/171
  * Tully et al. (2014) Tully, R. B., Courtois, H., Hoffman, Y., & Pomarède, D. 2014, Nature, 513, 71, doi: 10.1038/nature13674
  * Tully et al. (2016) Tully, R. B., Courtois, H. M., & Sorce, J. G. 2016, AJ, 152, 50, doi: 10.3847/0004-6256/152/2/50
  * Tully & Fisher (1977) Tully, R. B., & Fisher, J. R. 1977, A&A, 54, 661
  * Tully & Fisher (1987) —. 1987, Nearby galaxies Atlas (Cambridge: University Press, 1987)
  * Tully et al. (2019) Tully, R. B., Pomarède, D., Graziani, R., et al. 2019, ApJ, 880, 24, doi: 10.3847/1538-4357/ab2597
  * Tully et al. (1992) Tully, R. B., Scaramella, R., Vettolani, G., & Zamorani, G. 1992, ApJ, 388, 9, doi: 10.1086/171124
  * Tully & Shaya (1984) Tully, R. B., & Shaya, E. J. 1984, ApJ, 281, 31, doi: 10.1086/162073
  * van der Marel et al. (2012) van der Marel, R. P., Fardal, M., Besla, G., et al. 2012, ApJ, 753, 8, doi: 10.1088/0004-637X/753/1/8
  * Vogelsberger et al. (2014) Vogelsberger, M., Genel, S., Springel, V., et al. 2014, MNRAS, 444, 1518, doi: 10.1093/mnras/stu1536
  * Warren et al. (2006) Warren, M. S., Abazajian, K., Holz, D. E., & Teodoro, L. 2006, ApJ, 646, 881, doi: 10.1086/504962
  * Watson et al. (2013) Watson, W. A., Iliev, I. T., D’Aloisio, A., et al. 2013, MNRAS, 433, 1230, doi: 10.1093/mnras/stt791