Interacting galaxies in the IllustrisTNG simulations - VI: Reconstructed orbits, close encounters and mergers

Cosmological simulations have been used to study interacting galaxies as a function of galaxy pair separation, enabling comparisons with observational studies of galaxy pairs. The study of interacting galaxies as a function of time (i.e. merger stage) has mostly been limited to high resolution merger simulations, due to the poor time sampling available in cosmological simulations. Building on an earlier study of galaxy pairs in the IllustrisTNG cosmological simulations, we reconstruct the orbits of galaxy pairs involving massive galaxies ( 𝑀 ∗ > 10 10 𝑀 ⊙ ) at redshifts of 0 ⩽ 𝑧 < 1, using a novel kinematic interpolation scheme to model the orbits in between the IllustrisTNG snapshots (which are separated by 162 Myr on average). We assess the accuracy of these interpolations using a pre-existing suite of merger simulations, and find that kinematic interpolations provide a remarkable improvement in accuracy compared with interpolations that use only radial separations or 3D positions. We find that nearly 90 per cent of the closest pairs ( 𝑟 < 25 kpc) have had a pericentre encounter within the past Gyr. Many of these close pairs are found on rapidly shrinking orbits, and roughly 85 per cent of these pairs will merge within 1 Gyr. However, approximately 3 per cent of these close pairs appear to be flyby systems that will never merge. These reconstructed orbits will be used in future studies to investigate how and when galaxy properties change during close encounters and mergers between galaxies in IllustrisTNG.


INTRODUCTION
Early studies of galaxy properties showed that strongly interacting galaxies have morphologies, colours and star formation rates that are significantly different from those of relatively isolated galaxies (Arp 1966;Toomre & Toomre 1972;Larson & Tinsley 1978;Kennicutt et al. 1987).Large statistical studies of galaxy pairs have subsequently shown clear correlations between pair separation and numerous galaxy properties, suggesting that the degree to which galaxy properties are perturbed may be closely tied to the proximity of the encounters.For example, the star formation rates of galaxies in close pairs (with typical projected separations less than about 30 kpc) are found to be enhanced by a factor of roughly 2-3 with respect to relatively isolated galaxies (Ellison et al. 2008;Scudder et al. 2012;Barrera-Ballesteros et al. 2015;Stierwalt et al. 2015;Pan et al. 2018;Garduño et al. 2021;Steffen et al. 2021;Shah et al. 2022), with smaller enhancements extending out to separations of about 150 kpc (Patton et al. 2013(Patton et al. , 2020;;Brown et al. 2023).Galaxies in close pairs have also been found to have lower gas-phase metallicities (Ellison et ★ E-mail: dpatton@trentu.caal. 2008;Rupke, Kewley, & Chien 2010;Scudder et al. 2012;Bustamante et al. 2020;Garduño et al. 2021), higher asymmetries (Patton et al. 2005;De Propris et al. 2007;Casteels et al. 2014;Patton et al. 2016), and higher active galactic nuclei (AGN) fractions (Ellison et al. 2011(Ellison et al. , 2013;;Silverman et al. 2011;Satyapal et al. 2014;Ellison et al. 2019;Steffen et al. 2023;Bickley et al. 2024) than matched control samples of galaxies without close companions.All of these differences have been shown to diminish and then disappear at larger pair separations.
Simulations of merging galaxies have been used to interpret correlations between galaxy properties and pair separation, in attempts to identify the physical mechanisms by which these perturbed galaxy properties may arise.In general, these correlations can be explained by a model in which gravitational and hydrodynamical interactions trigger morphological disturbances and cause the infall of gas into the central regions of galaxies, leading to a dilution of gas-phase metallicity, enhanced star formation, and elevated accretion of gas on the central supermassive black hole (Mihos & Hernquist 1996;Di Matteo, Springel, & Hernquist 2005;Hopkins et al. 2008;Torrey et al. 2012;Renaud et al. 2014;Moreno et al. 2019).These simulations allow researchers to track changes in galaxy properties as a function of time, revealing when and how galaxy properties are altered by interactions and mergers.However, merger simulations are typically carried out on idealized pairs of galaxies on pre-defined orbits in otherwise isolated environments.In reality, the correlations seen in observed samples of galaxy pairs arise within a wide range of environments (Sol Alonso et al. 2006;Ellison et al. 2010;Lin et al. 2010), with a diversity of progenitor galaxies, galaxy orientations and orbital configurations.
In recent years, the availability of cosmological simulations has facilitated the theoretical study of more diverse samples of interacting galaxies, albeit at much lower resolution than merger simulations (Sparre & Springel 2016;Bustamante et al. 2018;Blumenthal et al. 2020;Hani et al. 2020;McAlpine et al. 2020;Patton et al. 2020;Bottrell et al. 2022;Byrne-Mamahit et al. 2023).Cosmological simulations offer an unprecedented opportunity for researchers to bridge the gap between high resolution simulations of galaxy mergers and large observational samples of galaxy pairs.In principle, cosmological simulations allow one to analyze correlations in galaxy properties as a function of pair separation and also to analyze changes in galaxy properties as a function of time, using the same set of simulated galaxies.In practice, however, the time sampling of cosmological simulations is sufficiently poor (typically > 150 Myr) that the orbits of interacting galaxies are crudely approximated at best.
The first paper in this series (Patton et al. 2020) reports on the identification of galaxy pairs within the IllustrisTNG cosmological simulations (Nelson et al. 2019a).In that study, the closest companion of each massive galaxy (with stellar masses of  * > 10 10  ⊙ ) was identified, with companions required to have at least 10 percent of the stellar mass of the galaxy in question.By comparing the specific star formation rates (star formation rates per unit stellar mass; hereafter sSFRs) of these galaxies with a control sample of relatively isolated galaxies matched on stellar mass, redshift and environment, Patton et al. (2020) found that galaxies with close companions (3D separations ∼ 15 kpc) have sSFRs that are enhanced by a factor of about two, with statistically significant enhancements detected out to separations of 280 kpc.By repeating their analysis in projected space, Patton et al. (2020) showed that these results are consistent with the correlations between sSFR and projected separation in galaxy pairs from the Sloan Digital Sky Survey (Abazajian et al. 2009).In a subsequent investigation, Brown et al. (2023) found that this consistency holds true when dividing the pairs sample based on whether the closest companion is star forming or passive.Additional papers in this series have reported on the star formation rates and AGN activity of paired and post-merger galaxies in IllustrisTNG (Hani et al. 2020;Quai et al. 2021;Byrne-Mamahit et al. 2023, 2024).
In this study, we aim to set the stage for the investigation of interacting galaxy properties as a function of time (rather than as a function of pair separation) by reconstructing the orbits of the Illus-trisTNG galaxy pairs identified by Patton et al. (2020).Our initial goals include the identification and characterization of the significant past and future close encounters that occur between the galaxies in each pair.For pairs that subsequently merge, we also determine the time of the merger.These reconstructed orbits can be used to determine how often close encounters and mergers occur between galaxies and their closest companions.In addition, these orbits can be used to determine if the differences between paired galaxies and their controls are indeed triggered by close encounters, as is typically assumed in studies of observed galaxy pairs.This will provide a direct link between studies of galaxy pairs in observations and simulations, facilitating the interpretation of galaxy pairs detected in observational studies.
In the following section, we outline our procedure for reconstruct-ing the orbits of interacting galaxies in IllustrisTNG, including our novel approach for interpolating the orbits in between simulation snapshots.In Section 3, we use a suite of merger simulations to investigate the accuracy of our interpolations.In Section 4, we describe the diversity of orbits in our sample, and quantify the time and separation of each galaxy pair's adjacent pericentres and apocentres, along with the time until each pair merges.In Section 5, we carry out a statistical analysis of the close encounters and mergers that are associated with the closest companion pairs in our sample.In the final section, we summarize our conclusions and describe future work that will benefit from the methods and/or reconstructed orbits presented in this study.

The IllustrisTNG Simulations
The IllustrisTNG project1 consists of a suite of cosmological gravomagnetohydrodynamical simulations that model dark matter, gas, stars and black holes within cubic cosmological volumes of approximate dimension 50 Mpc (TNG50), 100 Mpc (TNG100) and 300 Mpc (TNG300), adopting a cosmology consistent with the Planck 2015 results (Planck Collaboration et al. 2016).These simulations produce realistic galaxies in a wide range of environments, and they are well-matched to observations (Pillepich et al. 2018a).The public release of these simulations is described in Nelson et al. (2019a) and Nelson et al. (2019b).Additional details of the simulations can be found in Nelson et al. (2018), Springel et al. (2018), Pillepich et al. (2018b), Naiman et al. (2018) and Marinacci et al. (2018).
Our initial goal in this study is to reconstruct the orbits of the closest-companion galaxy pairs identified in the TNG100-1 simulation by Patton et al. (2020).The TNG100-1 simulation is the highest resolution run of the IllustrisTNG 100 Mpc volume.Given our desire to track interactions between large numbers of reasonably well-resolved galaxies spanning a wide range of environments, this simulation offers the best compromise between resolution and volume of the available IllustrisTNG simulations.

The Patton et al. (2020) TNG100-1 Closest Companion Sample
The Patton et al. (2020) pairs sample was assembled by extending the observational closest companion methodology of Patton et al. (2016) to cosmological simulations.Patton et al. (2020) found the closest companion for each massive ( * > 10 10  ⊙ ) galaxy at  < 1 (snapshots 50-99) in the TNG100-1 simulation, using the 3D separation (hereafter ) between the centre of potential of the two galaxies to make this determination.Companions were required to have stellar masses that are at least 10 percent as massive as the galaxy in question.With this definition, a companion may be more or less massive than the galaxy in question.
A key limitation in the study of close galaxy pairs in cosmological simulations is the numerical stripping that can occur when two galaxies have significant overlap in their stellar mass distributions.In these situations, the IllustrisTNG Subfind algorithm (Springel et al. 2001;Dolag et al. 2009) has a tendency to assign the outermost stellar mass particles of lower mass companions to their more massive neighbour (Rodriguez-Gomez et al. 2015).Patton et al. (2020) demonstrate that this effect can be minimized (although not eliminated) by replacing A histogram of the 3D distance to the closest companion for massive galaxies ( * > 10 10  ⊙ ) at  < 1 in TNG100-1 is shown, using the closest companion catalog of Patton et al. (2020).Each bin is 25 kpc wide.Closest companions are required to have at least 10 percent as much stellar mass as the galaxy in question.Galaxies that have no such companions within 2 Mpc are not shown.This sample is likely to be somewhat incomplete at the smallest separations (especially  < 25 kpc) due to the inability of the Subfind algorithm to disentangle highly overlapping galaxy pairs.each galaxy's current stellar mass with the maximum stellar mass that it has had within the past 0.5 Gyr.This revised determination of  * is used only in cases where the separation of the two galaxies is less than twice the sum of their stellar half mass radii.Patton et al. (2020) applied these stellar mass and closest companion criteria to all massive galaxies in IllustrisTNG at  < 1, yielding closest companion information for an average of 5976 massive galaxies at each snapshot.In Figure 1, we present a histogram of the 3D distance to the closest companion for massive galaxies ( * > 10 10  ⊙ ) in TNG100-1, as determined by Patton et al. (2020).This figure includes the 90.2 percent of the sample for which the galaxy's closest companion is closer than 2 Mpc (the remainder of the sample is not shown).The majority of galaxies do not have a companion that would qualify as "close" using the criteria of most studies of interacting galaxies (close pair samples are often restricted to separations of less than 30 kpc, and few studies of interacting galaxies extend beyond 100 kpc).In general, the number of galaxies increases with decreasing 3D distance to the closest companion.However, this trend reverses direction below 50 kpc, with fewer close pairs than would be expected based on inward extrapolation of the trend at larger separations.This turnover in the number of pairs may be due to the inability of the Subfind algorithm to correctly disentangle strongly overlapping galaxy pairs, despite the efforts taken by Patton et al. (2020) to minimize the effects of numerical stripping on their sample.

Galaxy Velocities
We will use the relative velocities of galaxies in IllustrisTNG pairs when reconstructing the orbits of interacting galaxies.In order to be useful, these velocities should accurately capture the relative motion of the centre of each galaxy in a pair at any given snapshot, with each galaxy's centre being determined by the location of the particle with the lowest gravitational potential energy.
The subhalo velocities that are provided in the IllustrisTNG public data release (Nelson et al. 2019a) are computed using the massweighted velocities of all particles in the subhalo.These centre-ofmass velocities are determined mainly by the dark matter in these galaxies.Unfortunately, the velocity (and position) of a galaxy's centre of mass can be significantly different from the velocity (and position) of the centre of the galaxy (i.e. the centre of potential).This distinction is particularly important when the galaxy is interacting with or accreting a nearby galaxy (Gómez et al. 2015;Vasiliev, Belokurov, & Evans 2022;Salomon, Libeskind, & Hoffman 2023).
We therefore determine our own centre-of-potential velocities for each massive galaxy in TNG100-1, computing the mass-weighted average velocity in the ,  and  directions for all stellar particles that are located within the galaxy's stellar half mass radius.We have confirmed that these velocities provide a much better match to the actual motion of each galaxy than the centre-of-mass velocities available from IllustrisTNG by comparing these velocities with the change in position of each galaxy over time from one snapshot to the next.

Assembly of Orbit Snapshots
We use the IllustrisTNG SubLink merger trees (Rodriguez-Gomez et al. 2015) to trace the orbital history and future of every galaxy pair in our closest companion sample.The SubLink merger trees connect Subfind subhaloes across snapshots, identifying the progenitors and descendants of each galaxy.
We step back in the merger tree for each galaxy by using its SubLink FirstProgenitor, which is the progenitor with the "most massive history" (Nelson et al. 2015).A galaxy's FirstProgenitor is usually located at the preceding snapshot.However, in cases where two galaxies become temporarily blended together at one snapshot, the SubLink algorithm may intentionally skip a snapshot, resulting in the FirstProgenitor occurring two snapshots earlier (Rodriguez-Gomez et al. 2015;Nelson et al. 2015).We follow each galaxy's sequence of FirstProgenitors (known as the Main Progenitor Branch) back in time until we reach a redshift of 1 (SnapNum=50) and/or until the stellar mass of either galaxy drops below 10 8  ⊙ (below this threshold, there are fewer than 70 stellar particles per galaxy, making estimates of stellar mass, position and/or velocity unreliable).
Similarly, we step forward in the merger tree for each galaxy by using its SubLink Descendant.Again, while the Descendant is usually one snapshot forward in time, snapshot skipping will occasionally cause the Descendant to be two snapshots forward in time.Where possible, we follow each galaxy's sequence of Descendants forward in time until a redshift of 0 (SnapNum=99).
In cases where the galaxies in the pair merge with one another, the two galaxies will have the same Descendant at all times following the merger.Following Hani et al. (2020), we take the time of the merger (hereafter  merge ) to be the time of the first snapshot at which the two galaxies have the same Descendant, while recognizing that the merger may in fact have occurred at any time between the final pair snapshot and  merge .Given that the average time between TNG100-1 snapshots at 0 ⩽  < 1 is 162 Myr (Patton et al. 2020), this means that the actual merger for any given pair will presumably have occurred approximately 80 ± 80 Myr earlier than  merge .
In Figure 2, we illustrate the results of this approach by plotting the three-dimensional separation () versus time relative to the present day (hereafter ) for a sample orbit in TNG100-1, covering the nearly 8 Gyr timespan corresponding to our redshift range of 0 ⩽  < 1.This pair has a separation of about 300 kpc at  ∼ 1, and merges nearly 7 Gyr (and 42 snapshots) later.For most of this time, this is not a closest companion system; in fact, the first time that this companion becomes the closest companion is at a lookback time of 2.0 Gyr, when the pair separation is 15 kpc, and the pair is 840 Myr away from merging.There are two places in the orbit where the separation increases, indicating that this pair undergoes at least two pericentre encounters before it merges.
While the series of snapshots depicted in Figure 2 gives us a rough idea of when these two galaxies begin interacting and when they merge, the relatively sparse time sampling of the IllustrisTNG snapshots severely limits our ability to identify and characterize the encounters between these galaxies.First of all, while it is clear that at least two encounters occur before the merger, we cannot say precisely when these encounters occur, nor do we know how close the galaxies were to one another at the pericentres of these encounters.In addition, it seems likely that we are missing additional close encounters in the latter stages of this merging sequence, as encounters typically occur more frequently as a merging pair approaches coalescence.These limitations are a common feature in studies of interacting galaxies in cosmological simulations (e.g., Blumenthal et al. 2020;Sotillo-Ramos et al. 2022;Łokas 2023).
In the following section, we address these limitations by introducing and comparing three different methods of interpolating orbits in between simulation snapshots.Our objective is to improve our ability to identify and characterize close encounters between interacting galaxies in IllustrisTNG.

Interpolation Between Snapshots
The following section describes three interpolation schemes for reconstructing the orbital path of a given galaxy pair between two consecutive snapshots. 2We define the time of the initial and final snapshots to be   and   respectively, and we define the corresponding 3D separations to be   and   .

1D Interpolation Between Snapshots
The simplest approach to interpolating between two consecutive snapshots is a linear interpolation in radial separation ().We refer to this approach as 1D interpolation.When applied to a plot of pair separation versus time (such as the orbit depicted in Figure 2), this is equivalent to simply drawing a straight line between each adjacent pair of snapshots.

3D Interpolation Between Snapshots
We can improve on simple 1D interpolation by taking full advantage of the three dimensional positional information that is available for each member of the galaxy pair at each snapshot.The simplest assumption is that the separation between the two galaxies changes linearly with respect to time in each of the ,  and  directions.We refer to this approach as 3D interpolation.Our expectation is that 3D interpolation should be an improvement over 1D interpolation, due to the incorporation of the additional 3D information.However, we will test this hypothesis in Section 3.

6D Interpolation Between Snapshots
Our 1D and 3D interpolation schemes account for the positions of both galaxies in a given pair, but they do not take the relative velocity of the system into account.We now introduce a kinematic interpolation scheme, which uses the separation and relative velocity of the pair as boundary conditions for integrating between two consecutive snapshots.We refer to this as 6D interpolation, since it requires 3D positions and 3D velocities.The purpose of this interpolation scheme is to find equations for Δ(), Δ(), and Δ() (and therefore  ()) that satisfy the boundary conditions.We compute a cubic Hermite spline in each coordinate direction, finding a unique third-order polynomial3 with respect to time .For completeness, we now outline this procedure in detail for Δ().
The pair separation in the -direction, Δ(), can be written as where   ,   ,   and   are constants.We have four constraints on this equation (Δ  , Δ  , Δ   and Δ   ), which allows us to solve for the four constants, as follows: Analogous equations for Δ() and Δ() have the same form, and are not shown.We can then compute  () by combining Δ(), Δ() and Δ() in quadrature.
This kinematic interpolation scheme is illustrated in Figure 3, which shows the interpolation in each of Δ, Δ and Δ, along with the resulting interpolation in , for the same sample pair shown in Figure 2, narrowing our focus to the final 3 Gyr of the merging sequence.This plot shows that Δ(), Δ(), Δ() and  () all vary smoothly with time, with the interpolations being continuous at each snapshot (a consequence of the boundary conditions in both separation and relative velocity).

1D, 3D and 6D Interpolation of a TNG100-1 Orbit
We have described three different interpolation schemes (1D, 3D and 6D) for predicting the orbital path between each pair of snapshots.In Figure 4, we apply all three interpolation schemes to the IllustrisTNG orbit that was depicted in Figure 2. In the early stages of this merger sequence, and surrounding the first apocentre, all three interpolation schemes yield similar results.On the other hand, these interpolation schemes yield strikingly different predictions near each apparent pericentre, as well as during the late stages of this merging sequence.
Which of these interpolations is best?Qualitatively speaking, the 6D interpolations appear to be the most realistic, with a smooth upand-down curve that is reminiscent of analogous "bounce" plots from idealized merger simulations (e.g.2022 ;Vasiliev, Belokurov, & Evans 2022).The smoothness of these interpolations as they pass through each snapshot is a consequence of our velocity boundary conditions, which force (), () and () to be continuous at each snapshot.This behaviour is what we would expect in the presence of smoothly varying gravitational forces.On the other hand, the 3D and 1D interpolations (which do not have velocity boundary conditions) are continuous but not smooth as they pass through each snapshot.
While there are good reasons to think that our kinematic (6D) interpolations are a major improvement over their 1D or 3D counterparts, we would like to confirm this by assessing the accuracy of all three interpolation schemes.More specifically, we would like to quantify the uncertainties in our predictions of the time and separation of close encounters found using these interpolations, as this is central to our goal of characterizing the encounters that occur between pairs of galaxies in IllustrisTNG.We therefore turn to a preexisting suite of idealized merger simulations to assess the accuracy of our interpolations, as described in the following section.

The Patton et al. (2013) Merger Suite
To investigate the accuracy of our IllustrisTNG interpolations, we apply our interpolation techniques to a suite of galaxy mergers that was introduced by Patton et al. ( 2013) and further analyzed by Moreno et al. (2015).For the purposes of this validation exercise, the key benefit of these simulations is the superior time sampling that is available (10 Myr in these merger simulations, versus an average of 162 Myr in IllustrisTNG).These merger simulations use the GADGET-2 N-Body/SPH simulation code of Springel (2005) to model the interaction and merger of two galaxies (with initial stellar masses of 1.4 × 10 10  ⊙ and 5.7 × 10 9  ⊙ ) placed on a variety of realistic merging orbits.Five choices of orbital eccentricities (0.85, 0.9, 0.95, 1.0 and 1.05) and impact parameters (2, 4, 8, 12 and 16 kpc) yield a set of 25 orbits that were chosen to be consistent with the orbital element distribution functions derived from cosmological simulations by Wetzel (2011).In addition, three sets of merger disc orientations were chosen, using the ,  and  orientations from Robertson et al. (2006).This yields a set of 75 merger simulations that exhibit a range of initial encounter separations, subsequent apocentre separations, and merger timescales, crudely mimicking the variety of merging orbits that are found in cosmological simulations.The output from these simulations was saved every 10 Myr, providing time sampling that is about 16 times better than for the IllustrisTNG snapshots.Please see Patton et al. (2013) for additional details about this suite of merger simulations.

Reconstruction of the Fiducial Orbit
We begin by examining the merger suite orbit that has intermediate values of impact parameter (8 kpc), eccentricity (0.95) and orientation (  ), calculating the pair separation at every 10 Myr snapshot throughout the interaction and merger sequence.We take this to represent the true orbit for this fiducial simulation.This orbit contains four close encounters, with the merger occurring approximately 2 Gyr after the first encounter.
Next, we create a mock IllustrisTNG orbit for this same fiducial simulation by sampling the orbit once every 160 Myr (i.e., at every 16th output in the merger simulation).We then apply our three interpolation schemes (1D, 3D and 6D) to this set of mock IllustrisTNG snapshots, and compare the reconstructed orbit with the true orbit.We pay particular attention to the accuracy of the interpolations during the first close encounter, which occurs at  = 0.50 Gyr and  = 12.3 kpc in the true orbit, where  is the time relative to the start of the simulation, and  is the 3D separation of the galaxies.
In Figure 5, we plot the 1D, 3D and 6D interpolations for this fiducial orbit, along with the true orbit.In the upper panel of Figure 5, we see that the 1D interpolation yields a single encounter that corresponds to the first pericentre of the true orbit.However, the separation at pericentre is equal to 33 kpc, which is nearly 3 times larger that the true pericentre separation.The 1D interpolation also does a poor job of identifying the time at which this pericentre occurs, with the time of the first pericentre found to be 70 Myr later than the true value.The 1D interpolation then completely misses all subsequent encounters, depicting a monotonic decline in separation until the system merges.
In the middle panel of Figure 5, the 3D interpolation is seen to provide a much better fit to the true orbit than the 1D interpolation, capturing the first two encounters reasonably well.The timing of the first encounter is consistent (within the 10 Myr time sampling) with the true orbit, and the separation is found to occur at 1.5 kpc (i.e.∼ 10 kpc smaller than the true separation).
In the lower panel of Figure 5, the 6D interpolation is seen to provide a substantial improvement over the 3D interpolation, and a vast improvement over the 1D interpolation.The 6D interpolation provides an excellent fit to the first encounter, correctly identifying the time, and estimating the separation to be 8.1 kpc (i.e.4.2 kpc smaller than the true separation).In fact, the 6D interpolation provides an excellent fit to all but the final 0.15 Gyr of the orbit.

Generalizing Over All Time Offsets
For a set of evenly spaced sparsely-sampled snapshots, the quality of the interpolation at a given close encounter is very sensitive to the proximity of the snapshots to the encounter.In particular, if a snapshot happens to fall close to the time of the pericentre, the resulting interpolation will be much more accurate at (and close to) that pericentre, for all three interpolation techniques.
In Figure 5, we intentionally offset our sparsely-sampled snapshots so as to ensure that none of the snapshots fell close to the first or second pericentres, thereby avoiding the possibility of highly accurate interpolations driven by good luck alone.We now generalize this approach by reconstructing the fiducial orbit 16 times, using all 16 possible offsets in time.This allows us to assess the accuracy of the interpolations around each pericentre in an unbiased manner, by including the best case scenarios (when the snapshot falls very close to the pericentre), the worst case scenarios (when the pericentre lies midway between the adjacent snapshots), and everything in between.
The results of this comparison are shown in Figure 6.In the upper panel of this plot, the places where the interpolation provides a poor fit to the true orbit appear as a black webbed pattern.The 1D interpolations perform worst in the vicinity of the first and second pericentres, as well as during the final stages of the merger.These 1D interpolations overestimate the separation at each pericentre, since the encounters identified from 1D interpolations are the local minima found within the available set of snapshots (i.e. in all but the ideal- ized case of the snapshot landing at the pericentre, the interpolation minimum will be larger than the true minimum).
In the middle panel of Figure 6, the 3D interpolations are seen to provide better fits to the true orbit than the 1D interpolations.The 3D interpolations almost always underestimate the pair separation in the vicinity of the first encounters and in the late stages of the merger sequence.The 3D interpolations consistently capture the second encounter, unlike the 1D interpolations.
In the lower panel of Figure 6, the 6D interpolations provide a much better fit to the true orbit than even the 3D interpolations.As with Figure 5, we again see that these interpolations usually underestimate the separation at each pericentre.Remarkably, most of these 6D interpolations also provide a good fit to the third close encounter.

Interpolation Accuracy at the First Pericentre
One of the primary goals of this study is to identify and characterize close encounters within the orbits of IllustrisTNG galaxy pairs.While Figure 6 provides useful insights in the accuracy of our interpolations for the fiducial orbit from Patton et al. (2013), we now extend our comparison to the full suite of 75 orbits from this merger suite, focusing on the accuracy of the interpolations at the first pericentre.In order to increase the precision of this analysis, we improve the time sampling of the true orbits by applying our 6D interpolation technique to the merger simulation snapshots (which are separated by 10 Myr).This allows us to reliably estimate the time of the true pericentre to the nearest Myr (interpolations are much more accurate within a small time interval).For each of the 75 orbits, we then create mock IllustrisTNG orbits (i.e., with 160 Myr time sampling) for all 16 time samplings.After reconstructing each of these 1200 mock orbits using 1D, 3D and 6D interpolations, we determine the time and separation of the first pericentre (i.e. the first local minimum in the interpolation), and then calculate the offset of these values from the true time and separation of the first pericentre.
In Figure 7, we plot the offset in pericentre separation ( peri −  true peri ) vs. the offset in pericentre time ( peri −  true peri ) for all 1200 reconstructed orbits.Above and to the right of this scatter plot, we include histograms of time and separation offsets respectively.In addition, we quantify the time and separation offsets for all three interpolation techniques by computing the median offsets for each, along with the standard deviation, reporting these results in Table 1.1D interpolations exhibit a wide spread of offsets in time and separation4 , with pericentre separations always overestimated (as expected).Compared with 1D interpolations, 3D interpolations have much smaller offsets, yielding roughly a factor of five improvement in time and a factor of three improvement in separation, with separation offsets that are almost always underestimated.6D interpolations are exceptionally good at identifying the time of the pericentre, with no significant median offset, and a spread of only 3.3 Myr.This is ∼ 15 times better than the 1D interpolations and ∼ 3 times better than the 3D interpolations.6D interpolations also excel at identifying the pericentre separation, with a small (negative) median offset of 1.4 kpc and a spread of only 1.8 kpc.

Success at Detecting Pericentres
Broadly speaking, Figure 6 indicates that all three interpolation schemes become increasingly inaccurate during the final stages of the merger, eventually leading to the failure to detect successive pericentres and apocentres, even for our 6D interpolations.These trends are consistent with what we would expect when we compare the timescales of the orbit to the time sampling of the snapshots.The orbital period (hereafter ) at any given stage in the merger sequence can be estimated by computing the time between successive pericentres or apocentres.In the final stages of the merger sequence,  becomes smaller than the average time between snapshots (∼ 160 Myr), posing serious challenges for any interpolation scheme.
In order to identify where our 6D interpolations are no longer able to detect every pericentre, we carried out a detailed inspection of our 1200 mock IllustrisTNG orbits (with time sampling of 160 Myr), comparing our reconstructed orbits with the true orbits.In cases where no more than one pericentre lies between an adjacent pair of snapshots, we find that the reconstructed orbits are always successful at detecting the pericentre.On the other hand, if there are two pericentres in between an adjacent pair of snapshots, the 6D interpolations often fail to detect both pericentres, with the likelihood of success depending on where the snapshots happen to lie in relation to the pericentres.These findings indicate that 6D interpolations are reliable at detecting pericentres as long as the orbital period is greater than the time between snapshots.
Of course, the challenge with interpreting orbits that have been reconstructed from sparely sampled snapshots is that we cannot accurately determine the orbital period during the late stages of the merger sequence.In particular, if the interpolation misses a close encounter, the orbital period may appear to be significantly larger than it is.To address this issue, we exploit the well-known relationship between orbital period and semi-major axis (hereafter ) that is enshrined in Kepler's third law.We investigate this relationship in the late stages of the merger sequence by determining  and  following the second, third and fourth pericentres in the Patton et al. (2013) merger suite.In each case, we compute the orbital period by doubling the time between the given pericentre and the subsequent apocentre, and we find the length of the semi-major axis by taking the average of the pair separation at the given pericentre and the subsequent apocentre.
We first apply this methodology to the true orbits in the merger suite (6D interpolations between 10 Myr snapshots), plotting the results as black stars in Figure 8.A clear correlation5 between  and  is seen, with orbital periods of ∼ 10-400 Myr and semi-major axes of ∼ 0.5-20 kpc.We note that the orbits that have  > 160 Myr all have  > 7 kpc.This suggests that if we were able to restrict our analysis to orbit segments with  > 7 kpc, we could avoid orbital periods of less than 160 Myr.
We next apply this methodology to our 1200 mock IllustrisTNG orbits, which sample the Patton et al. (2013) merger suite using 160 Myr snapshots.We plot the corresponding values of  and  in Figure 8, colour-coding the data based on whether the true orbital period is greater than 160 Myr (green circles) or less than 160 Myr (brown circles).The sparse time sampling of the reconstructed orbits induces a substantial spread in the estimates of both  and , particularly in the lower part of the plot, where the estimated orbital periods are much larger than their true orbital periods.However, in all cases where the estimated semi-major axis is greater than 10 kpc, the true (and estimated) orbital period is greater than 160 Myr.This suggests the following rule of thumb: 6D interpolations based on snapshots separated by 160 Myr should reliably detect all pericentres within snapshot intervals where the average pair separation is at least 10 kpc.
This rule of thumb confirms our impression from Figure 6 that 6D interpolations typically fail only in the very final stages of the merger sequence, when the average separation indeed becomes smaller than ∼ 10 kpc.In practice, the quality of the interpolations will depend on the details of a given orbit, including the dynamical masses of the galaxies (which affect the slope of the - correlation) and the specific spacing of the IllustrisTNG snapshots (which varies from 123 to 236 Myr at  < 1).We note also that this rule of thumb scales roughly linearly with the available time sampling (e.g., with 80 Myr snapshots, we have confirmed that pericentres can be reliably detected if  > 5 kpc).
In summary, our 6D interpolation technique provides a dramatic improvement in accuracy over 3D and especially 1D interpolation between snapshots, and is successful at detecting encounters as long as the average pair separation is ≳ 10 kpc.It therefore seems likely that 6D interpolations between IllustrisTNG snapshots should be capable of identifying most of the meaningful encounters between galaxies, while yielding estimates of the time and separation of these pericentres that are much more accurate than for 1D or even 3D interpolations.

DESCRIPTION AND CHARACTERIZATION OF RECONSTRUCTED ORBITS
We have shown that our 6D interpolation scheme is successful at identifying many close encounters in merger simulations with time sampling of 160 Myr.We now apply this technique to the orbits of galaxy pairs in IllustrisTNG, which have an average time between snapshots of 162 Myr at  < 1.We apply this technique to all TNG100-1 closest companion galaxy pairs with  < 2 Mpc, creating reconstructed orbits in every case where the galaxy pair (or corresponding post-merger) is identified at two or more snapshots. 6This procedure yields a sample of 32,319 unique reconstructed orbits.

Multiple Views of a Sample Reconstructed Orbit in IllustrisTNG
Each reconstructed orbit provides a record of the 3D separation () of a galaxy pair at each available IllustrisTNG snapshot, as well as the interpolated separation at each Myr in between snapshots.In addition to  (), the interpolations provide estimates of the separation in each of the ,  and  directions (see Equation 1).Finally, by differentiating Equation 1 with respect to time, we can determine the velocity in each of the ,  and  directions, and we can combine these quantities in quadrature to find the relative velocity of each pair at every Myr along the reconstructed orbit.Taken together, these quantities provide a wealth of information about the shape of the orbit and of the relative motion of each galaxy pair.
In Figure 9, we provide a visualization of this information by plotting multiple views of the TNG100-1 orbit that was initially portrayed in Figure 2. In the upper left hand panel, 3D separation is plotted versus time, extending the 6D interpolation from Figure 4 to the full time interval shown in Figure 2. In the upper right panel, we plot Δ versus Δ, showing the trajectory of this orbit within the  −  plane.The orbital path that is seen is approximately elliptical in shape, with the size of the ellipse shrinking during each close encounter.In the lower left panel, 3D velocity is plotted versus time, showing that the velocity peaks at around the time of each close encounter.In the lower right panel, 3D velocity is plotted vs. 3D separation, showing the interplay between these two quantities throughout the merging sequence.

Identification of Pericentres and Apocentres
One of the goals of this study is to identify and characterize the close encounters that are associated with each galaxy pair.In principle, it is straightforward to extract this information from our reconstructed orbits, with each local minimum (maximum) corresponding to a pericentre (apocentre).In practice, however, some orbits are much more complex than the one depicted in Figure 9, with small variations in separation producing subtle local minima and local maxima that are unlikely to be associated with any meaningful changes in galaxy properties.We therefore choose to restrict our analysis by requiring each pericentre to have a 3D separation that is at least 10 percent smaller than the 3D separation at the adjacent (preceding and following) apocentres.The locations of the corresponding pericentres and apocentres are identified on the reconstructed orbits that follow later in this manuscript.
Having identified pericentres and apocentres within each reconstructed orbit, we then situate each closest companion pair with respect to adjacent encounters within its orbit.Of paramount importance is the identification of the most recent pericentre encounter between the two galaxies, as that encounter may have had a strong influence on the properties of the galaxies in the pair.Using the reconstructed orbit, we compute the time since the most recent pericentre (hereafter Δ peri − ) and the separation of the pair at that pericentre (hereafter  peri − ).
Given the possibility that a given pair may be more strongly influenced by an imminent close encounter than any past encounters, we also compute the time until the next pericentre (hereafter Δ peri + ), along with the separation at that pericentre (hereafter  peri + ).Finally, we compute the corresponding times and separations of the adjacent apocentres (hereafter Δ apo − ,  apo − , Δ apo + and  apo + ).

The Diversity of Reconstructed Orbits in TNG100-1
The TNG100-1 orbit portrayed in Figure 9, like many other reconstructed orbits in our sample, is reminiscent of those seen in high resolution simulations of idealized galaxy mergers, whereby two galaxies merge in an environment that is devoid of other galaxies or large scale structure.In these scenarios, the galaxies in a pair undergo a series of progressively closer encounters, typically ending with a merger within 1-2 Gyr (Di Matteo et al. 2007;Cox et al. 2008;Lotz et al. 2008;Hopkins et al. 2013;Moreno et al. 2021).In reality, however, interacting galaxies do not normally evolve in isolation, and neighbouring galaxies or structures may exert a strong dynamical influence on a given galaxy pair (Moreno et al. 2013;An et al. 2019;Contreras-Santos et al. 2023).These more complex interactions between galaxies arise naturally in IllustrisTNG and other cosmological simulations, and can be seen in many of the reconstructed orbits in our sample.
In Figure 10, we illustrate the diversity of closest-companion orbits in IllustrisTNG by plotting 3D separation versus time for a representative set of 9 reconstructed orbits.The top row consists of three pairs that merge, including a fairly typical merging pair (Orbit a), a pair that merges rapidly after its initial infall (Orbit b), and a pair with a complex pre-merger sequence (Orbit c).The middle row depicts three pairs that have experienced multiple encounters but which appear unlikely to merge soon (if ever).These include a long-lived pair whose orbit is decaying slowly (Orbit d), a long-lived pair whose orbit does not appear to be decaying (Orbit e), and a pair that is undergoing strong interactions with other galaxies in its vicinity (Orbit f).The bottom row of Figure 10 shows a very close flyby encounter (Orbit g), a much wider flyby encounter (Orbit h), and an orbit that has a pair of encounters that are separated by more than 7 Gyr (Orbit i).

Single Companions Versus Multiple Companions
Throughout this analysis, we have described the orbits of closestcompanion galaxy pairs which -by definition -contain exactly two galaxies at all times.While some of the orbits in Figure 10 are described as being influenced by additional neighbouring galaxies, this may nevertheless have given the reader the mistaken impression that we are modelling the orbits of isolated galaxy pairs.In reality, galaxies often have many companions at a range of separations.For example, in Figure 4 of Patton et al. (2020), the average number of companions within 2 Mpc ( 2 ) ranges from 7 to 22 for galaxies in our TNG100-1 sample that have a closest companion within 1 Mpc.In addition, the "closest companion" designation may switch from one galaxy to another at different snapshots.This can occur if there are two or more companions with orbits whose radial extents overlap with one another, if the change in mass of a neighbouring galaxy moves it above or below the 10 per cent threshold used to identify companions (see Section 2.2), or if the closest companion merges with the galaxy in question, thereby passing on the title of closest companion to the next closest companion.
Figure 11 provides an illustration of the latter scenario, showcasing a galaxy that has five different closest companions over the redshift range of 0 ⩽  < 1.The reconstructed orbits of all of these companions are shown, with black circles used to identify the closest companion at each snapshot.Four of these companions merge with the given galaxy after undergoing at least two close encounters.After each of these mergers, the closest companion designation switches to a different companion that was further out and which was already on a decaying orbit.

CLOSE ENCOUNTERS AND MERGERS
In this section, we move from the examination of individual reconstructed orbits to a statistical analysis of the orbit sample as a whole.Our focus in this section is to study the prevalence and nature of close encounters and mergers in the Patton et al. (2020) sample of closest-companion galaxy pairs in IllustrisTNG.

Close Encounters
In assessing the extent to which a pair of galaxies may have been altered by a previous interaction, two key questions arise: how long has it been since the galaxies had their last pericentre encounter (if ever), and how close together (in space) were the galaxies at that time?If the encounter was sufficiently close, it is reasonable to expect that the galaxies may have undergone significant disruptions as a result.And if the encounter was relatively recent, one might expect that signs of any such disruptions might still be detectable.Pericentres and apocentres are identified using + symbols.For the four orbits that end in a merger, the time of the merger is depicted with a black triangle.At every snapshot in time, the closest companion is highlighted using a black circle.

The Times and Separations of Previous Encounters
Having computed the time since the most recent pericentre (Δ peri − ) and the pair separation at that pericentre ( peri − ) for every pair in our sample (see Section 4.2), we can tackle both of the questions posed in the preceding paragraph.In Figure 12, we plot  peri − versus Δ peri − for the 93,926 galaxies that have experienced a closest companion encounter within 500 kpc during the past 3 Gyr.The majority of these encounters have occurred relatively recently (within the past Gyr) and at separations that are sufficiently close (within about 100 kpc) that the galaxies in the pair would likely have had substantially overlapping dark matter haloes during the encounter.This suggests that a substantial fraction of the closest companion pairs in our sample may have been altered by relatively recent encounters.

The Previous Pericentre Fraction ( 𝑓 peri )
While the distributions of  peri − and Δ peri − shown in Figure 12 are instructive, more information is needed if we are to quantify the prevalence of close encounters in our sample.In particular, we must account for the fact that there is a limit to how far into the past we can track the orbit of each galaxy pair, and this will limit our ability to detect encounters at any given time in the past.In order to quantify the prevalence of recent close encounters in our sample of closest companion pairs, we therefore set out to compute the fraction of galaxies that have had a previous pericentre encounter with their closest companion (hereafter  peri ) as a function of lookback time (hereafter   ).Given the redshift limits of our pairs sample (0 ⩽  < 1), we can track some orbits nearly 8 Gyr into the past (those at very low redshift), while we have no available orbital history for pairs at snapshot 50 ( ∼ 1).We address this issue by computing  peri at a range of   , including only those orbits that have not yet dropped out of the sample at any given lookback time.
In Figure 13, we plot  peri versus lookback time for different ranges of 3D closest companion separations ().Pairs with smaller separations have the highest average  peri at all times.This is to be expected, since pairs with small  are more likely to be found near the pericentre of their orbit, and pairs with small  are also more likely to be found on smaller orbits with shorter orbital timescales.The pericentre fraction rises rapidly with   for the closest pairs ( < 25 kpc), reaching 88.2 ± 0.4 per cent by a lookback time of 1 Gyr.In other words, nearly every close pair has had a recent pericentre encounter with its closest companion.It is also notable that  peri reaches a plateau of about 97 per cent at 3 <   < 7 Gyr.This suggests that the remaining 3 per cent of close pairs never have a discernable first pericentre (e.g. if they merge very soon after their first close approach on a nearly radial orbit).
At wider pair separations,  peri drops quickly as a function of increasing , at any given lookback time.For example, only 49 per cent of pairs with 100 <  < 200 kpc have had a pericentre encounter within the past Gyr, and only 10 per cent of pairs with 500 <  < 1000 kpc have had a pericentre encounter within the past Gyr.And at the widest separations that we track (1000-2000 kpc), fewer than 20 per cent of the pairs have had a pericentre encounter within the past 7 Gyr.In summary, larger pair separations are associated with a decreased likelihood of the pair having had a recent pericentre encounter.

Merging vs. Non-merging Pairs
Our sample of reconstructed orbits contains a wealth of information about the types of interactions that take place between galaxies in IllustrisTNG.As seen in Figure 10, some closest-companion pairs are on merging orbits, while others appear unlikely to merge soon (if ever).We can assess the relative proportions of merging and nonmerging pairs in the sample as a whole by comparing the pericentre separations of the encounters that precede and follow each galaxy pair.For a sample of merging pairs, we would expect that the separation at the next encounter ( peri + ) would typically be smaller than the separation at the preceding encounter ( peri − ).Conversely, for a sample of non-merging pairs, we would expect no significant difference between  peri − and  peri + .
In Figure 14, we explore the relationship between  peri − and  peri + by plotting each of these quantities versus pair separation (in the upper and middle panels respectively) and by plotting  peri + / peri − versus pair separation in the lower panel.We restrict this sample to pairs with separations less than 500 kpc that have at least one pericentre in their past and in their future, thereby including the same set of galaxy pairs in all three panels.
In the upper panel of Figure 14, we can see that  peri − increases as a function of .In addition, while  peri − is usually smaller than , there is a distinct subset of close pairs that have  peri − > , which is what we would expect in a shrinking and potentially merging orbit (see e.g."Orbit a" in Figure 10).In the middle panel, wider pairs are found to have (on average) larger  peri + .In addition,  peri + is rarely larger than , especially at small separations.In these cases, the orbits are increasing in size (see e.g."Orbit f" in Figure 10).
In the lower panel of Figure 14, we plot  peri + / peri − versus pair separation, and compare with the null hypothesis of  peri + / peri − ∼ 1 that would be expected in the absence of mergers.At small separations (roughly  < 50 kpc), the median  peri + / peri − is approximately In the upper panel, the separation at the previous pericentre ( peri − ) is plotted vs. the 3D distance to the closest companion ().The dashed line denotes  =  peri − , and the solid blue line represents the median  peri − (with the line width corresponding to the 1 uncertainty in the median).In the middle panel, the separation at the next pericentre ( peri + ) is plotted vs. , with the dashed line denoting  =  peri + and the solid blue line representing the median  peri + .In the lower panel,  peri + / peri − is plotted vs. , with the dashed line denoting  peri + =  peri − (i.e.no change in minimum separation from one pericentre to the next), and the solid blue line representing the median  peri + / peri − .0.4; that is, the median pericentre decreases by about 60 per cent from one encounter to the next.This implies that many of these close pairs are on rapidly shrinking merging orbits.At larger separations, the median  peri + / peri − increases, nearly reaching the null hypothesis by  ∼ 500 kpc.This suggests that closest companion pairs with separations larger than about 0.5 Mpc are dominated by galaxies on non-merging or slowly decaying orbits.

Mergers
For the reconstructed orbits in which the pairs merge, we define the time of the merger ( merge ) to be the time of the first snapshot when the two galaxies have the same descendant (see § 2.4).By comparing this value with the time at which the pair is identified, we can compute the time until the pair merges.For pairs at the high redshift end of our sample, we can track pairs for nearly 8 Gyr into the future, providing plenty of opportunity to see if the pair merges.For pairs near  ∼ 0, there is instead only a short window of time within which a merger may occur.
Following the same methodology that was used in Section 5.1, we compute the fraction of pairs that will merge (hereafter  merge ) as a function of time, considering only those systems that have not yet dropped out of the sample at any given time.In Figure 15, we plot  merge versus time, for different subsets of 3D pair separation.We find that  merge rises rapidly for the closest pairs ( < 25 kpc), reaching 84.5 ± 0.4 per cent within 1 Gyr, and 92.6 ± 0.3 per cent within 2 Gyr.On longer timescales (> 5 Gyr),  merge appears to plateau at about 97 per cent, suggesting that roughly 3 per cent of the closest pairs are flyby systems that never merge.At larger pair separations, the merger fraction becomes progressively smaller at any given time.
In addition, there is an increasingly large time delay before  merge begins to rise for the wider pairs, as its may take billions of years for the galaxies to reach each other in pairs with large separations.
In Figure 16, we examine the relationship between pair separation and time relative to the merger (hereafter Δ merge ).Pairs that will merge in 1 Gyr have a median separation of 60 kpc, and 90 per cent of these pairs have separations of 18-136 kpc.On longer timescales, pairs that will merge in 5 Gyr have a median separation of 268 kpc, and 90 per cent of these pairs have separations of 65-694 kpc.A linear fit to the median separation versus Δ merge has a slope of −55 kpc/Gyr, giving us a sense of the typical time until the merger at a given separation (for merging pairs).

CONCLUSIONS AND FUTURE WORK
We have reconstructed the orbits of massive galaxies ( * > 10 10  ⊙ ) and their closest companions in the highest resolution run of the 100 Mpc IllustrisTNG cosmological simulation (TNG100-1).For each closest companion pair at  < 1 identified by Patton et al. (2020), we trace the orbit backwards in time (to a redshift of  ∼ 1) and forwards in time (to  = 0), using the 50 snapshots that are available.We interpolate between each pair of snapshots (which are 162 Myr apart on average) using a novel kinematic interpolation technique, and we use the merger simulations of Patton et al. (2013) to show that these 6D interpolations are much more effective at identifying and characterizing close encounters than 3D interpolations (which use separations in each of the ,  and  directions) or 1D interpolations (which use only radial separations).After analysing the resulting sample of reconstructed orbits, our main conclusions are as follows: (i) IllustrisTNG closest companion pairs include a diverse mix of rapidly merging pairs, slowly merging pairs, long lived satellites, flyby encounters, pairs in groups and clusters, etc.
(ii) almost 90 per cent of the closest pairs ( < 25 kpc) have undergone a pericentre encounter within the past Gyr (iii) galaxies in the closest pairs are often found on rapidly shrinking orbits, and roughly 85 per cent of these pairs will merge within 1 Gyr (iv) roughly 3 per cent of the closest pairs appear to be flyby systems that will never merge (v) the median separation of merging pairs decreases by ∼ 55 kpc/Gyr during the final 5 Gyr preceding the merger The results and trends reported here are likely to depend on various properties of the galaxy pairs, including the stellar mass and mass ratio of the galaxies, the environments that the pairs reside in, and the redshifts of the pairs.While our sample is large enough to investigate these dependencies, this analysis is beyond the scope of this study and is left to a future paper.
The reconstructed orbits presented in this work can also be used to study how galaxy properties change during close encounters, enabling us to uncover the physical processes that give rise to the correlations between galaxy properties and pair separation that have been reported in observed samples of galaxy pairs.By following merging galaxy pairs from their earliest close encounters through to their mergers and beyond, it will be possible to draw more direct connections between studies of pre-merger and post-merger galaxies in cosmological simulations.In addition, the properties of the orbits themselves (such as eccentricity and orbital period) can be used to analyze the dynamics of interacting galaxies in a realistic cosmological context, potentially relating the characteristics of these pre-merger orbits to the properties of the post-merger galaxies.
Finally, the kinematic (6D) interpolation scheme introduced in this study can also be applied to other cosmological simulations, allowing researchers to reconstruct the orbits of galaxies in between the available snapshots.While we have designed this approach with interacting galaxies in mind, these methods should be equally useful for researchers wishing to reconstruct the orbits of galaxies in groups and clusters of galaxies.
Figure1.A histogram of the 3D distance to the closest companion for massive galaxies ( * > 10 10  ⊙ ) at  < 1 in TNG100-1 is shown, using the closest companion catalog ofPatton et al. (2020).Each bin is 25 kpc wide.Closest companions are required to have at least 10 percent as much stellar mass as the galaxy in question.Galaxies that have no such companions within 2 Mpc are not shown.This sample is likely to be somewhat incomplete at the smallest separations (especially  < 25 kpc) due to the inability of the Subfind algorithm to disentangle highly overlapping galaxy pairs.

Figure 2 .
Figure 2. 3D separation () is plotted versus time relative to the present day () for a sample orbit in TNG100-1.Red circles denote the separation at each snapshot, and a black triangle indicates the time at which this pair merges ( merge ).

Figure 3 .
Figure 3. Separation is plotted versus time for a sample orbit in TNG100-1.The upper three panels show how the separations along the ,  and  axes change with time.The lowermost panel shows how the 3D separation changes with time.Within each panel, red circles denote the separation at each snapshot, the blue solid line depicts the kinematic (6D) interpolation between snapshots, the horizontal dashed line represents a separation that is equal to zero, and the black triangle denotes the merger.

Figure 4 .
Figure 4. 3D separation () is plotted versus time () for a sample orbit in TNG100-1.Red circles denote the separation at each snapshot, and a black triangle indicates the time at which this pair merges ( merge ).The black line depicts a linear (1D) interpolation of  with respect to , and yields two close encounters.The orange line depicts a 3D interpolation (positional interpolation in each of ,  and ), yielding four close encounters.The blue line depicts a 6D interpolation (kinematic interpolation in each of ,  and ), yielding four close encounters.

Figure 5 .
Figure 5. Three interpolation schemes are shown for the fiducial orbit from Patton et al. (2013), using a single offset in time.In all three panels, the red line shows the true orbit (10 Myr time sampling) and the red circles denote every 16th snapshot (with 160 Myr time sampling).The black line in the upper panel show a 1D interpolation between every 16th snapshot.The orange line in the middle panel show a 3D interpolation between every 16th snapshot.The blue line in the lower panel show a 6D (kinematic) interpolation between every 16th snapshot.

Figure 6 .
Figure 6.Three interpolation schemes are shown for the fiducial orbit from Patton et al. (2013), using all 16 offsets in time.In all three panels, the red line shows the true orbit.Each black line in the upper panel show a 1D interpolation between every 16th snapshot.Each orange line in the middle panel show a 3D interpolation between every 16th snapshot.Each blue line in the lower panel show a 6D (kinematic) interpolation between every 16th snapshot.

Figure 7 .
Figure 7.The accuracy of the 1D, 3D and 6D interpolation schemes are compared by plotting the offsets in time and separation at the first pericentre for all 16 time samplings of all 75 orbits fromPatton et al. (2013).A perfect interpolation would be located at the origin of this scatter plot.Colour coding is used to distinguish between 1D interpolations (black), 3D interpolations (orange) and 6D interpolations (blue).A histogram of  peri −  true peri is shown above the scatter plot, and a histogram of  peri −  true peri is shown to the right of the scatter plot.

Figure 8 .
Figure 8. Semi-major axis () is plotted versus orbital period () following the second, third and fourth encounters in the 75 orbits from the Patton et al. (2013) merger suite.Black stars depict the "true" values of  and , determined using 6D interpolations between snapshots separated by 10 Myr.Coloured circles show the corresponding values of  and  computed using 6D interpolations between snapshots separated by 160 Myr, with green circles for orbits with true  > 160 Myr and brown circles for orbits with true  < 160 Myr.

Figure 9 .
Figure9.A sample TNG orbit is depicted.3D separation is plotted versus time in the upper left panel.3D relative velocity is plotted versus time in the lower left panel.Relative position within the  −  plane is plotted in the upper right panel.3D relative velocity is plotted vs. 3D separation in the lower right panel.Within each panel, red circles denote the orbital properties at each snapshot, while the blue solid line depicts the kinematic interpoloation between snapshots.In all panels, horizontal and vertical dashed lines denote separations or velocities that are equal to zero.

Figure 10 .
Figure10.A sampling of orbits illustrates the diversity of interactions that are present in our sample of reconstructed orbits.In each panel, filled red circles denote snapshots, the solid blue line shows the 6D interpolation, + symbols locate pericentres and apocentres, and a black triangle denotes the time of the merger (if one occurs).The top row (Orbits a-c) depicts pairs that merge, with a variety of pre-merger orbits.The middle row includes pairs that are unlikely to merge soon (if ever) despite having experienced multiple encounters (Orbits d-f).The bottom row shows pairs with flyby or widely spaced encounters (Orbits g-i).

Figure 11 .
Figure 11.An example of a galaxy that has five different closest companions between z=1 and z=0.The orbits of all five companions are shown, using a different colour for each orbit.Within each orbit, the snapshots are depicted with filled circles, while the 6D interpolation is shown with a solid line.Pericentres and apocentres are identified using + symbols.For the four orbits that end in a merger, the time of the merger is depicted with a black triangle.At every snapshot in time, the closest companion is highlighted using a black circle.

Figure 12 .
Figure 12.The 3D separation at the previous pericentre ( peri − ) is plotted versus the time since the previous pericentre (Δ peri − ) in the main panel of this figure.The histogram at the top shows the distribution of Δ peri − (summed over the given range in  peri − ), while the histogram on the right shows the distribution of  peri − (summed over the given range in Δ peri − ).Closest companion pairs which have  peri − > 500 kpc or Δ peri − > 3 Gyr are relatively rare and are not shown.

Figure 13 .
Figure 13.The fraction of galaxies that have had a previous encounter with their closest companion (  peri ) is plotted versus lookback time (  ) for different ranges of 3D closest companion separations.The shaded regions depict 1 binomial errors on  peri .

Figure
Figure 14.In the upper panel, the separation at the previous pericentre ( peri − ) is plotted vs. the 3D distance to the closest companion ().The dashed line denotes  =  peri − , and the solid blue line represents the median  peri − (with the line width corresponding to the 1 uncertainty in the median).In the middle panel, the separation at the next pericentre ( peri + ) is plotted vs. , with the dashed line denoting  =  peri + and the solid blue line representing the median  peri + .In the lower panel,  peri + / peri − is plotted vs. , with the dashed line denoting  peri + =  peri − (i.e.no change in minimum separation from one pericentre to the next), and the solid blue line representing the median  peri + / peri − .

Figure 15 .
Figure15.The merger fraction (  merge ) is plotted versus time for subsets of the sample divided according to the 3D distance to the closest companion, as specified in the legend.The shaded regions depict 1 binomial errors on the merger fraction.The majority of galaxies with very close companions merge quickly.At larger separations, mergers are increasingly rare and occur much later.

Figure 16 .
Figure 16.Pair separation () is plotted versus time relative to the merger (Δ merge ) for closest companion pairs that merge by  = 0.The thick red line denotes median , and the thickness of the line depicts the 1 uncertainty on the median.The blue shaded region encloses the 95th and 5th percentile in .The black dashed line corresponds to a linear fit to the median, with a slope of −55 kpc/Gyr.

Table 1 .
Patton et al. (2013) time and separation at the first pericentre, computed using 16 time samplings of the 75 orbits ofPatton et al. (2013).The quoted uncertainties are standard deviations.