Multiple Beads-on-a-string: Dark Matter-Deficient Galaxy Formation in a Mini-bullet Satellite-satellite Galaxy Collision

Dark matter-deficient galaxies (DMDGs) discovered in the survey of ultra-diffuse galaxies (UDGs), in apparent conflict with standard CDM, may be produced by high-velocity galaxy-galaxy collisions, the $\textit{Mini-bullet}$ scenario. Recent observations of an aligned trail of $7-11$ UDGs near NGC1052, including DMDGs DF2 and DF4, suggesting a common formation event, $\sim8.9\pm1.5$ Gyr ago, provide a test. Hydro/N-body simulations, supplemented by galaxy orbit integrations, demonstrate that satellite-satellite collisions outside the host-galaxy virial radius can reproduce the observed UDGs in the NGC1052 group. A trail of $\sim10$ DMDGs is shown to form, including two massive ones that replicate the observed motions of DF2 and DF4. The linear relation, $v=Ax+v_{0}$, conjectured previously to relate positions ($x$) and velocities ($v$) of the aligned DMDGs as a signature of the collision event, is approximately obeyed, but individual DMDGs can deviate significantly from it. The progenitors whose collision spawned the trail of DMDGs survive the collision without, themselves, becoming DMDGs. We predict one progenitor is located at the end of the trail, testable by observing the difference between its stars, formed pre-collision, from those of the DMDGs, formed post-collision. By contrast, stellar ages and metallicities of the DMDGs are nearly identical. We further offer a hint that the tidal field of host NGC1052 may contribute to making DMDGs diffuse. $\Lambda$CDM simulation in a 100 cMpc box finds our required initial conditions $\sim10$ times at $z<3$. These results indicate current observations are consistent with the $\textit{Mini-bullet}$ scenario.


INTRODUCTION
In recent years, several galaxies have been observed to contain a lower amount of dark matter than predicted by galaxy formation theory in the standard Cold Dark Matter ("CDM") model.The latter posits that halo formation occurs in the pressure-free, collisionless dark matter prior to the gravitational infall of the baryonic component.On very small mass scales, below the baryonic Jeans mass of the pregalactic medium, gas pressure in the baryons can resist gravity, making the baryon mass fractions of the lowest-mass halos below the cosmic mean baryon-to-dark matter density ratio.For objects well above this baryonic Jeans-filter scale, the infall of baryons is supersonic and pressure forces are unimportant, so the baryons collapse along with the dark matter, Corresponding author: Joohyun Lee joohyun.lee@austin.utexas.eduand the baryon mass fraction inside virialized halos is close to that cosmic mean density ratio.When gaseous baryons are heated by feedback processes inside (e.g.SNe, AGNs) and/or outside (e.g.reionization) of the halo to which they would have been bound, pressure forces can suppress their infall or reverse it, resulting in a baryon-to-dark ratio well below the cosmic mean.However, in all these cases, the halos that form are dark matter-dominated.It was notable, therefore, when van Dokkum et al. (2018bvan Dokkum et al. ( , 2019) ) reported the existence of two ultra-diffuse galaxies (UDGs), NGC1052-DF2 and NGC1052-DF4 (hereafter DF2 and DF4, respectively), that are located in close proximity to the massive elliptical galaxy NGC1052 and exhibit a deficiency in dark matter, rather than the dark matter-dominance described above 1 .Subsequently, dark matter-deficient galaxies (DMDGs) have been identified across various environments and mass scales.These include the local group and isolated low-mass galaxies (Guo et al. 2020), a distant low-mass galaxy (Mancera Piña et al. 2022), and even a massive early-type galaxy in a cluster environment (Comerón et al. 2023).
The formation model that explains DF2 and DF4's dark matter deficiency should also address their exceptionally luminous globular cluster (GC) population (van Dokkum et al. 2018a(van Dokkum et al. , 2019)).To account for both phenomena at the same time, Silk (2019) proposed a "Mini-Bullet (cluster)" event or a "Bullet dwarf" scenario, in which a high-velocity (≳300 km s −1 ) collision of low-mass (dwarf) galaxies dissociates collisionless dark matter from baryons.As the name suggests, this scenario was inspired by the famous example of separation of dark matter and baryons observed in the Bullet Cluster, which has been explained by the collision of two cluster-scale halos at a high velocity, greater than either of their virial velocities, in which the collisionless nature of CDM and stars allows those components of each colliding halo to pass through each other, while the baryonic intracluster gas in each is prevented from doing so by its fluid behavior, resulting in a bow shock seen in X-ray emission (Tucker et al. 1998;Liang et al. 2000;Markevitch et al. 2002;Clowe et al. 2006).In the analogous Mini-Bullet collision, strong shock compression is induced, which in turn triggers star formation and the formation of massive star clusters.The large mass and narrow mass range of the observed star clusters are explained in this model by large gas surface densities that lead to a large lower limit to the initial cluster mass function while large galactic shear limits their mass range from above (Trujillo-Gomez et al. 2021).The galaxy collision produces these necessary conditions by strong radiative shocks that compress the gas and large-scale motions.
The separation of dark matter and baryons on various scales has been extensively studied on various scales, ranging from the GC scale (Kim et al. 2018;Ma et al. 2020b;Madau et al. 2020) to the galaxy cluster scale (Springel & Farrar 2007;Milosavljević et al. 2007;Mastropietro & Burkert 2008;McDonald et al. 2022), within the framework of the ΛCDM cosmology, using theoretical modeling and simulations.With regard to the DMDGs, we previously demonstrated this Mini-Bullet scenario using idealized galaxy collision simulations in Shin et al. (2020, hereafter Paper I), where we used two different N-body hydrodynamics codes with distinct numerical schemes to constrain the collision parameter space.Furthermore, we showed that massive star cluster formation is indeed triggered by high-velocity galaxy collision and that the star cluster properties from the simulation are roughly in line with the observations (Lee et al. 2021, hereafter Paper II).
The Mini-Bullet model is not the only one advanced so far to explain the DMDGs.One of the most frequently studied alternative mechanisms for the formation of these unique systems is tidal interaction that transfers dark matter to the more massive system (tidal stripping) (Ogiya 2018;Macciò et al. 2021;Jackson et al. 2021;Ogiya et al. 2022a;Moreno et al. 2022;Montero-Dorta et al. 2023;Mitrašinović et al. 2023;Katayama & Nagamine 2023).In the case of DF2 and DF4, however, although there are measurements of their tidal distortion, the results do not require that the galaxies' dark matter was removed by tidal interaction (Keim et al. 2022).Müller et al. (2019) also argued that there is no sign of stellar streams induced by tidal interaction near DF2 and DF4, which is also claimed by Montes et al. (2021) with respect to DF2.On the other hand, Montes et al. (2020) claimed that DF4 is undergoing tidal disruption.In any case, the tidal interaction scenario cannot explain the exceptionally bright GC population as a natural outcome of the scenario.
Another scenario is the tidal dwarf galaxy ("TDG") formation mechanism, in which DMDGs formed from efficiently cooled and fragmented gas after it was ejected during a strong tidal encounter with a disk galaxy (Recchi et al. 2007;Duc et al. 2014;van Dokkum et al. 2019;Haslbauer et al. 2019;Fensch et al. 2019).This mechanism, too, has, so far, not been shown to involve the simultaneous formation of bright GCs, however.
Another idea, suggested by Trujillo-Gomez et al. (2022), explains the dark matter deficiency as a consequence of the formation of those bright GCs and their back-reaction on the gaseous galactic baryons, which in turn modified the dark matter distribution.They argued that powerful stellar feedback from massive GC populations can induce a gravitationally-coupled expansion of the dark matter content, reducing its contribution to the dynamical mass of the galaxy (also see Li et al. 2023, for similar work on the response of dark matter to gas ejection).
Recently, van Dokkum et al. (2022a) presented a new clue that supports the "Mini-bullet" scenario in the case of DF2 and DF4 by measuring their line-of-sight velocities.The authors conclude that both UDGs were formed from a single event that occurred about 8 Gyr ago, likely a high-velocity galaxy collision, which has been shown to be capable of producing the observed lack of dark matter and bright GC populations (Paper I; Paper II).Using the HST observation of DF2 and DF4 and a catalog of low-surface brightness galaxies in the NGC1052 group studied by Román et al. (2021), the authors identified an alignment of 7 − 11 UDGs as a possible "trail of DMDGs" (see Figure 1 in van Dokkum et al. 2022a).This argument is further substantiated by follow-up observation performed by the same group, which revealed that the GCs of DF2 and DF4 have the same color (van Dokkum et al. 2022b).Additionally, age measurements of the GCs and stel-lar bodies by Fensch et al. (2019) (DF2 only) and Buzzo et al. (2022) (DF2 and DF4) yielded consistent results of ∼ 8 Gyr for the age of the GCs and stellar bodies.
More recently, Buzzo et al. (2023) studied the large-scale structure of GCs in the NGC1052 group and found that the GC distribution is consistent with scenarios involving a single galaxy-galaxy interaction event and subsequent coeval formation of the GCs and the DMDGs, which includes the tidal dwarf galaxy scenario and the Mini-bullet scenario.To distinguish the two scenarios, a possible "smoking gun" signature of the Mini-bullet scenario was suggested by Gannon et al. (2023), a linear relationship between the line-of-sight velocities of the aligned DMDGs produced by the collision and their distance from their common point of origin in the collision; the further a galaxy was from this point, the larger must its launch velocity have been to reach that distance in the same elapsed time.The authors noted that the known velocities and projected distances from NGC1052 of DF2 and DF4 were consistent with such a linear relationship.To test this further, Gannon et al. (2023) used spectroscopy to measure the stellar age, stellar metallicity, and line-of-sight velocity of one more DMDG in the trail, NGC1052-DF9 (DF9).They concluded that the age and metallicity of the galaxy are similar to those of DF2 and DF4, consistent with the Mini-bullet scenario, but the observed line-of-sight velocity of DF9 deviates significantly from the expected linear relationship.
In response to the observations, Ogiya et al. (2022b) claimed the Mini-bullet scenario may face challenges due to the strong tidal forces exerted by the host galaxy, NGC1052, which can strip GCs from the formed DMDGs shortly after their formation.This argument is supported by two main factors.Firstly, the spatial distribution of the observed GCs in DMDGs DF2 and DF4 is extended, and taking into account their orbital decay due to dynamical friction, their formation epoch-distribution should have been even more extended than their current distribution (Dutta Chowdhury et al. 2019, 2020).Secondly, the galaxy-galaxy collision occurs at or within the virial radius of NGC1052.The authors argue that the combined effect of the extended distribution of GCs and the strong tidal field exerted by NGC1052 makes the Mini-bullet scenario implausible.
In this paper, we shall explore the ability of the Mini-bullet scenario to account for the enigmatic characteristics of the UDGs in the NGC1052 group, their dark matter deficiency, and alignment.Our investigation will be carried out through a series of gravitohydrodynamic simulations and galaxy orbit integrations using the ENZO code and the Rebound code, respectively.The initial conditions of the simulations will be designed to match the observed physical properties, including stellar masses and kinematics, and alignment of the NGC1052 group UDGs.Our results will demonstrate that appropriate initial structural and orbital parameters of the colliding satellite progenitor galaxies can produce a "trail of DMDGs" that includes two massive DMDGs with M ⋆ > 10 8 M ⊙ corresponding to DF2 and DF4, whose motions agree with the observed values.We will show that while the positions and velocities of the DMDGs on the trail generally follow a linear relationship, as previously suggested to be a signature of their collision origin, there can be deviations in the positions and velocities of individual DMDGs from that simple relation.We will find that the stellar ages and metallicities of the DMDGs are nearly identical, but we will also examine the scatter in their values.We will compare the simulated DMDGs with observed UDGs and discuss which physical processes need to be taken into account.We will also quantify the occurrence of such Mini-bullet events in the Universe using a large simulation of galaxy formation from cosmological initial conditions, TNG100-1.
This paper is structured as follows.In Section 2, we describe our effort to constrain the initial conditions by using idealized galaxy-galaxy collision simulations and backward (i.e.time-reversed) orbit integration.Section 3.1 presents the simulation results, including the stellar masses and orbits (positions and velocities) of the DMDGs that formed.In Section 3.2, we discuss the stellar properties, stellar metallicities, ages, and sizes of the product DMDGs.We compare these results with previous observational and theoretical work in Section 4.1.Section 4.2 discusses the statistical likelihood of the Mini-bullet satellite-satellite galaxy collision in a large cosmological simulation ILLUSTRISTNG.Our summary and conclusions are presented in Section 5.

SIMULATIONS
We present a three-step methodology aimed at aligning hydrodynamic simulations with the observational findings of van Dokkum et al. (2022a), focusing on the formation of multiple DMDGs through a single Mini-bullet collision between two progenitor galaxies orbiting around the massive host halo of NGC1052.Our primary objectives are to match (i) the stellar mass of DF2 and DF4, the two most massive DMDGs among the NGC1052 group UDGs, (ii) the line-ofsight velocity difference of DF2 and DF4, (iii) the positions of DF2 and DF4, and (iv) the number of resultant DMDGs.
In the first step, we conduct idealized galaxy-galaxy collision simulation experiments similar to what we have done in Paper I and Paper II to explore the parameter space of progenitor galaxy properties, such as the dark matter halo mass (M DM ), gas mass (M gas ), gas distribution, and collision configuration, including relative collision velocity (v col ) and pericentric distance (r min ), that can produce ∼ 10 aligned DMDGs after the collision.As the second step, utilizing the information obtained from the previous step regarding collision configurations capable of generating aligned DMDGs and their associated positions and velocities, we conduct backward orbit integrations to determine the initial conditions for the hydrodynamic simulations based on the observation of the NGC1052 group UDGs.Finally, in the third step, we perform hydrodynamic simulations using the established initial conditions to examine the feasibility of the Mini-bullet scenario in producing a trail of DMDGs and predict their characteristics.In the following sections, we elaborate on each step in detail.

Confining parameter space-idealized galaxy-galaxy collision simulations
In the idealized galaxy-galaxy collision simulations, we use the publicly available adaptive mesh refinement (AMR) code ENZO (Bryan et al. 2014;Brummel-Smith et al. 2019) with its hydrodynamics solver ZEUS (Stone & Norman 1992a,b).The GRACKLE library (Smith et al. 2017) is used to compute radiative gas cooling and heating assuming the homogeneous ultraviolet (UV) background of Haardt & Madau (2012) at z = 0, by interpolation of the lookup table generated from the CLOUDY code (Ferland et al. 2017) 2 .We refine the simulation box of (2.621 Mpc) 3 down to a spatial resolution of ∆x = 10 pc (level l max = 12), which is eight times coarser than the simulations carried out in Paper II.It is necessary here to reduce the computational cost of simulating a larger box over a longer time than before, so we must relax the resolution, aware that it cannot resolve star cluster formation as was a primary goal of Paper II, with its spatial resolution of ∆x = 1.25 pc, but is sufficient to resolve the galaxy formation properties we need here.Our refinement strategy is super-Lagrangian, meaning if a cell contains more mass than the mass threshold (i.e. its mass density exceeds a density threshold), that cell splits into eight child cells.At refinement level l, for gas with star formation threshold gas number density n th , M l ref, gas = 2 −0.333(l−12) × M 12 ref, gas , where M 12 ref, gas = 10000 M ⊙ = n th ∆x 3 ≃ 2.5 M Jeans (T = 100 K, n = 400 cm −3 ), and for (dark matter and stellar) particles, M l ref, part = 2 −0.107(l−12) × M 12 ref, part , where M 12 ref, part = 16000 M ⊙ .We note that, while refinement proceeds down to a cell size that is small enough not to contain more than a fixed physical mass of either gas or particles, this refines the force length resolution of the gravity and gas pressure forces but does not refine the particle mass resolution.Dark matter and stellar particles (after releasing stellar feedback) have a fixed mass, at all times and do not refine.
To model feedback-regulated star formation, we adopt subgrid algorithms for under-resolved small-scale physical processes, just as we did in Paper II, but with some adjustments for coarser resolution and an assumption of higher thermal energy released by SNe associated with the massive stars.A parcel of gas is determined to form stars, according to the approach in Cen & Ostriker (1992), and the outcome of that star formation is assumed to follow the star-forming molecular cloud model (SFMC; for details, see Kim et al. 2013Kim et al. , 2019)).In brief, an SFMC particle is created when the following criteria are met: (1) the density of a gas cell exceeds n th = 400 cm −3 , (2) the gas flow is converging, (3) the cooling time of the cell is shorter than its dynamical time (t dyn ), and (4) the mass in that cell is enough to create an SFMC particle heavier than m SFMC = 5 × 10 3 M ⊙ (which leads to a "permanent" star particle mass -i.e. total mass that follows the initial mass function of the stars that form, of m ⋆, new = 10 3 M ⊙ ).The SFMC particle returns 80% of its original mass to the gas in a time equal to 12t dyn , according to our assumed star formation efficiency for converting gasto-star mass (Krumholz & Tan 2007), along with supernova thermal feedback of 10 51 erg per 50 M ⊙ of permanent star mass that peaks at 1t dyn and 2% of metal yield (see also Kim et al. 2011).
We follow the method used in Paper I and Paper II to initialize two progenitor galaxies by utilizing the DICE code (Perret 2016).We place two identical progenitor galaxies 60 kpc apart and set their relative velocity to be 300 − 600 km s −1 .Since the observed line-of-sight velocity difference of DF2 and DF4 is 358 km s −1 (van Dokkum et al. 2022a), relative collision velocities need to be higher than that.The pericentric distances r min is set to be 2 kpc.
While the parameters for these idealized galaxy-galaxy collisions in Step 1 are similar to those adopted in Paper I, there are several important distinctions between the simulations here in Step 1 and those in Paper I. First, the progenitors here are taken to be spherical halos with gaseous baryons, which self-consistently form stars before their collision, while in Paper I, the intention was to model the collision of present-day galaxies with well-established disks.Second, the spatial resolution here is much higher than in Paper I, which makes a difference in the ability of the collision simulations to produce the post-collision stellar systems.To model the NGC1052 galaxy group, it is necessary to demonstrate that a single collision with realistic parameters can naturally produce a trail of ∼10 DMDGs.The coarser spatial resolution of simulations in Paper I, of 80 pc (as opposed to this paper's 5 pc resolution) prevented us from forming ∼ 10 DMDGs there, however.In those lower-resolution simulations, fewer objects of higher mass were formed, in general, and it was necessary to tune the choice of parameters just so as to maximize this number.Even so, only up to 6 DMDGs were formed.In the current paper, we believe we have converged with spatial resolution, as demonstrated by the comparison of the runs with 5 and 10 pc resolution, respectively.NOTE-(1) run name, (2) relative velocity of the two progenitors at 60 kpc distance, (3) pericentric distance (i.e.distance at closest approach), (4) scale radius of the PIS gas density profile, (5) the total mass of a progenitor, (6) gas fraction f gas = M gas /M total , (7) stellar mass of the most massive DMDG formed, (8) the largest relative velocity difference between the DMDGs, (9) the number of formed DMDGs, (10) time since the pericentric approach of the two progenitor galaxies when we end the simulation.
Figure 1.Illustrative Mini-Bullet collision, before and after.Snapshots of the time-history of the ENZO simulation of an idealized collision of two identical gas-rich dwarf galaxies, each with M total = 1.89 × 10 10 M ⊙ , with a collision velocity of 500 km s −1 (with black arrows indicating the progenitors' moving directions; Run 3 in Table 1), at t = −110 (initial time-slice of simulation), 90, 390, 690 Myr.t = 0 is set to the moment when the two dwarf galaxies are at pericenter (i.e.closest approach).Surface densities of dark matter (top row), gas (middle row), and only those stars that formed after the start of the simulation (110 Myrs before the orbits of the colliding haloes reached pericenter) (bottom row) are presented.All projections are conducted in 100 kpc-depth layers.After the galaxy collision, seven DMDGs are formed and survive (t = 690 Myr; fourth column).
And the formation of ∼ 10 DMDGs did not require such finetuning as before, either.
Paper II, on the other hand, had even higher spatial resolution than the simulations here, but only by applying the AMR to a limited "zoom-in" region surrounding the most massive DMDG formed by an idealized galaxy-galaxy collision.As such, it did not address the questions at issue here, of producing a trail of ∼10 DMDGs.
The progenitors are initialized with only gas and dark matter.Their dark matter halos have M 200 = 1.66 × 10 10 M ⊙ , J 200 = 1.03 × 10 12 M ⊙ kpc km s −1 , and R 200 = 51.7 kpc, and follow the NFW profiles (Navarro et al. 1997), with con-centration parameter c = 133 .For the gas density profiles, instead of initializing the progenitor galaxies with exponential disks as in Paper I and Paper II, this time, we adopt the pseudo-isothermal (PIS) profile without any rotation, with a scale radius of R s, gas = 2 kpc, which is more commonly observed and used to simulate low-mass satellite galaxies, to initialize the gas density profile of the progenitors (e.g., Kurapati et al. 2020).The total mass is 1.89 × 10 10 M ⊙ , and the gas fraction is f gas = M gas /M total = 12.4%, with no stars at the beginning. 4The gas is initially set to a temperature of 10 4 K and metallicity of Z = 0.1 Z ⊙ = 0.002041 to match the metallicity of a low-mass galaxy at z ∼ 1 − 2. We  The surface density of the gas (top row) and density-weighted average gas temperature with velocity vectors overplotted (bottom row) are presented.Gas velocity vectors are projected onto the image plane, which is the yz plane, in which the center of the mass of the colliding progenitor system was at t = 0 Myr.The velocity vectors are proportional to the length of the velocity vector arrows and the largest arrow in each panel corresponds to ∼200 km s −1 .Gas is projected in 10 kpc-depth layers.
use 2 × 10 6 dark matter particles to model each progenitor, resulting in mass resolution of m DM = 8.29 × 10 3 M ⊙ .
We perform a suite of four galaxy-galaxy collision simulations, each with a different set of initial configurations.In Table 1, we summarize this set of structural parameters we test and the stellar masses of the most massive DMDG that result, the velocity difference of the most and second massive DMDGs (corresponds to DF2 and DF4), and the number of DMDGs produced at t end , after the collisions.Among the tested results, it is notable that there are collision configurations that result in a large velocity difference (∼ 300 km s −1 ) between the first and second most massive DMDGs that form which is similar to that between DF2 and DF4 (van Dokkum et al. 2022a).We note that, in these idealized collision simulations, the post-collision separation velocity does not increase with time, while in the more realistic simulations we will describe below, in which the progenitors are also satellites of the massive galaxy NGC1052, their postcollision separation velocity can increase with time.As a result, it is possible that their immediate post-collision separation velocity was somewhat less than observed now for DF2 and DF4.
The time history of the mini-Bullet system from pre-to post-collision is illustrated by Figures 1, 2, and 3, which show snapshots of the results of our ENZO simulation for Run 3 in Table 1.At t = −110 Myr, two progenitor galaxies are initialized with a separation of 60 kpc, approaching at 500 km s −1 .Figure 1 zooms out for time-slices long before, during, and long after the progenitor collision, to show how the Mini-bullet collision dissociates collisionless components, dark matter and stars, of these progenitors from their gas (t = 90 Myr in Figure 1).In Figure 2, we zoom in, both in time and space, to show the immediate effects of the collision on the two progenitors and the distinctions between those effects on their collisionless components (dark matter and stars), which can pass through each other, versus their collisional component, their gaseous baryonic interstellar media (ISM), which, as a fluid, cannot.Instead, the ISM of each progenitor is halted from its supersonic approach to its counterpart in the other progenitor, by a strong shock that heats the gas to T > 10 7 K, in which the sound speed is comparable to the collision velocity (t = −30 and 0 Myr in Figures 2 and 3).There are two shocks, one on each side of their collision center, one to decelerate each incoming progenitor's ISM, producing a layer of shock-heated gas between them, whose mass grows as the collision proceeds to overtake more and more of each progenitor's ISM.The shock-heated gas begins cooling radiatively as soon as it is heated, and the shock becomes a radiative one, for which the cooling time is short compared with the time for the shock velocity to change (t = 30 Myr in Figures 2 and 3).In such shocks, the gas can cool isobarically, to well below 10 4 K, increasing its density in inverse proportion to the temperature, so its density increases as the square of the shock Mach num-ber.Since the gas in the ISM of each progenitor is centrallyconcentrated, the cooling time is shortest for shocked gas that was initially closer to the centers of the progenitors, so this is where the shocks become radiative first, and the cooling gas is subject to dynamical and gravitational instabilities.For this off-axis collision, dense gas clumps are formed due to radiative cooling near the line that connects the progenitors (y = 0 line).This results in a combination of cold gas which collapses towards their centers while surrounding shells of still-hot shocked gas are driven to expand away from the line of collision, by pressure forces, as shown in Figure 3, from the view along the line of collision (x-axis).After star formation in the collision-induced dense gas clumps, we identify seven DMDGs that are roughly aligned between the progenitors, forming an "S-shaped" trail (t = 390 Myr in Figure 1).Gas as the fuel of star formation is almost exhausted before this time, preventing star formation in the DMDGs on the trail.The trail of DMDGs is tidally stretched over time by the gravitational field of the separating progenitors at each end, which retain their dark matter and stars, making the "Sshaped" trail longer as the separation between galaxies grows (t = 690 Myr in Figure 1).

Producing "Tailor-Made" initial conditions-backward orbit integration
Equipped with this knowledge from our suite of idealized galaxy-galaxy collision simulations regarding which parameters are suitable for producing DMDGs like DF2 and DF4, we now move to establish a realization of colliding progenitor pairs which are satellites of a massive host galaxy.This realization should produce two DMDGs, each with M ⋆ ≳ 10 8 M ⊙ (for DF2 and DF4) and several UDGs along a line after ∼ 8 Gyr, to match the observations (Román et al. 2021;van Dokkum et al. 2022a).These observations will include the radial separations and positions on the sky of DF2 and DF4, their radial velocity separation, and their radial velocities relative to that of NGC1052.
The unknown initial conditions that will lead to these final conditions for DF2 and DF4 are the pre-collision locations and velocity vectors of the collision progenitors, relative to each other and to NGC1052.To derive these we will integrate the orbits of DF2 and DF4 backward in time to locate the collision that produced them and their velocity vectors at that time.By doing a large enough sample of such orbit integrations, we were able to identify those for which the orbits of DF2 and DF4 once crossed in the past, indicating the time and place of the collision we postulate to have formed them.From this subset of the orbit integrations that led to collisions in the past, we refined our sample further, to those for which the lookback time of the collision was large enough to explain the observationally inferred ages of the GCs of DF2 and DF4, i.e. ∼ 8 Gyr (Ma et al. 2020a;Fensch et al. 2019;van Dokkum et al. 2018c).
We noticed that, for this refined set of cases, the relative velocities of DF2 and DF4 immediately after the collision that produced them, were typically lower than the presentlyobserved radial separation velocities, by ≳ 100 km s −1 , reflecting the relative acceleration of their post-collision orbits over time.This means that the immediate post-collision separation velocities were consistent with the collisions in Runs 2 and 3 in Table 1, thereby confirming that such post-collision outcomes were, indeed, possible from suitably high-velocity galaxy-galaxy collisions.
The next step was to select three of these orbit integrations that yielded the collisions at the right lookback time, and determine the pre-collision parameters of the progenitor galaxies that would produce these immediate post-collision separations and velocities for DF2 and DF4.In principle, if we had a large enough sample of cases in Table 1, we might have had a close enough match to the orbit integrations that we could use them to identify the pre-collision configuration that would lead to the post-collision configuration derived from a given orbit integration.In practice, however, it is computationally prohibitive to run so many cases of idealized collisions for this purpose.Instead, we determined the pre-collision configuration as follows.
The exact pre-collision configurations of the two progenitors are decided by setting relative collision velocities v col to two times the relative velocities of DF2 and DF4 ∆v DF2−DF4 , resulting in the two progenitors collide with ∼ 500 km s −1 , and pericentric distances of r min ∼ 2 kpc.This follows the observation from the simulations presented in Section 2.1 that the velocity difference of the first and second most massive DMDGs formed in ∆v DMDG, max2 is roughly a half of v col of the progenitors ∆v DF2−DF4 = ∆v DMDG, max2 ∼ 0.5v col (with ignorance of the long-term velocity change due to the presence of the host gravity; refer to Section 2.1 and Table 1).
Finally, with the immediate pre-collision configurations derived in this way for each of the three selected cases above, another backward time integration was required for each, this time of the orbits of their two progenitor galaxies, destined to meet at their collision moment, to place them at large separations at earlier times, from which we could perform full gravitohydrodynamic simulations of their collision events, from start to finish.Now we will describe these steps in more detail.
We use the Integrator with Adaptive Step-size control, 15th order (IAS15), a gravitational dynamics integrator, implemented in the Rebound orbit integration code (Rein & Liu 2012;Rein & Spiegel 2015) to backtrace the orbits of DF2 and DF4 based on the currently observed quantities and find the collision point of the progenitor galaxies, where DF2 and DF4 should have started to form at the same time.Then, us-ing the confined collision parameters that can produce multiple aligned DMDGs (Section 2.1), we place progenitor galaxies to make initial conditions for the hydrodynamic simulation (Section 2.3).Finding the collision point is not trivial, because depending on the initial locations and motions, the back-traced DF2 and DF4 may not collide.Therefore, we should find the initial conditions from the parameter space that consists of the range of current locations and motions of DF2 and DF45 .
In Figure 4, we present an example of the three initial conditions based on the orbit integration calculation results.We place a dark matter-only massive host (black circle) with a halo mass of M DM = 6.12 × 10 12 M ⊙ to represent the host, NGC1052, at the center (Forbes et al. 2019).The orbits of DF2 and DF4 with the stellar mass of M ⋆ = 2 × 108 M ⊙ need to be time-reversed integrated to find the collision point.The radial distances of DF2 and DF4 are constrained with a certain amount of uncertainties (22.1 ± 1.2 Mpc for DF2 and 20.0 ± 1.6 Mpc for DF4) (van Dokkum et al. 2018b(van Dokkum et al. , 2019(van Dokkum et al. , 2022a) ) 6 , but the distance to NGC1052 is not well measured.Thus, in this study, we will assume the distance to NGC1052 as 19 Mpc consistent with previous studies using the surface brightness fluctuation method (Tonry et al. 2001;Blakeslee et al. 2001) and the Virgo-infall corrected radial velocity (Gil de Paz et al. 2007) 7 .The transverse distances of DF2 and DF4 from NGC1052 can be simply calculated from their angular separations (Román et al. 2021): −0.081 Mpc (DF2) and +0.168 Mpc (DF4) assuming the observed distance to DF2 and DF4.For the radial distances, we set errors to be 0.25×(observation error) and for the transverse distances, we adopt 0.002 Mpc as errors.The three galaxies and the observer are not exactly on the same plane but we ignore the deviation and assume that the transverse (proper) and radial (line-of-sight) distances can be used as the x-axis and the yaxis in our orbit integration as if they are on the same plane.The radial velocities of DF2 and DF4 relative to NGC1052 are also measured with unspecified errors (van Dokkum et al. 2022a): +315 km s −1 (DF2) and −43 km s −1 (DF4).Since the observation errors for radial velocities are not specified, we set their errors to be 10 km s −1 .For those known param-eters (radial and transverse distances and radial velocities), we let them vary within observed values ± errors.However, since the transverse motion of DF2 and DF4 cannot be observed, we vary the transverse velocities within fixed ranges: [10 km s −1 , 170 km s −1 ] for the progenitor 1 (bound satellite, blue dashed in Figure 4) and [0 km s −1 , 80 km s −1 ] for the progenitor 2 (unbound satellite, red dashed in Figure 4) with an interval of ∆v = 2.5 km s −1 .To confine the parameter space, we uniformly sample the parameter space, testing > 1, 500, 000 cases, and find ∼ 3000 cases where DF2 and DF4 started from places closer than 10 kpc, the place of a Mini-bullet galaxy collision.
After finding the collision point, we place two progenitor galaxies (blue and red circles on dashed lines) with dark matter halo mass of M DM = 1.42 × 10 10 M ⊙ and M gas = 2.47 × 10 9 M ⊙ ( f gas = 0.148) with no stars 8 , which is similar to the mass we tested, using information from the idealized two-galaxy collision simulations in Section 2.1.One of the progenitors (dashed blue) is gravitationally bound to the host and has a highly eccentric orbit.The other progenitor (dashed red) is not bound to the host but comes from outside of the NGC1052 system with a high velocity of > 400 km s −1 relative to the host, which is higher than the host virial velocity Using these orbit integration results, we set up initial conditions consisting of the massive host and two progenitors at 0.3 Gyr before their Mini-bullet collision, to make them relax after the artificial starbursts that occur after initialization.Their gas density distributions are given by the PIS profile, while their dark matter densities follow an NFW profile.Using the DICE code (Perret 2016), the density profile of NGC1052-like host is sampled with 10 6 dark matter particles (m DM, host = 6.1 × 10 6 M ⊙ ) and the colliding progenitors' NFW profile is realized with 10 7 dark matter particles (m DM, prog = 1.42 × 10 3 M ⊙ ) each.Since the lookback time of 8 Gyr corresponds to z ∼ 1, we will now set the progenitors' halo concentration to be c = 7, about a factor of two lower than c = 13, the value which was adopted in the idealized simulations in Table 1.This is because, for a halo with a given mass, the concentration parameter is lower for a halo that formed earlier.We summarize sets of collision configurations and structural parameters of the progenitors we test in Table 2. Henceforth, we shall refer to all of the initial conditions described here in Section 2.3 as "Tailor-Made" and the cases in Table 2 as "Tailor-Made" cases (e.g.TM1 shall refer to Tailor-Made case 1).Among those initial conditions, we specifically focus on TM1, TM2, and TM3 for our analysis presented in Section 3.1 and 3.2, as representative runs that start Table 2. Initial conditions for simulations of satellite-satellite galaxy collisions around a massive host with M ⋆ = 6.12 × 10 12 M ⊙ .The initial conditions are "tailor-made" in attempt to produce collisions whose outcome matches observations of the NGC1052 group (van Dokkum et al. 2022a), including backward orbit integration described in Section 2.2.10) orbits of DF2 and DF4 in backward orbit integration calculation described in Section 2.2.The difference in the orbits originates from different starting positions and velocities of DF2 and DF4 in the time-reversed orbit integrations, as listed in columns ( 6)-( 9), accounting for observational error and uncertainties in transverse velocities of DF2 and DF4, which cannot be observed.Runs (TM4−8) share the starting positions and velocities of DF2 and DF4 that define Orbit 1 with TM1.To be clear, all these "starting" values for backward time integration are actually "final" values when time is expressed as going forward to the present epoch.

Run name
from different pre-collision orbital parameters of the progenitor galaxies and yield results that reproduce the observations.(Note: TM5 is a variation on TM1, with different separation distance at the closest approach for the colliding progenitors, but with a similar enough outcome to TM1's that we do not give a more detailed description here.

Mini-bullet satellite-satellite galaxy collision and the formation of multiple dark matter-deficient galaxies
With these Tailor-Made initial conditions, summarized in Table 2, we again use the ENZO code to simulate the highvelocity collision of two satellite galaxies and the orbital evolution of the progenitors and resultant DMDGs after the collision.Hydrodynamics and star formation physics are the same as the idealized simulations in Section 2.1.The box size is increased to (5.243 Mpc) 3 to fully capture the orbits of the product DMDGs, but the nested refinement region is set to resolve only the progenitors and product DMDGs, not the massive host galaxy.Static refinement up to level 3 (∆x = 10.24 kpc) is applied outside the nested refinement regions.We change the refinement scheme due to the refinement level change (∆x = 5 pc at level l max = 14) and refine more for the particles to resolve star-dominated structures better: at refinement level l, for gas, M l ref, gas = 2 −0.264(l−14) × M 14 ref, gas , where M 14 ref, gas = 4000 M ⊙ ≃ 2.5 M 100 K Jeans , and for par-ticles, M l ref, part = 2 −0.273(l−14) × M 14 ref, part , where M 14 ref, part = 4000 M ⊙ .Following the change of cell spatial resolution and refinement criteria, we change the star formation density threshold to n th = 1.6 × 10 3 cm −3 and the SFMC particle mass threshold to m SFMC = 2.5 × 10 3 M ⊙ (permanent star particle mass of m ⋆, new = 5 × 10 2 M ⊙ ).The simulations are run for 2.3 Gyr, 0.3 Gyr before the collision, and 2 Gyr after the collision at which we analyze the properties of the product DMDGs.Due to the significant computational time needed to run the full ∼ 8 Gyr of simulation, we further supplement the simulations with orbit integration code Rebound for the remaining ∼ 6 Gyr to track the orbital evolution of the DMDGs and colliding progenitors during later stages.

RESULTS
In this section, we present the results of the idealized satellite-satellite galaxy collision simulation introduced in Section 2.3 and discuss how the product DMDGs can be compared to the observation and what the physical properties of the DMDGs are.
3.1.Orbits of the produced galaxies: a trail of dark matter-deficient galaxies and progenitors Figures 5 and 6 depicts a time sequence of the Mini-bullet collision of two satellite galaxies orbiting around a massive host, which corresponds to NGC1052 (see Section 2.2 and  2) stellar mass of the most and second-massive DMDGs ("-" = no DMDG formed), (3) line-of-sight velocity difference of the most and second massive DMDGs at t = 2 Gyr ("-" = no DMDG formed), ( 4) line-of-sight velocity difference of the most and secondmost massive DMDGs obtained from orbit integration until the distance between the two most massive DMDGs that correspond to DF2 and DF4 reaches 2.1 Mpc apart ("N/A" corresponds to the runs that do not involve orbit integration due to their lack of the number of DMDGs), ( 5) the number of formed DMDGs with M ⋆ > 10 7 M ⊙ , (6) time since the pericentric approach of the two progenitor disks when the ENZO simulation ends, (7) time since the pericentric approach of the two progenitor disks when the orbit integration ends, i.e., when the distance between the two most massive DMDGs is 2.1 Mpc, (8) whether the results from TM ("Tailor-Made") produce enough DMDGs to match the number of UDGs observed in the NGC1052 group.2.3 for the exact host properties), and the resulting formation of spatially aligned DMDGs in the TM1 run in Table 2 and  Table 3.The time of the two colliding satellite galaxies at the pericentric approach is defined as t = 0.At t = 2 Gyr, which is the last snapshot of the simulation, the progenitors are approximately 1 Mpc apart, and the sequence of DMDGs spans over 600 kpc.
The DMDGs are identified by the HOP halo finder (Eisenstein & Hut 1998).We employ the boundary of each DMDG to be the tidal radius, or the Jacobi radius, the distance at which the gravitational force of the massive host becomes equivalent to the DMDG's self-gravitation, defined by Binney & Tremaine (2008): where r host is the distance of the DMDG from the massive host, M DMDG is the mass of the DMDG, and M host = 6.12 × 10 12 M ⊙ is the mass of the massive host.
The middle row of Figure 5 illustrates the gas distribution during and after the galaxy collision.Consistent with our findings in previous studies (Paper I; Paper II), the baryonic component of the progenitor galaxies is stripped from their dark matter halos and undergoes severe shock compression during the collision process (as shown in the second, t = 0 Gyr panel of the upper right panels).This results in a burst of star formation inside the compressed dense gas clouds (see Figures 1 and 3 in Paper II, for more details).Meanwhile, the stripped gas is tidally elongated, and multiple self-gravitating gas clumps form along the line connecting the two progenitors, generating stars inside these clumps.As a result, after .Snapshots from the ENZO simulation of two colliding satellite galaxies of NGC1052 from TM1 run ("Tailor-Made" fiducial run) initial conditions at t = −0.1 Gyr, t = 0 Gyr, t = 1 Gyr, and t = 2 Gyr, centered on the two colliding satellite progenitor galaxies in the reference frame of their center of mass.t = 0 is set to the moment when the two progenitors are at a pericentric approach.Surface densities of dark matter (top row), gas (middle row), and stars (bottom row) in 100 kpc-depth layers are presented.Black arrows in the top-left panel indicate the approach of the progenitor galaxies toward each other before they collide.In the top row, white arrows on each plot point to the direction of the host galaxy NGC1052, centered at the distances labeled, and the length scale in each panel is indicated by the distance rulers of 40 kpc overplotted there.A few hundred Myr after the collision, almost no gas is left.At t = 1 Gyr and t = 2 Gyr, a sequence of DMDGs can be observed.approximately 1 Gyr, a number of similar stellar structures are identifiable (as shown in the bottom panel of the third column in Figure 5, t = 1 Gyr).At this point, the fuel for star formation is depleted due to both the star formation burst (major) and gravitational infall to the host (minor), causing the stellar mass and population of the DMDGs to remain almost constant.Consequently, between t = 1 Gyr and t = 2 Gyr, the gravitational field of the host and the progenitors is the primary driver of the dynamical evolution of the DMDGs and the progenitors.
As summarized in Table 3, the results from the simulations successfully satisfy two of the goals listed in Section 2: the stellar mass of DF2 and DF4 and the number of product DMDGs.We find ∼ 10 (20; 30) self-gravitating stellar structures with M ⋆ > 10 7 M ⊙ (M ⋆ > 10 6 M ⊙ ; M ⋆ > 10 4 M ⊙ ) in the simulations we analyze in detail: TM1 (Tailor-made1), TM2, and TM3 (refer to Table 3 for the exact setups).In our  3 for more detail.
analysis, we restrict our focus to self-gravitating star clumps with M ⋆ > 10 7 M ⊙ and designate them as DMDGs.Table 3 also shows the velocity differences between the two most massive DMDGs, which can be compared to the observed line-of-sight velocity difference of DF2 and DF4 at the end time of the simulations (t = 2 Gyr for TM1, TM2, and TM3 run; corresponding to ∼ 6 Gyr ago from the present), in the TM1 run and other simulations we test.The velocity differences are smaller, but orbital motions for ∼ 6 Gyr need to be taken into account.
Notably, the formation of DMDGs with M ⋆ > 10 8 M ⊙ and multiple DMDGs is highly sensitive to the structural parameters of the progenitors.Too small r min and large R s, gas result in the gathering of gas near the collision point, yielding too massive DMDGs and a small velocity difference between the most and second massive DMDGs.Eventually, the DMDGs gravitationally pull each other and no trail of DMDGs is produced.Furthermore, the concentration of progenitor dark matter halo also affects the DMDG formation.More DMDGs are formed in the runs with lower concentra-tion parameters, meaning that the strong tidal field during the Mini-bullet collision can suppress gas and newly-born stars from gathering and forming self-gravitating structures.
In order to compare our simulation results with observations of the NGC1052 group, which by design, are expected to exhibit good alignment, we further examine the orbital evolution of the DMDGs and their progenitor satellite galaxies by integrating their orbits with the Rebound code. Figure 7 displays the orbits of the product DMDGs and the progenitors, tracked in the hydrodynamic simulations from t = −0.3Gyr to t = 2 Gyr (solid lines) and in the orbit integration from t = 2 Gyr to t ∼ 8 Gyr (dashed lines).At the latter time, the line-of-sight distance (y-coordinates in Figure 4, 5, and 6 and x-coordinates in Figure 1 of van Dokkum et al. 2022a) between the two most massive DMDGs, which correspond to DF2 and DF4, is 2.1 Mpc apart, is set to the observed amount in van Dokkum et al. (2022a).As summarized in Table 3, in the TM1 (TM2; TM3) run, this occurs at t = 8.375 (6.65; 8.502) Gyr.The line-of-sight velocity difference of the two most massive DMDGs after the orbit integration ∆v DMDG, max2, now = 324 in the TM1 run (554, 343 km s −1 in TM2, TM3 run) is roughly in line with the observed velocity difference of DF2 and DF4, 358 km s −1 (van Dokkum et al. 2022a), complying with the remaining goals of matching the velocity difference and positions of DF2 and DF4 listed in Section 2.
Conserving the alignment of the DMDGs at t = 2 Gyr depicted in Figures 5 and 6, our orbit integration analysis reveals that the DMDGs remain spatially aligned after orbit integration, consistent with the observed peculiarity of the UDGs in the NGC1052 group (Román et al. 2021;van Dokkum et al. 2022a).In several runs with different initial orbital and structural parameters of the colliding satellites, we consistently confirm that with appropriate initial conditions, multiple (∼ 10) aligned DMDGs with a considerable amount of stars are formed from a single Mini-bullet collision.However, note that the alignment of the DMDGs is not exactly on a line.There are deviations in DMDG distribution from the line, at most ∼ 100 kpc at t = 2 Gyr and ∼ 300 kpc at t ∼ 8 Gyr.Estimating these deviations in simulation could provide valuable insight for future observations in identifying which UDG is on the sequence of DMDGs and shares the common origin with DF2 and DF4, potentially providing a way to disprove or substantiate the Mini-bullet scenario.
We present the spatial displacements of the DMDGs from the "trail" and the deviations in their line-of-sight velocities in Table 4, along with their positions at the end of the orbit integration (Figure 7).The displacement from the sequence of DMDGs, denoted as ∆d, is defined as the distance from the DMDG to the line y = Ax + y 0 that best fits all the DMDG positions.Line-of-sight velocities (v y ) of DMDGs other than simulated DF2 and DF4 candidates are estimated using either projected positions (x) or line-of-sight positions (y) through a simple linear relationship v y, x-pred = Ax + v y, 0 (v y, y-pred = Ay + v y, 0 ) derived from the positions and velocities of DF2 and DF4 candidates.Estimating velocities based on the linear relationship derived from x and y results in significant differences between the estimated line-of-sight velocities v y, x-pred and v y, y-pred , with v y, y-pred being more accurate.This is expected since the sequence of DMDGs is elongated along the y-axis and relative deviations in x are more significant than in y.
In Figure 8, we summarize the correlation between v y and x, and v y and y of the product DMDGs, including simulated DF2 and DF4 candidates and the progenitors, along with linear regression fits.We perform least-square regression linear fits to the correlations using the data of the product DMDGs shown in Table 4 (not including the data of the progenitors), overplotting the fits with dashed lines.While both correlations follow linear relationships fairly well, the correlation between v y and y is tighter than that of v y and x.The scatter in the correlation between v y and x is crucial in interpreting the trail of UDGs in the NGC1052 group.For instance, even if a UDG is located at the end of the trail when projected on the sky (i.e. with the greatest projected distance from NGC1052), a different UDG, located at a smaller projected distance, can have the most extreme line-of-sight velocity.
The implications of these findings for future observations are that determining the membership of a UDG in the NGC1052 group within the sequence of DMDGs is a challenging task.The ideal approach involves precise measurements of line-of-sight velocities, distances, and projected distances of the aligned UDGs.However, achieving high precision in line-of-sight distance measurements is extremely difficult even with deep imaging with a significant amount of telescope time.Thus, the first step would involve establishing a correlation between projected distances and line-ofsight velocities.
In a study by Gannon et al. (2023), the line-of-sight velocity of DF9 was measured and found to be inconsistent with the expected linear relation between projected distances and line-of-sight velocities, if it was a member of the trail along with DF2 and DF4.This caused the authors to question whether DF9 was correctly identified as a member of the trail or shared a common origin with the others in the trail, whether 3D geometry in projection might alter the interpretation relative to the simple linear relationship if it was a member, or, finally, whether the idea of a common origin for all the galaxies in the trail was even correct.As the simulations here of the Mini-bullet scenario predict, however, individual DMDGs are expected to exhibit deviations in their positions and velocities from that simple linear relation, so these observations of DF9 are not inconsistent with that, so far.By gathering data from future deep observations of mul-  tiple UDGs on the trail, it will be possible to observe more of the details of the galaxies in the apparent trail, to better compare with the Mini-bullet formation scenario of the aligned UDGs in the NGC1052 group.
Another interesting prediction we make is the correlation between the predicted deviation of line-of-sight velocity and transverse displacement of a DMDG from the best-fitting line of the common trail of DMDGs.This is shown in Figure 9, by defining the deviations as the difference between predicted v y using projected positions (x) and the actual v y (v y − v y, x-pred ).
In the TM1 run, 2 Gyr after the collision, the stellar mass of the most (second) massive DMDG is 2.6 × 10 8 M ⊙ (1.5 × 10 8 M ⊙ ), showing a good agreement to the observed stellar mass, ∼ 2.0 × 10 8 M ⊙ for DF2 and ∼ 1.8 × 10 8 for DF4 (van Dokkum et al. 2018b(van Dokkum et al. , 2019)).Moreover, their line-ofsight velocity difference at the end of the orbit integration (t = 8.375 Gyr for the TM1 run) is 324 km s −1 , which is similar to the observed value of 358 km s −1 (v DF2 = 315 km s −1 and v DF4 = −43 km s −1 relative to the host, NGC1052; van Dokkum et al. 2022a).These results, along with those from TM2 and TM3 runs, are summarized in Table 4.The stellar masses show little variation, but the line-of-sight velocity differences can vary several times depending on the initial orbital parameters of the progenitors and the subsequent orbits of the product DMDGs.This is because, as the orbit of the second massive DMDG (corresponding to DF4) passes closer to the host, it experiences stronger gravity, resulting in a more negative line-of-sight velocity.
van Dokkum et al. (2022a) postulated that NGC1052-DF7 (hereafter DF7) and RCP32, located at the farthest end of the UDG trail, represent the remnants of the progenitor galaxies.Notably, a follow-up study showed that the identified GCs near RCP32 exhibit stellar populations that are more akin to those of the host galaxy NGC1052, in contrast to the GCs associated with DF2 and DF4 (Buzzo et al. 2023).This finding provides support for the notion that RCP32 may be the preserved remnant of a post-collision satellite galaxy.In Figure 7, the trajectories of the progenitors in our simulations and orbit integration are illustrated.One of the progenitors is located at the farthest from the host, roughly on the trail of DMDGs, suggesting that the progenitor might correspond to RCP32.On the other hand, the other progenitor passes through and orbits the host halo rather than remaining in the sequence of the DMDGs and the other progenitor.This suggests the possibility that one of the progenitor galaxies may  Correlation between line-of-sight velocity and both line-of-sight and projected distance from the host galaxy.Line-of-sight velocity relative to the host galaxy (v y ) vs. projected distance (x) from the host galaxy (left column).Line-of-sight velocity relative to the host galaxy vs. line-of-sight distance (y) from the host galaxy (right column).Black dashed lines are linear regression fits of the correlations using the data of the product DMDGs excluding the progenitors.The progenitors are marked with squares, with Progenitor 1 (the one farther from the host) and Progenitor 2 (the one closer to the host).The correlation between v y and y is tighter than that of v y and x.See Table 4 for the exact numbers used to plot this figure.4 for the exact numbers used to plot this figure.
have experienced strong tidal forces exerted by the host, potentially rendering it challenging to be detected in the present day.

Properties of product DMDGs: stellar metallicities, ages, and sizes
Now that we have confirmed that the spatial distribution, line-of-sight velocities, and stellar masses of the product DMDGs in our simulation are consistent with those observed by van Dokkum et al. (2022a), we turn our attention to their detailed physical properties.Figure 10 presents the stellar masses, average metallicities, and average formation time of the stars of the formed DMDGs at t = 2 Gyr, when the hydrodynamic simulation ends.The two most massive DMDGs, which correspond to DF2 and DF4, are marked with red and magenta, respectively, to distinguish them from other less massive DMDGs.On the left panel, we show the stellar masses and the average metallicities with the error bars indicating the standard deviations of the metallicities.Overall, the average metallicities of the DMDGs exhibit a remarkable similarity with a deviation of [M/H] ∼ 0.2, regardless of their stellar mass.The right panel shows the metallicities and formation time of the stars, or the age of the stars at the end of the orbit integration on the x-axis above, which corresponds to the observed age at the present epoch.It is notable that massive DMDGs form stars for a longer period than other less massive DMDGs, along with the accretion of nearby gas.
Even though the initial metallicity of the gas starts with Z = 0.1 Z ⊙ , the metallicities range from a few Z ⊙ to ∼ 10 Z ⊙ , which is a lot higher than the observed metallicities of DF2 ([M/H] = −1.07± 0.12; Fensch et al. 2019).We argue that this discrepancy originates from two effects: pre-enrichment of the gas in the progenitor galaxies before the collision and the stellar feedback model we employ in our simulation.As described in Section 2.1, the stellar feedback is implemented with thermal feedback of supernova.However, simplifying various modes of stellar feedback, such as stellar winds and radiation feedback that affects nearby interstellar medium over greater distances, can result in overly efficient self-enrichment of the stars that form subsequently.This does not mean that the metal formation is overestimated but the distribution of metal is not efficient.Especially in the case of bursty star formation resulting from shock compression of gas due to high-velocity galaxy collision, metal distribution should be carefully modeled (for example, see Han et al. 2022, for the case of GC formation in cloud-cloud collisions).Including more realistic stellar feedback modes will be the subject of future studies, and our results should be interpreted as controlled simulations showing that the formed DMDGs have almost the same metallicities and ages, which can be tested in future observations.
Then in Figure 11, the temporal evolutions of the 3D stellar half-mass radii (R 1/2 ) of the DMDGs from the TM1, TM2, and TM3 runs are presented.The tracking of the DMDG is based on the star particle IDs at t = 2 Gyr, and the center of the DMDG is defined as the center of mass of the identified stars at that time, at which we measure R 1/2 .It is worth noting that under the assumption of spherical symmetry and density profile, R 1/2 can be converted to the effective radius, or 2D projected half-light radius, R eff , by R eff = R 1/2 /A, where A ∼ 1.3 − 1.35 (Wolf et al. 2010).The color-coded lines in Figure 11 indicate the 3D half-mass radii, with the distances from the host galaxy at t = 2 Gyr, r host,t=2 Gyr , determining the color scheme.Age ( t end t form ) (Gyr) 5.9 6.0 6.1 6.2 6.3 7.8 7.9 8.0 8.1 8.2   Despite the suite of hydrodynamic simulations, TM1, TM2, and TM3, start from nearly identical initial conditions, the sizes of the product DMDGs are notably different.This means that the DMDG sizes are highly sensitive to various complicated physical processes involved, especially the turbulent behavior of gas driven by stellar feedback processes during the DMDG formation.Since the simulations we perform are limited by the spatial resolution of 5 pc and simple stellar feedback physics adopted, the goal of realistically modeling the turbulent gas emerging from cloud scale and the influence of the tidal field on the distribution of gas and stars at the same time is extremely hard to achieve.Therefore, drawing a definitive conclusion on the sizes of DMDGs in Mini-bullet scenario from the simulations presented in this work is challenging.
Instead, we focus on the effect of simple physical processes including the tidal field of the host and the merging of stellar structures long after the initial burst of star formation at t ≲ 0.4 Gyr.In general, the larger and more distant the DMDG is from the host, the smaller its size becomes, and vice versa.These results suggest that the primary factor influencing the sizes of the formed DMDGs is their susceptibility to the tidal forces exerted by the host galaxy, which is determined by a combination of the distance from the host and the mass of the DMDG.The post-formation evolution of size is influenced by processes such as tidal interaction with the host and merging of star clumps.Removal of stars on the outskirts of a DMDG by tidal interaction results in a gradual decrease in size, while merging events are evident in Figure 11 as sudden changes in size9 .

Differences between the observation and simulation
We have so far demonstrated that simulations starting from initial conditions designed to lead to a satellite galaxy-galaxy collision whose products match the spatial and kinematic properties and stellar masses of the NGC1052 group UDGs can indeed produce the observed trail of UDGs and their alignment, along with the dark matter deficiency of DF2 and DF4.However, upon more detailed inspection, the simulated UDGs do not exactly reproduce all the observed characteristics of these objects.Indeed, this is not surprising, since some details of the simulation outcome are dependent on numerical limitations and choices, such as spatial and mass resolution and the subgrid prescription for star formation and feedback physics.In this section, we focus on how the simulated outcomes differ from observational results, especially the sizes and metallicities of the stellar components and GCs hosted by the produced galaxies.We also discuss what needs to be taken into account for the comparison beyond our simulations by discussing other previous work on this subject.

Size of DMDGs
As presented in Section 3.2 and shown in Figure 11, the half-mass radii of stellar components of the produced galaxies tend to be a few times smaller than the observed sizes of DF2 and DF4.One of the factors that affect their sizes is the tidal field of the host.Although the tidal force can reduce the size of a DMDG by removing stars from its outskirts, as appears in Figure 11 as a gradual decrease, heating of the stars in the central region of the galaxy can expand the DMDG (see, e.g., Gnedin et al. 1999;Gnedin 2003).For instance, in the context of transforming normal low-mass (dwarf) galaxies into UDGs by tidal heating, Jones et al. (2021) claimed to have found possible evidence of this process in UDGs NGC2708-Dw1 and NGC5631-Dw1.Also, while Liao et al. (2019);Tremmel et al. (2020); Liao et al. (2019) demonstrates this in cosmological simulations, Carleton et al. (2019) uses a semi-analytic model applied to dark matter-only simulation and shows the expansion of low-mass galaxies due to tidal heating.Product DMDGs that are close to the host including the simulated DF4 candidate will experience tidal heating even after t = 2 Gyr, the end time of our hydrodynamic simulation, limiting the assessment of the tidal heating effect in our work.Therefore, modeling of tidal heating of dark matter-less galaxies that occurs beyond our simulation is necessary to gauge stellar component size expansion.
On the other hand, apart from long-term size evolution on the orbits, strong supernova feedback from the formation of massive GCs during the initial star formation epoch (t ≲ 0.4 Gyr) can effectively expand stars in DMDGs (Trujillo-Gomez et al. 2022).This is particularly efficient in the case of a galaxy without a co-centered dark matter halo.Since the question of how powerful the feedback should be to expand a DMDG produced in the Mini-bullet scenario into UDG needs to be studied further quantitatively, a self-consistent high-resolution simulation that resolves massive star cluster formation and its feedback is necessary to understand the expansion process.

Metallicity of DMDGs
Another difference between the simulated DMDGs from observation is that the stellar metallicities are super-solar in the simulated galaxies, whilst the observation of DF2 revealed sub-solar metallicity of [M/H] ∼ −1.The overshooting of metallicity could also be observed in simulated DMDGs and star clusters in Paper II, mainly caused by the simple prescription there for supernova metal ejection and an over-efficient self-enrichment -the process that locks up metals in successive generations of star formation rather than dispersing them into the ISM-in massive star clusters.
However, it is well known that metal dispersal within a galaxy is dependent on the computational recipe for metal mixing (see e.g., Shin et al. 2021).Furthermore, Han et al. (2022) found in their radiation-hydrodynamic simulations that radiation feedback from GC stars can reduce the selfenrichment of second-generation stars by preventing immediate accretion of supernova-enriched gas and delaying subsequent star formation, which allows the released metals to mix sufficiently with interstellar gas before the next episode of star formation.The results of this work suggest that the inclusion of radiative feedback from massive stars is crucial in studying heavy element mixing in hydrodynamic simulations.To model the evolving chemical enrichment of gas and stars more realistically inside the star clusters formed in the galaxy-galaxy collisions simulated here, therefore, it will be necessary to add the radiative transfer of ionizing UV photons that photoionize H and He, and of optical/infrared photons that exert radiation pressure on the surrounding gas, mediated by dust.We will consider that in future work.4.1.3.Survivability of globular clusters in the host tidal field Ogiya et al. (2022b) argued that it is challenging for the Mini-bullet scenario to explain both the number and spatial extent of massive GCs in DF2.Considering dynamical friction, the GC distribution should have been more extended just after their formation than it is now.Thus, the GCs are more susceptible to the strong tidal force from NGC1052, demanding a larger number of GCs at the time of the formation of DF2, meaning that most of the stellar mass in DF2 should have been in GCs, which seems unlikely.
This argument relies on the assumption that the location of the Mini-bullet collision was inside the virial radius of NGC1052, however, in order to have sufficient strength of the tidal force.In fact, in the constrained initial conditions designed here to lead to a collision whose end-result matches the current positions and velocities of DF2 and DF4 after their orbital evolution for ∼ 8 Gyr, that collision should have taken place much further out, at ∼ 2× the virial radius of NGC1052 (see Section 2.2 and Figure 4).Otherwise, if the collision point is too close to NGC1052, the alignment of DMDGs does not occur, because product DMDGs close to the host will not be part of the sequence of DMDGs because they will be accelerated by gravity and end up on the other side of the host.Combining these, we claim that in the case of the NGC1052 group system, the Mini-bullet collision should have occurred far from the host galaxy, enabling both the formation of a trail of DMDGs and extended GC distribution in the most massive DMDGs.We support this claim with further quantitative analysis.Figure 12 displays the tidal radii R tidal , stellar masses M ⋆ , and the distances of the product DMDGs from the host at t = 2 Gyr, r host,t=2 Gyr , for each simulation run, TM1, TM2, and TM3.As a result of being farther from the host, the tidal radii R tidal of the simulated DF2 and DF4 candidates are larger than 20 kpc, which is > 2 times the current spatial extent of the DF2 GCs and large enough to retain extended GCs.4.2.Statistical likelihood of Mini-bullet satellite-satellite collision in a large-volume ΛCDM cosmological simulation Quantifying the probability of the Mini-bullet event is important, as such an event is expected to be rare, considering the high v col and the small r min of the collision.Moreover, as discussed in Section 4.1.3,the collision should occur outside of the virial radius of the host halo, possibly reducing the likelihood of its making the occurrence even further.To address this, we will now use a large-volume cosmological simulation of galaxy and large-scale structure formation in a ΛCDM universe to study the frequency of such Mini-bullet collisions in what follows 10 .
As discussed in the previous sections, in the Minibullet scenario, two low-mass satellite galaxies with M 200 ∼ 10 9−10 M ⊙ collide with v col ∼ 500 km s −1 and r min ≲ R s, gas .This condition is expected to be rare in the Universe, considering that the progenitor satellite galaxies collide outside the host galaxy's virial radius (Figure 4) with a relative velocity well in excess of that host's virial velocity but with such a small impact parameter that their centers approach within a small fraction of their own virial radii.Moreover, as described in Section 2.1, if one of the colliding satellites is unbound to the system and the other satellite is bound to the system, the likelihood of the collision will be lower.Therefore, quantifying the statistical likelihood of the Mini-bullet scenario is essential in studying how plausible the scenario is.
Towards this end, we investigate the frequency of occurrence of such satellite-satellite galaxy collision events in a large-volume cosmological N-body/hydrodynamics simulation of galaxy formation in ΛCDM, in a box 100 cMpc on a side, TNG100-1, from the simulation suite known collectively as ILLUSTRISTNG (Pillepich et al. 2018a,b; Nelson et al. 2019; Marinacci et al. 2018; Nelson et al. 2018;   10 After this work was completed, we learned of the recent paper Otaki & Mori (2023), which also considered the likelihood of Mini-bullet collisions a ΛCDM universe, by a semi-analytical approach.Their focus was on collisions between subhalos inside the virial radius of a host galaxy.As we have described here, however, the collisions between satellite dwarf galaxies responsible for producing DMDGs like those observed by van Dokkum in NGC1052 took place far outside the virial radius of the host galaxy NGC1052.Springel et al. 2018;Naiman et al. 2018).ILLUSTRISTNG is the recent set of cosmological simulations by the ILLUSTRIS project, based on the moving-mesh code AREPO (Springel 2010).The cosmological parameters adopted in the ILLUS-TRISTNG are from the Planck 2015 results: Ω Λ,0 = 0.6911, Ω m,0 = 0.3089, Ω b,0 = 0.0486, σ 8 = 0.8159, n s = 0.9667, and h = 0.6774 (Planck Collaboration et al. 2016).TNG100-1, which has the highest resolution among the other TNG100 simulations, is composed of the 1820 3 baryon particles with 1.4 × 10 6 M ⊙ and 1820 3 dark matter particles with 7.5 × 10 6 M ⊙ .The particle-level data are stored in 100 snapshots from z = 20.05 to z = 0.In the halo catalogs, the halos identified with the friends-of-friends (FoF) algorithm (Davis et al. 1985) and the subhalos identified with the SUBFIND algorithm (Springel et al. 2001;Dolag et al. 2009) are stored.The ILLUSTRISTNG team also provides merger tree data generated by the SUBLINK algorithm (Rodriguez-Gomez et al. 2015) that traces the descendants and the progenitors of the subhalos.
In Paper I, we searched for collision-induced DMDGs directly in TNG100-1 but noted that the numerical resolution of TNG100-1 was insufficient to form DMDGs. To see DMDGs form, that simulation would have required higher mass and length resolution in both the dark matter and baryonic gas, in order to follow the shock compression of gas leading to star formation and then track the motion of those stars after they formed.Thus, in this section, we take a different approach, aiming to find, not the collision-induced DMDGs themselves, but rather the number of pre-collision galaxy pairs in TNG100-1 that meet the conditions established here for being possible progenitors of the Mini-bullet satellitesatellite galaxy collisions that form DMDGs.A similar approach was also performed in (Paper I), with criteria applied to the two-galaxy pairs found in TNG100-1.By contrast with our previous study (Paper I), however, we focus here, instead, on 3-body host-satellite-satellite systems and examine galaxy collisions in more detail by studying the location of the collision point and the gravitational boundedness of the colliding galaxies, to capture systems like those we showed here can produce the aligned galaxies of the NGC1052 system (van Dokkum et al. 2022a).
The exact search criteria are as follows.Conditions (i)-(iv) given below are applied to find the collision events that could have formed DMDGs and are based on the conditions we found to be capable of producing collision-induced DMDGs in Section 2.1 and our previous study (Paper I).
• (i) The colliding galaxies collide v col > 300 km s −1 , (ii) the gas fraction of the galaxies, f gas = M gas /M 200 , is greater than 0.05, (iii) the initial mass ratio of the colliding galaxies satisfies 1/3 < M 1 /M 2 < 3, and (iv) the initial gas mass ratio satisfies 1/3 < M 1, gas /M 2, gas < 3. Right: collisions that occur outside of R 500 of the collision pair's host halo.The collision pairs satisfy the criteria for the relative collision velocity, pericentric distance, gas fraction, total mass ratio, gas mass ratio, and total mass of colliding galaxies and the host.See Section 4.2 for a detailed description of the criteria.
The following conditions (v) and (vi) are added to find the collision events similar to the collision pairs that are studied in Section 3.1.
• (v) The total mass of each colliding progenitor galaxy, M 1 and M 2 , satisfies 10 9 M ⊙ < M 1 , M 2 < 10 11 M ⊙ , and (vi) the collision is associated with a host halo and the total mass of the host halo is M 200, host > 10 11 M ⊙ .
The following conditions are newly applied in this paper.
• (vii) The distance of the collision from the host halo is R > R 200 or R > R 500 .The Mini-bullet satellite-satellite galaxy collision likely to produce the NGC1052 system should occur far from the host galaxy.
• (viii) r min < (R eff, 1 + R eff, 2 )/2, where R eff, 1 and R eff, 2 are the half-mass radii of the colliding galaxies.To complement the temporal resolution of the TNG100-1 simulation, r min are calculated using the Rebound orbit integration code.Starting from the last snapshot before the collision, we take halos that are within 500 kpc from the colliding pair and mass of M 200 > 10 9 M ⊙ , calculating the halos' orbits with an assumption of their being point masses.
• (ix) We select pairs in which at least one of the colliding galaxies is gravitationally unbound to the host halo, meaning the mechanical energy E mec = E K + E P > 0, where E K is the kinetic energy and E P the gravitational potential energy.E P is calculated assuming that the host M 200, host is at the center of mass of the FOF halo.
Figure 13 displays the distribution of the collision pairs that satisfy the search criteria (i) − (viii) in the parameter space of v col and r min occurring outside of R 200 or R 500 of the host halo.With all the search criteria (i) − (ix) applied to 95 halo catalogs from z = 3 (10) to 0.01, we find 11 (36) collision pairs happening outside of R 200 of the host halo and being unbound to the center of mass of the host halo 11 .z = 3 corresponds to the lookback time of 11.6 Gyr, roughly corresponds to the observed stellar age of DF2 and DF4 +2× (measurement error).On the other hand, z = 10 is employed as a starting point of the search range z = 10 to 0.01 to exclude too frequent galaxy collisions at the early epoch of galaxy formation z > 10.Other numbers for less strict conditions are summarized in Table 5.In conclusion, at z < 3, in the simulated ∼ (100 Mpc) 3 -sized Universe of TNG100-1, we attest that there are ∼ 10 Mini-bullet satellite-satellite galaxy collisions that are expected to produce a sequence of DMDGs even with the most conservative criteria.Since these numbers are subject to the accuracy of the orbit integration, we examine and discuss this in Appendix (Section A).

CONCLUSION AND SUMMARY
In this paper, using gravitohydrodynamic simulations starting from the carefully designed initial conditions and orbit integrations applied to the simulations, we have demonstrated that the Mini-bullet satellite-satellite galaxy collision scenario is capable of explaining the observed properties of the unusually aligned galaxies in the NGC1052 group.
Informed by results from idealized two-body hydrodynamic simulations (Section 2.1), we set the initial orbital and structural parameters of the collision progenitors, one being a bound satellite galaxy and the other being an unbound satellite galaxy of a massive host galaxy that corresponds to NGC1052 (Section 2.2).We simulate satellite-satellite galaxy collision in the three-body host-satellite-satellite system with a suite of high-resolution (5pc) ENZO simulations (2.3).The simulations are augmented with the orbit integration code Rebound to study the orbital evolution of dark matter-deficient galaxies (DMDGs) produced in the satellitesatellite galaxy collision for ∼ 8 Gyr (3.1).
Our main results are as follows: • We demonstrate that the Mini-bullet satellite-satellite galaxy collision that can reproduce the DMDGs in the NGC1052 group observed by van Dokkum et al.
(2022a) should have taken place outside the virial radius of the host galaxy NGC1052, i.e more than ∼ 2× the virial radius, in fact.As a result, the globular clusters (GCs) formed during this galaxy collision are able to survive tidal disruption by the host galaxy.(Section 2.2 and Section 4.1.3) • We show that the Mini-bullet satellite-satellite galaxy collision event is capable of producing a "trail of 11 Note that this number is smaller than the number we got in Paper I, 248.This is mainly due to the addition of conditions (v) − (ix).Counting the number of collisions that happen inside of host halo R 200 , this number will be much larger.
• We find that while the positions (x x x) and velocities (v v v) of the aligned DMDGs approximately conform to a linear relationship, v v v = Ax x x + v v v 0 , individual DMDGs can significantly deviate from that simple relation (Section 3.1).
• We forecast that one of the collision progenitor galaxies is likely to be located at the end of the trail of DMDGs and might be able to be confirmed to be distinct from other aligned UDGs in future observations (Section 3.1).
• We investigate the stellar ages and metallicities of the product DMDGs and find that they are nearly identical within the standard deviations (Section 3.2).
• We study the size evolution of the product DMDGs and conclude that the simulated DMDGs are a few times smaller than the observed size of DF2 and DF4.The tidal field of the host galaxy plays a major role in making a DMDG larger, but further study is needed to capture complicated physical processes during the formation of DMDGs (Section 3.2).
We also discuss differences between the observations and our simulation results, concluding they are not a significant challenge to the mini-Bullet scenario (Section 4.1).The observed DMDGs, for example, are somewhat larger and more diffuse than the simulated ones.However, this can be explained as a possible outcome of tidal interaction with the host galaxy and nearby galaxies during the formation of DMDGs, which our simulations did not take into account.With higher numerical resolution, the simulations would also have resolved the formation and feedback of massive star clusters, which would serve to puff the DMDGs up (Section 4.1.1).Inclusion of this enhanced stellar feedback would also serve to expel some of the metals before they are recycled into new stars, thereby reducing the metallicity of the simulated DMDGs, found here to exceed the observed values somewhat (Section 4.1.2).The massive GCs observed in the DMDGs DF2 and DF4 suggest that such massive star cluster formation did take place.
To quantify the statistical likelihood of the Mini-bullet satellite-satellite galaxy collision event in the Universe, we inspect a large cosmological simulation TNG100-1.We confirm that 11 galaxy collision pairs at 0.01 < z < 3 satisfy the search criteria which we carefully choose to match the characteristics of the initial conditions that resulted in the formation of multiple aligned DMDGs similar to the NGC1052 group galaxies in the hydrodynamic simulations (Section 4.2).

APPENDIX A. ACCURACY OF ESTIMATING PERICENTRIC DISTANCE WITH ORBIT INTEGRATION CODE IN TNG100-1
As briefly mentioned in Section 4.2, we use the Rebound orbit integration code to compute r min of the colliding galaxy pairs in TNG100-1.To improve the accuracy of the orbit integration, we take into account the gravitational force of not only the pair of colliding galaxies but also other surrounding galaxies together if they are within 500 kpc from the colliding pair and have the total mass of M 200 > 10 9 M ⊙ .Given that the galaxies are approaching with a high velocity (≳ 300 km s −1 ), the orbit is computed for 200 Myr, within which the pericentric approach (collision) happens.We assume that the galaxies behave as point mass particles and ignore the impact of the Hubble flow on the orbits of the galaxies because the orbits are calculated for a very short time (200 Myr).
Since we are interested in the distances between the colliding galaxies, not the actual position, we compare the distances between the colliding galaxies in TNG100-1 snapshots with the distances calculated using the Rebound code to investigate the accuracy of the orbit integration.In TNG100-1, there are only 100 snapshots so we cannot take the continuous trajectories.We choose the snapshot right before the pericentric approach and compare the distances between the colliding galaxies at that moment.We study 1804 collision events in the redshift range of 0.01 < z < 10 (the first row of Table 5), in which galaxies collide with r min < 15 kpc and outside of R 500 of the host halo.In 868 events among the 1804 collision pairs, galaxy collision happens earlier than the output time step, so we use the remaining 936 collision events for the comparison.
The results of the comparison are shown in Figure 14.The distances between colliding galaxies before the collision in Rebound orbit integration and NG100-1 are presented.We can see that the majority of points are located slightly above the y = x line.This means that the predicted distances of colliding galaxies in orbit integration calculation are smaller than the distances in TNG100-1.This is due to the point mass assumption of galaxies.Also, some points are located far away from the y = x line.The main reason for this is the intervention of surrounding galaxies: if one of the surrounding galaxies approaches one of the colliding galaxies very closely, it can disturb the orbit of the colliding galaxy.However, the prediction from the Rebound orbit integration is generally consistent with TNG100-1, justifying the accuracy of the orbit integration.Thus, we argue that the existence of ∼ 10 Mini-bullet satellite-satellite galaxy collisions Figure 14.Scatter plot of the comparison of the distances between colliding galaxies in Rebound orbit integration code and TNG100-1.From TNG100-1, the snapshot right before the pericentric approach is taken, and the distances are compared at that moment from Rebound orbit integration calculation.The dashed red line is the y = x line.
at 0.01 < z < 3, in a ∼ 100 3 Mpc volume, satisfying the search criteria (i) − (ix) described in Section 4.2, is robust enough that we can believe the plausibility of the Mini-bullet event in the Universe.

Figure 2 .
Figure 2. Illustrative Mini-Bullet Simulation: collision in-progress (collision-plane view).Same as Figure 1 but zoomed-in in time and space to show snapshots of the galaxies in the midst of their collision, t = −30, 0, 30, 60 Myr.The surface density of dark matter (first row), surface density of gas with velocity vectors overplotted (second row), density-weighted average gas temperature (third row), and surface density stars -just those formed after the start of the simulation -(fourth row) are presented.Gas velocity vectors are projected onto the image plane, which is the xy plane to which the centers of mass of the colliding-galaxy orbits are confined.The velocity vectors are proportional to the length of the velocity vector arrows and the largest arrow corresponds to ∼300 km s −1 .Black arrows in the top-left panel show the progenitors' moving directions.Dark matter and stars are projected in 50 kpc-depth layers and gas is projected in 10 kpc-depth layers.

Figure 3 .
Figure3.Illustrative Mini-Bullet Simulation: collision in-progress (viewed perpendicular to collision-plane).Same simulation snapshots as Figure2, at t = −30, 0, 30, 60 Myr, but viewed in projection along the x-axis in the yz-plane -i.e. the view close to the line of collision.The surface density of the gas (top row) and density-weighted average gas temperature with velocity vectors overplotted (bottom row) are presented.Gas velocity vectors are projected onto the image plane, which is the yz plane, in which the center of the mass of the colliding progenitor system was at t = 0 Myr.The velocity vectors are proportional to the length of the velocity vector arrows and the largest arrow in each panel corresponds to ∼200 km s −1 .Gas is projected in 10 kpc-depth layers.

Figure 4 .
Figure 4.Backward Orbit Integration to Produce a "Tailor-Made" Initial Condition.An example of backtraced orbits of DF2 (solid red) and DF4 (solid blue), computed orbits of the progenitors (dashed red and blue) that are expected to produce a trail of DMDGs including DF2 and DF4 in a Mini-bullet collision, and their expected orbits after the collision.The purple cross indicates the location of the Mini-bullet collision.The black circle is the virial radius (R 200 = 400 kpc) of the host, NGC1052.This example corresponds to TM3 run in Table2See Section 2.2 for more information.

Figure 5 .
Figure5.Simulation of Mini-Bullet Satellite-Satellite Galaxy Collision from "Tailor-Made" Initial Conditions: Making the DMDGs/UDGs near NGC1052.Snapshots from the ENZO simulation of two colliding satellite galaxies of NGC1052 from TM1 run ("Tailor-Made" fiducial run) initial conditions at t = −0.1 Gyr, t = 0 Gyr, t = 1 Gyr, and t = 2 Gyr, centered on the two colliding satellite progenitor galaxies in the reference frame of their center of mass.t = 0 is set to the moment when the two progenitors are at a pericentric approach.Surface densities of dark matter (top row), gas (middle row), and stars (bottom row) in 100 kpc-depth layers are presented.Black arrows in the top-left panel indicate the approach of the progenitor galaxies toward each other before they collide.In the top row, white arrows on each plot point to the direction of the host galaxy NGC1052, centered at the distances labeled, and the length scale in each panel is indicated by the distance rulers of 40 kpc overplotted there.A few hundred Myr after the collision, almost no gas is left.At t = 1 Gyr and t = 2 Gyr, a sequence of DMDGs can be observed.

Figure 6 .
Figure 6.Simulation of Mini-Bullet Satellite-Satellite Galaxy Collision from "Tailor-Made" Initial Conditions: Segregating Dark Matter, Stars, and Gas.Same as Figure 5, t = −0.1 Gyr, t = 0 Gyr, t = 1 Gyr, and t = 2 Gyr snapshots, but shows stars overplotted on dark matter particles and gas.Dark matter particles are plotted as blue (left column).Star particles formed after the simulation starts are overplotted red (left column) and white (right column).Black arrows in the top-left panel indicate the approach of the progenitor galaxies toward each other before they collide.Distance rulers of 80 kpc are overplotted on the left column.Gas is projected in 100 kpc-depth layers with a color bar showing gas surface density (Σ gas ; right column).

Figure 7 .
Figure 7. Orbits the progenitor galaxies and the product DMDGs.Solid lines indicate the orbits in the ENZO simulations and dashed lines are computed orbits from the Rebound orbit integration code after the simulation end-time.The thick blue lines denote the trajectories of the progenitors (black galaxy icons).The thick red and purple lines are the orbits of the two most massive DMDGs, which correspond to DF2 and DF4, respectively (realistic galaxy icons).The orbits of other DMDGs are color-coded by their stellar mass indicated by the color bar on the right (small blue galaxy icons).The black circle indicates the virial radius of the host galaxy, NGC1052 (R 200 = 400 kpc and M 200 = 6.124 × 10 12 M ⊙ ).Left: Orbits in the TM1 (Tailor-made) run for t = −0.3− 8.375 Gyr.Middle: Orbits in the TM1 run for t = −0.3− 6.65 Gyr.Right: Orbits in the TM1 run for t = −0.3− 8.502 Gyr.The end time of orbit integration is when simulated DF2 and DF4 candidates are 2.1 Mpc apart in line-of-sight direction (y-axis).See Section 3.1 and Table3for more detail.
NOTE-(1) run name, (2) projected positions (x) and line-of-sight positions (y) relative to the host, line-of-sight velocities relative to the host (v y ), transverse displacement of a DMDG from the best-fitting line of the common trail of DMDGs (∆d; − signs indicate if the y coordinate of a DMDG is above (+) or below (−) the sequence in xy-plane), differences of predicted velocities v y, x-pred (v y, y-pred ) from a simple linear relationship between x (y) positions and actual velocities v y of simulated DF2 and DF4 candidates v y − v y, x-pred (v y − v y, y-pred ),(3) corresponding quantities of the simulated DF2 candidate, (4) those of the simulated DF4 candidate, (5) those of other aligned DMDGs.

Figure 8 .
Figure8.Correlation between line-of-sight velocity and both line-of-sight and projected distance from the host galaxy.Line-of-sight velocity relative to the host galaxy (v y ) vs. projected distance (x) from the host galaxy (left column).Line-of-sight velocity relative to the host galaxy vs. line-of-sight distance (y) from the host galaxy (right column).Black dashed lines are linear regression fits of the correlations using the data of the product DMDGs excluding the progenitors.The progenitors are marked with squares, with Progenitor 1 (the one farther from the host) and Progenitor 2 (the one closer to the host).The correlation between v y and y is tighter than that of v y and x.See Table4for the exact numbers used to plot this figure.

LEEFigure 9 .
Figure9.Correlation between the predicted deviation of line-of-sight velocity and displacement of a DMDG from the trail of DMDGs.Deviation of predicted line-of-sight velocity (v y ) based on the linear relationship between v y and projected distance (x) from and its simulated value (v y − v y, x-pred ) vs. transverse displacement from the line of DMDGs (∆d).Signs in ∆d indicate if the y coordinate of a DMDG is above (+) or below (−) the sequence of the DMDGs fitted with a simple linear relationship between x and y.See Table4for the exact numbers used to plot this figure.

Figure 10 .
Figure 10.Stellar Mass and Metallicity of the Collision-Produced DMDGs.Metallicity vs. stellar mass (left column).Metallicity vs. formation time (and age) (right column).Error bars correspond to the standard deviation of each physical quantity of the DMDG stars.Simulated DF2 and DF4 candidates are indicated with red and magenta dots respectively.

Figure 11 .
Figure 11.Stellar System Sizes of the Collision-Produced DMDGs vs. UDG Size.Time evolution of three-dimensional half-mass radius R 1/2 of the product DMDGs.After t = 1 Gyr, the amount of star formation of the product DMDGs is almost negligible.DMDGs with M ⋆ > 10 7 M ⊙ are color-coded by their distance to the host galaxy at t = 2 Gyr, the end of the simulations.Simulated DF2 and DF4 candidates are indicated with triangle markers.The observed stellar mass and metallicity of DF2 are marked as a red triangle on the left panels.Gray shade areas indicate the UDG criterion, R eff ≥ 1.5 kpc (R 1/2 = AR eff = 1.95 kpc assuming A = 1.3).The two most massive DMDGs, denoted by larger markers (representing simulated DF2 and DF4 candidates), exhibit up to a factor of ∼ 10 times more compact sizes compared to the observed R eff of DF2 (2.2 kpc; van Dokkum et al. 2018b) and DF4 (1.6 kpc; van Dokkum et al. 2019), except the simulated DF4 candidate in TM1 run, which has R 1/2 ∼ 3 kpc.The fact that this DMDG has the size of a UDG is a novel finding and enhances the plausibility of the Mini-bullet scenario.Prior to this work, in our previous works (Paper I; Paper II), it has not been demonstrated that the produce can be diffuse.The size difference in simulation runs is not negligible; the physical origin of the difference and the factors that affect the sizes of DMDGs will be discussed in detail in Section 4.1.1.Despite the suite of hydrodynamic simulations, TM1, TM2, and TM3, start from nearly identical initial conditions, the sizes of the product DMDGs are notably different.This means that the DMDG sizes are highly sensitive to various complicated physical processes involved, especially the turbulent behavior of gas driven by stellar feedback processes during the DMDG formation.Since the simulations we perform are limited by the spatial resolution of 5 pc and simple stellar feedback physics adopted, the goal of realistically modeling the turbulent gas emerging from cloud scale and the influence of the tidal field on the distribution of gas and stars at the same time is extremely hard to achieve.Therefore, drawing a definitive conclusion on the sizes of DMDGs in Mini-bullet scenario from the simulations presented in this work is challenging.Instead, we focus on the effect of simple physical processes including the tidal field of the host and the merging of stellar structures long after the initial burst of star formation at t ≲ 0.4 Gyr.In general, the larger and more distant

Figure 12 .
Figure 12.Star-Cluster Retention Radii of the Collision-Produced DMDGs.The tidal radii R tidal of the product DMDGs and the stellar masses M ⋆ of the DMDGs (left).R tidal and the distances from the host at t = 2 Gyr, r host,t=2 Gyr (right).Simulated DF2 and DF4 candidates are indicated with red and magenta triangles respectively.Other DMDGs with M ⋆ > 10 7 M ⊙ are color-coded by their stellar masses as presented in the color bar on the right.Top: 9 DMDGs produced in the TM1 run.Middle: 11 DMDGs produced in the TM2 run.Botttom: 9 DMDGs produced in the TM3 run.

Figure 13 .
Figure 13.Mini-bullet Progenitor Pairs With NGC1052-like host in a (100 cMpc) 3 Cosmological Simulation of ΛCDM.2D histogram of the collision pairs found in the TNG100-1 simulation volume.Left: collisions that occur outside of R 200 of the host halo of the collision pair.Right: collisions that occur outside of R 500 of the collision pair's host halo.The collision pairs satisfy the criteria for the relative collision velocity, pericentric distance, gas fraction, total mass ratio, gas mass ratio, and total mass of colliding galaxies and the host.See Section 4.2 for a detailed description of the criteria.

Table 1 .
A suite of 10 pc-resolution idealized dwarf galaxy-dwarf galaxy collision pair simulations listed with their initial configurations.

Table 3 .
Results from a suite of simulations of satellite-satellite galaxy collisions around a massive host.See Section 3.1 for a description and discussion of the properties of the product DMDGs.Run nameM ⋆, DMDG, max2 ∆v DMDG, max2 ∆v DMDG, max2, now N DMDG t end, sim t end, orbit

Table 2
See Section 2.2 for more information.

Table 4 .
Positions and line-of-sight velocities of DMDGs at the end time of orbit integration, presented in Figure7from TM1, TM2, and TM3 runs.See Section 3.1 for details.

Table 5 .
The number of colliding satellite-satellite galaxy pairs from z = 3 to 0.01 found in TNG100-1 that satisfy the search criteria.The numbers inside parentheses are the number of pairs from z = 10 to 0.01.See Section 4.2 for details.