Open Access
Issue
A&A
Volume 689, September 2024
Article Number A256
Number of page(s) 14
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202450346
Published online 17 September 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

Observations show that galaxy interactions and mergers play a significant role during the evolution of galaxies (Struck 1999; Smith et al. 2010). In addition to increasing the star formation rate (Kennicutt 1998; Georgakakis et al. 2000; Rossa et al. 2007), there are indications that mergers have an impact on the magnetic field strength in galaxies (Condon et al. 2002; Chyży & Beck 2004; Drzazga et al. 2011), suggesting an amplification of the field due to dynamo action. While the α-Ω dynamo induced by the differential rotation of disk galaxies is thought to play an important role in the generation of galactic magnetic fields from a weak primordial field (Brandenburg & Subramanian 2005; Park et al. 2013; Steinwandel et al. 2020; Wissing & Shen 2023), the α2 dynamo (caused by inhomogeneous shear flows with kinetic helicity) and the turbulent small-scale dynamo (induced by non-helical, isotropic turbulence on small scales) are expected to contribute significantly to the field amplification in merging galaxies (Arshakian et al. 2009; Schober et al. 2013; Bhat & Subramanian 2014; Moss et al. 2014). In contrast to isolated disk galaxies, where turbulence is assumed to be driven by supernova explosions (Hopkins et al. 2013) or gravitational instabilities (Krumholz & Burkhart 2016), shear and tidal effects can produce additional turbulence during galaxy interactions (Renaud et al. 2014; Whittingham et al. 2021; Jiménez & Lagos 2023).

The interaction process has been extensively studied in idealized numerical simulations without cosmological initial conditions. The basic idea is to embed disks composed of gas and possibly stars in spherical dark-matter halos and to set them on a collision course. An original study of such systems using a simple particle-based model was presented by Toomre & Toomre (1972) five decades ago. More recently, smoothed particle hydrodynamic simulations (Mihos & Hernquist 1996; Di Matteo et al. 2008; Blumenthal & Barnes 2018) and semi-Lagrangian simulations (Hayward et al. 2014; Moreno et al. 2019) focused on the dynamics of the interstellar medium and interaction-induced star formation. Merger simulations with magnetic fields were first performed by Kotarba et al. (2010), using a magnetohydrodynamic (MHD) extension of smoothed particle hydrodynamics (Dolag & Stasyszyn 2009). Their simulation of the Antennae galaxies suggests that mergers play a significant role in the amplification of galactic magnetic fields. Magnetic field amplification during minor mergers, where the secondary galaxy has only a small fraction of the mass of the primary galaxy, was investigated with a similar approach by Geng et al. (2012) and in cosmological simulations with the semi-Lagrangian Arepo code by Pakmor et al. (2014). Morphological changes due to mergers were studied for a sample of galaxies from the Horizon-AGN cosmological simulation (Martin et al. 2018). Moreover, the role of magnetic fields and the small-scale dynamo was analyzed in great detail in zoom-in simulations of major mergers from the Illustris simulation by Whittingham et al. (2021). However, in contrast to simulating pairs of galaxies, the resolution of cosmological simulations is limited and interaction parameters cannot be varied systematically.

While Lagrangian methods undoubtedly have merits for the treatment of moving galaxies, Rodenbeck & Schleicher (2016) demonstrate that merger simulations are feasible in a greatly simplified setup with purely gaseous and adiabatic disks without dark-matter halos using the MHD adaptive mesh refinement (AMR) code Enzo (Bryan et al. 2014). They applied the method from Wang et al. (2010) to initialize an isolated equilibrium disk in a static external potential that represents the gravity of dark matter with additional support from magnetic pressure. However, they excluded the external potential because it is not applicable to moving disks in a galaxy merger scenario, although the contribution of dark matter is dominant. Even so, this approach avoids the difficulties of implementing MHD into particle-based codes (Dolag & Stasyszyn 2009) and, owing to a well-defined resolution scale, offers the advantage of better controllability of the small-scale physics. For realistic simulations of interacting galaxies, however, N-body dynamics is essential. Hydrodynamic AMR simulations of a major merger similar to the Antennae system with live N-body halos were performed by Teyssier et al. (2010) and Renaud et al. (2015). At spatial resolutions in the parsec range, they were able to analyze the production of turbulence in the interstellar medium due to tidal interactions and the compression-induced enhancement of star formation.

In this paper, we combine the approaches of Teyssier et al. (2010) and Rodenbeck & Schleicher (2016). To generate initial conditions for spherical dark-matter halos, we follow Drakos et al. (2017). The gaseous disks are based on Wang et al. (2010). In essence, the gas density and rotation velocity of a disk in hydrostatic equilibrium are determined iteratively from the self-gravity of the gas and the smoothed potential of the dark-matter halo. The advantage of this method of initializing disks is that the effects of the initial relaxation are reduced and changes in the disk structure are mainly due to physical effects, such as tidal forces and fragmentation induced by cooling and gravity. The applied disk and halo models are described in more detail in the following section. In Sect. 3, the numerical setup is outlined and an overview of our simulation suite is provided. In this study we focus on magnetic field amplification due to tidal interactions between the disks and their halos. While radiative cooling and stellar feedback are essential ingredients for a realistic disk structure, the findings of Ntormousi et al. (2020) indicate that these processes are of minor importance for the generation of galactic magnetic fields. Since nonadiabatic physics increases computational requirements substantially for the targeted spatial resolution of about 10 pc, we use adiabatic gas dynamics in this work. The analysis of the simulation is presented in Sect. 4: After a phenomenological discussion of representative examples, the average magnetic fields in central and outer disk regions from different simulations and subsequent encounters of the interacting galaxies in each simulation are compared. Moreover, the question of whether magnetic field amplification during interaction phases can be related to dynamo action on the smallest numerically resolved length scales is investigated. In the conclusions, we sum up and discuss our results.

2. Disk and halo models

To initialize the gaseous disks, we applied the potential method of Wang et al. (2010) for a disk with exponential surface density:

Σ ( r ) : = ρ ( r , z ) d z = Σ 0 exp ( r / r d ) , $$ \begin{aligned} \Sigma (r) := \int _\infty ^\infty \rho (r,z)\mathrm{d} z = \Sigma _0\exp (-r/r_\mathrm{d} ), \end{aligned} $$(1)

where r and z are, respectively, the radius and height in a cylindrical coordinate system with the disk midplane at z = 0. The surface density at the center is given by Σ 0 = M d /2π r d 2 $ \Sigma_0 = M_\mathrm{d}/2\pi {r}_\mathrm{d}^2 $ for a disk of total mass Md and scale length rd. In many isolated disk simulations, an exponential form is assumed for the volume density ρ(r, z) of the gas. However, as argued by Wang et al. (2010), this implies a ring in the surface density.

The gas density is obtained by integrating the equation of hydrostatic equilibrium in vertical direction. For an initially isothermal disk composed of an ideal gas and a magnetic field with a constant ratio of thermal to magnetic pressure, β = p/pm, we have the total pressure

p + p m = ( γ 1 ) ρ e + B 2 8 π = ( 1 + ϵ m ) ρ c 0 2 , $$ \begin{aligned} p + p_{\mathrm{m} } = (\gamma -1)\rho e + \frac{B^2}{8\uppi } = \left(1 + \epsilon _\mathrm{m} \right)\rho c_0^2, \end{aligned} $$(2)

where e is the specific internal energy of gas with adiabatic index γ, ϵm = 1/β vanishes in the hydrodynamic limit, and the isothermal speed of sound is defined by

c 0 2 : = p ρ = k B T μ m H , $$ \begin{aligned} c_0^2 := \frac{p}{\rho } = \frac{k_\mathrm{B} T}{\mu m_\mathrm{H} }, \end{aligned} $$(3)

with Boltzmann’s constant kB, the gas temperature T, the mean atomic weight μ, and the hydrogen mass mH. In this case, the gas density can be expressed as (see Wang et al. 2010)

ρ ( r , z ) = ρ 0 ( r ) exp ( Φ z ( r , z ) ( 1 + ϵ m ) c 0 2 ) , $$ \begin{aligned} \rho (r,z) = \rho _{0}(r)\exp {\left(-\frac{\Phi _z(r,z)}{(1 + \epsilon _\mathrm{m} )c_0^2}\right)}, \end{aligned} $$(4)

where Φz(r, z) = Φ(r, z)−Φ(r, 0) is the vertical potential difference of the total gravitational potential

Φ ( r , z ) = Φ dm ( r , z ) + Φ s ( r , z ) + Φ g ( r , z ) $$ \begin{aligned} \Phi (r,z) = \Phi _{\mathrm{dm} }(r,z) + \Phi _{\mathrm{s} }(r,z) + \Phi _{\mathrm{g} }(r,z) \end{aligned} $$(5)

at height z relative to the midplane. The dark matter potential Φdm(r, z) follows from the halo model (see below) and the stellar potential Φs(r, z) is zero in this work.

From the Poisson equation, it follows that the gravitational potential of the gas, Φg(r, z), is given by

d 2 Φ g ( r , z ) d z 2 = 4 π G ρ 0 ( r ) exp ( Φ z ( r , z ) ( 1 + ϵ m ) c 0 2 ) . $$ \begin{aligned} \frac{\mathrm{d}^2\Phi _{\mathrm{g} }(r,z)}{\mathrm{d} z^2} = 4\uppi G\rho _0(r)\exp {\left(-\frac{\Phi _z(r,z)}{(1 + \epsilon _\mathrm{m} )c_0^2}\right)}. \end{aligned} $$(6)

Equation (1) is satisfied if the midplane density ρ0(r) = ρ(r, 0) is given by

ρ 0 ( r ) = Σ ( r ) exp ( Φ z ( r , z ) / [ ( 1 + ϵ m ) c 0 2 ] ) d z · $$ \begin{aligned} \rho _0(r) = \frac{\Sigma (r)}{ \int _\infty ^\infty \exp {\left(-\Phi _z(r,z)/ [(1 + \epsilon _\mathrm{m} )c_0^2]\right)}\mathrm{d} z}\cdot \end{aligned} $$(7)

The integral in the denominator can be interpreted as a disk scale height that varies with radius and is determined by the model. Equations (6) and (7) cannot be solved analytically. A numerical solution is readily obtained by using an iterative scheme in combination with a numerical solver for the differential equation (see Sect. 3). An overview of all model parameters is provided in Table 1.

Table 1.

Overview of the initial disk and halo model parameters.

For the dark-matter halo, we assumed the radial density profile (Hernquist 1990)

ρ dm ( r ) = M dm 2 π a r ( r + a ) 3 , $$ \begin{aligned} \rho _{\mathrm{dm} }(r) = \frac{M_{\mathrm{dm} }}{2\uppi }\frac{a}{r(r+a)^{3}}, \end{aligned} $$(8)

where a is the scale length of the halo and Mdm the total cumulative mass for r → ∞. Here, r is a spherical radial coordinate. As proposed by Springel et al. (2005), the scale length of the Hernquist halo is chosen such that it contains the same mass M200 within the virial radius r200 as a Navarro-Frenk-White halo of given scale length rs and concentration parameter c=r200/rs:

a = r s 2 [ ln ( 1 + c ) c / ( 1 + c ) ] . $$ \begin{aligned} a = r_{\mathrm{s} }\sqrt{2[\ln {(1+c)} - c/(1+c)]}. \end{aligned} $$(9)

The analytic form of the gravitational potential is

Φ dm ( r ) = G M dm r + a · $$ \begin{aligned} \Phi _{\mathrm{dm} }(r) = -\frac{G M_{\mathrm{dm} }}{r + a}\cdot \end{aligned} $$(10)

We used the method developed by Drakos et al. (2017) to compute random initial positions and velocities of the particles. The algorithm populates the phase space of an N-body system based on the distribution function of a halo in virial equilibrium.

The rotation curve of the disk is given by Wang et al. (2010), Rodenbeck & Schleicher (2016)

v rot 2 ( r ) = v dm 2 ( r ) + v thin 2 ( r ) + ( 1 + ϵ m ) c 0 2 ln ρ ln r | z = 0 , $$ \begin{aligned} v_{\mathrm{rot} }^{2}(r) = v_{\mathrm{dm} }^2(r) + v_{\mathrm{thin} }^{2}(r) + (1+\epsilon _{\mathrm{m} })c_0^2\frac{\partial \ln {\rho }}{\partial \ln {r}}\bigg \vert _{\mathrm{z=0} }, \end{aligned} $$(11)

where

v dm ( r ) = G M dm r r + a $$ \begin{aligned} v_{\mathrm{dm} }(r) = \frac{\sqrt{G M_{\mathrm{dm} } r}}{r + a} \end{aligned} $$(12)

for a Hernquist halo, vthin(r) is the contribution from an infinitesimally thin disk with exponential surface density defined by Eq. (1) (for details, see Binney & Tremaine 2008), and the last term is a correction due to the pressure gradient of the gas following from Eq. (4).

Finally, the initial magnetic field is assumed to be toroidal with constant β inside the disk (Rodenbeck & Schleicher 2016):

B ϕ = 8 π ρ c 0 2 β = 8 π ϵ m ρ c 0 2 , $$ \begin{aligned} B_\phi = \sqrt{\frac{8\uppi \rho c_0^2}{\beta }} = \sqrt{8\uppi \epsilon _{\mathrm{m} }\rho c_0^2}, \end{aligned} $$(13)

and Br = Bz = 0 in the disk’s cylindrical coordinate system. The disk interior is defined by ρ ≥ ρmin, where ρmin = 6.768 × 10−30 g cm−3 is the constant density of the ambient medium. While we accounted for magnetic pressure in the computation of the equilibrium density and rotation curve, a toroidal magnetic field gives rise to an additional curvature force that is directed radially inward. In the case of a constant plasma beta, this force gives rise to a constant negative contribution to the rotation curve. We neglected this contribution because it is lower by several orders of magnitude than the maximum rotation speed for the chosen disk parameters.

The combination of a rotating disk in hydrostatic equilibrium with an N-body halo and a constant background gas density introduces several effects impeding the stability of the disk:

  • Adiabatic disks have relatively large scale heights, particularly in the outer regions, while the model of Wang et al. (2010) is based on a thin-disk approximation. Near the center, the scale height is small, but the disk structure cannot be mapped accurately to the grid due to the limited numerical resolution. In addition, the rotation velocity is independent of the height z as long as the temperature remains constant, but adiabatically evolving perturbations give rise to z-dependent deviations from the initial rotation curve.

  • The density of the gaseous disk given by Eq. (4) decreases rapidly with the scale height. To avoid unphysically low densities and density jumps at the interface between disk and ambient medium, the density floor ρmin mentioned above is applied when setting up the disk. Since a constant density is not consistent with the equilibrium solution, gas is subsequently falling toward the disk. This could be alleviated by modeling a circumgalactic medium (see, for example, Steinwandel et al. 2019). However, the total mass in the ambient medium is very small compared to the disk mass.

  • Since the gravitational potential is determined by discrete particles, changes in the particle positions entail fluctuations of the potential, which in turn affect the gas.

  • Although the gaseous disk is computed from the joint potential of halo and disk, the distribution of dark matter particles is computed assuming the analytical potential of the halo and neglecting the gravity of the disk.

As a consequence, perturbations of the initial equilibrium disk are unavoidable. For the most part, the disk does not undergo abrupt changes, but gradually (over a few hundred megayears) evolves into a disk that exhibits a spiral-like structure.

3. Numerical methods and overview of simulations

We used a modified version of the publicly available cosmological AMR code Enzo (Bryan et al. 2014; Brummel-Smith et al. 2019)1. The code is MPI-parallelized, features N-body dynamics based on a second-order drift-kick-drift algorithm in combination with cloud-in-cell interpolation to compute the joint gravitational potential of gas and particles, and a variety of different finite volume solvers for gas dynamics and magnetic fields. In our simulations, we applied the monotonic upstream-centered scheme for conservation laws with a local Lax-Friedrichs solver for sufficient numerical robustness. Adiabatic gas dynamics is applied under the assumption that magnetic field dynamics during interactions is mainly driven by tidal forces. Since radiative cooling, star formation, and feedback processes are not included, the spatial resolution and time stepping requirements are significantly reduced. Even so, the required computational resources for an extensive parameter study were substantial, as detailed below.

Originally, the computation of initial conditions for a disk in hydrostatic equilibrium were part of the initialization in Enzo (Rodenbeck & Schleicher 2016). However, this was not practical for the particle initial conditions required for the live dark-matter halo. For this reason, we implemented a Python tool for computing particle initial conditions along the lines of ICICLE (Drakos et al. 2017). By applying object-oriented programming, we were readily able to incorporate equilibrium disks as a subclass of the halo class into our Python code. The algorithms outlined in Sect. 2 are implemented as object methods and make use of the highly optimized numerical methods for root finding and integration from the SciPy library. The resulting disk and particle data can be saved to files, which in turn can be read by Enzo to set up one or more disk galaxies embedded in live dark-matter halos2.

We restricted ourselves to interactions of Milky Way-like galaxies of equal mass. The parameters chosen for our simulations are summarized in Table 1. Each galaxy initially consists of a gas disk of 1010 M and a dark-matter halo of M200 = 1012M, where M200 is the mass enclosed within the virial radius r200 = 210 kpc. The radial scale length of the disks is rd = 3.5 kpc. The intial conditions were computed using the models outlined in Sect. 2. The initial separation dsep, 0 of the two galaxies is set to the cut-off radius r200 of a dark-matter halo. For a periodic box of 2 Mpc linear size, this choice ensures that the separation of the centers of mass is small compared to the distance from the periodic boundaries. Although the halos are overlapping, the initial separation of the gas disks is significantly larger than their size. For our production runs, we used a 5123 root grid and 8 levels of refinement by overdensity and shear, corresponding to a maximal spatial resolution of 15 pc. The application of the second refinement criterion is important to cover interaction zones with high resolution, especially between the galaxies and in the outer disk regions where the gas density is relatively low. However, this greatly increases the computational cost of the simulations3.

The parameter space of galaxy interactions is spanned by

  • the initial translational velocities of the two disks, V1 and V2,

  • the impact parameter, b, which is the normal distance between a straight line in the direction of the initial velocity of the secondary galaxy, V2, and the center of the primary galaxy,

  • and the disk orientations defined by the angular momentum vectors L1 and L2.

In realistic interaction scenarios, galaxies move on trajectories of various shapes and orientations within clusters and filaments. As long as we consider a system of two interacting galaxies in analogy to the two-body problem and ignore their environment, it is only the motion of the primary relative to the secondary that matters (Toomre & Toomre 1972; Moreno et al. 2015). The relative velocity is increasingly dominated by the gravitational potential as the galaxies approach each other. Consequently, the magnitudes of the initial velocities V1, 2 of the galaxies are not overly important and can be estimated from the halo mass and initial separation.

In our simulation suite, we examined a series of impact parameters combined with varying orientations of the galaxies. To specify the orientations of the galaxies, we used the inclination angles of each galaxy’s angular momentum vector with respect to the vector of the initial relative velocity of the interacting system, Vrel = V1 − V2:

i 1 , 2 = ( L 1 , 2 , V rel ) , $$ \begin{aligned} i_{1,2} = \angle (\boldsymbol{L}_{1,2}, \boldsymbol{V}_{\mathrm{rel} }), \end{aligned} $$(14)

The relative inclination of the disks is given by the difference |i2 − i1|.

In Table 2, we provide an overview of the parameters of all production runs. For referencing, the simulations are numbered sequentially (second column). In addition, simulations are grouped according to the disk inclinations (first column):

  1. Group 0: edge-on collisions with parallel rotation axes perpendicular to relative motion (i1 = i2 = 90°).

  2. Group I: same orientation of primary, secondary inclined by i2 = 135° (|i2 − i1| = 45°).

  3. Group II: same orientation of primary, secondary’s rotation axis aligned with relative motion (i1 = 90°, i2 = 0°, |i2 − i1| = 90°).

  4. Group III: face-on collisions with relative velocity parallel or antiparallel to rotation axes (i1 = i2 = 0°).

  5. Group IV: primary and secondary inclined by i1 = 67.5° and i2 = 112.5°, respectively (|i2 − i1| = 45°).

  6. Group V: primary and secondary inclined by i1 = 45° and i2 = 135°, respectively (|i2 − i1| = 90°).

Table 2.

Overview of varied simulation parameters and minimum distances between galaxy centers.

In each group, the impact parameter b is varied between the limiting scenarios of a central collision (b = 0) and an interaction where the galaxies and their halos are moderately distorted by tidal forces (b on the order of virial radius). In the following, we use the angle defined by sin αb = b/dsep, 0 to specify the impact parameter. For the absolute value of the initial relative velocity, we assumed |Vrel|=v200 in all cases.

Owing to the required computational time and substantial storage requirements, only a few dozen scenarios with more or less arbitrary parameter choices were feasible. Thus, our simulation suite covers some typical scenarios, such as edge-on (group 0) and face-on (group III) collisions of aligned disks, and a selection of intermediate cases with nonzero relative inclinations of the disks. In groups IV and V, we investigated scenarios where the secondary galaxy moves along an inbound trajectory that is inclined with respect to the disk plane of the primary (similar to satellites of the Milky Way, such as the Magellanic Clouds). If feasible, the simulations were run until the post-merger stage. In group 0, this was achieved for Sim. No. 1–5 within 2.5 Gyr, but not for Sim. No. 6. In the other groups, reductions of the fraction of the time step implied by the Courant-Friedrichs-Lewy criterion required terminating the simulations already after 1.5 Gyr in order to save computational time. However, all simulations captured at least the first interaction between the galaxies.

4. Results

4.1. Phenomenology of interacting disks

We begin with a phenomenological discussion of representative examples. In Fig. 1, density slices from an edge-on collision of two galaxies with αb = 20° (Sim. No. 3 in Table 2) are overlaid with a visualization of magnetic field lines using the line integral convolution method4. The first three panels in Fig. 1 show the disks shortly before the collision, the time of closest approach (pericentric passage), and 100 Myr later. The separation of the galaxy centers is plotted as a function of time in Fig. 2. After reversing their motion, the galaxies are pulled toward each other once more and finally merge into a remnant of elliptical shape (Fig. 1, lower right). This is indicated by the time evolution of the shape parameters s and q plotted in Fig. 3. The two parameters specify the ratios of the three principal axes (thin disk: s = 1, q ∼ 0; sphere: q = s = 1). As can be seen, the initially disk-like galaxies are strongly distorted during the interaction stages and evolve into an ellipsoid in the post-merger phase.

thumbnail Fig. 1.

Slices of the gas density in combination with a line integral convolution of the magnetic field in an edge-on collision of two disk galaxies (Sim. No. 3 in Table 2). The first three panels show stages around the first encounter (t = 762.5 Myr). At the end of the simulation (t = 2.5 Gyr), the two galaxies have merged.

thumbnail Fig. 2.

Time evolution of the separation, d, between the centers of interacting galaxies for five cases. Sim. No. 1 and 15 are central collisions (αb = 0°) at different initial inclinations (see Table 2). For Sim. No. 3 (αb = 20°), different interaction stages are visualized in Fig. 1. Sim. No. 6 and 18 are examples for grazing interactions at αb = 45°. The dashed gray line marks a separation of 10 kpc, corresponding to the diameter of the center regions of the initial disks.

The gas density indicated by the color map in Fig. 1 shows that a substantial amount of gas is ripped off the initial disks by tidal forces. The magnetic field structure reveals that the initially ordered toroidal fields inside the two disks evolves into highly turbulent structures expanding into the surroundings of the galaxies. After the first encounter, this can be mainly observed in the region in between the two disks (the so-called bridge) and in the tidal tails. The final remnant exhibits an unordered field that extends far into low-density regions. This suggests that interactions play a role in the magnetization of the circumgalactic medium, although the magnetic field in a realistic cosmological environment is considerably larger than the weak background field of 10−25 G assumed in our simulations.

For nearly central collisions with a minimal separation below 10 kpc, the disks are strongly perturbed after the first interaction (Fig. 4, upper panels). With increasing impact parameter, only the outer regions of the disks interact and less turbulence is produced. A bridge connecting the disks and tidal tails become the most prominent interaction features in the case of a grazing encounter (Fig. 4, lower left). Although the magnetic field remains mostly ordered at this stage, the disks eventually undergo substantial changes and finally merge, as indicated by the time evolution of the separation and shape parameters plotted in Figs. 1 and 3, respectively. The distortion of the disks is reduced further if angular momentum is aligned with velocity (i1 = i2 = 0°), that is to say, if the disk planes of both galaxies are oriented perpendicular to the velocity vector Vrel (Fig. 4, lower right). This type of collision is called face-on (group III in Table 2). Consequently, not only the impact parameter, but also the disk orientation has a significant influence on the outcome of galaxy interactions.

thumbnail Fig. 3.

Evolution of galaxy shape parameters s = b/a (left) and q = c/a (right) of Sim. No. 1, 3, and 6, which represent a sequence of simulations with increasing impact parameter for fixed initial inclinations.

thumbnail Fig. 4.

Slices of the gas density in combination with a line integral convolution of the magnetic field showing post-interaction stages for four different initial orientations. From left to right and top to bottom: Sim. No. 1, 3, 6, and 18 (see Table 2). While the disk inclinations are the same in the first three simulations, Sim. No. 6 and 18 use the impact parameter αb = 45°.

Another example of a face-on collision from group III is shown in Fig. 5. In this simulation, the impact parameter is zero (Sim. No. 15). The galaxies have parallel angular momentum vectors (i.e., the rotation of the two disks is prograde). The morphology after the first encounter is markedly different from the central edge-on collision shown in Fig. 4 (upper left). Instead of tidal tails and a highly turbulent bridge, ring-like perturbations in the outer disks are accompanied by braid-like outward gas streams. As time progresses, more gas is scattered outward and most of the outer disk disperses (Fig. 5, lower left). The remnant after about 1.5 Gyr is highly irregular and does not resemble a disk any longer (Fig. 5, lower right). In the following sections, we investigate whether the morphological differences between different types of interactions are reflected in magnetic field statistics.

thumbnail Fig. 5.

Similar stages as in Fig. 1 for a simulation with impact parameter αb = 0° and inclination angles i1 = i2 = 0° (Sim. No. 15).

4.2. Average magnetic field evolution

As a quantitative measure, we averaged the magnetic field strength in cylindrical regions enclosing the initial gaseous disks (see also Rodenbeck & Schleicher 2016). These regions have a radius of 20 kpc and a height of 4 kpc. By additionally averaging over an inner disk of 5 kpc radial extent, we can distinguish between a central region (r ≤ 5 kpc) and an off-center region (5 kpc < r ≤ 20 kpc).

Figure 6 shows the mean magnetic field strength in the primary galaxy as a function of time for different subsets of simulations. The amplification of the magnetic field by interactions is indicated by the height of peaks. Edge-on collisions (group 0) with impact parameter varying between αb = 0° and 45° are shown in the left plots. While the collision is central for αb = 0°, the impact parameter b is about 150 kpc for αb = 45°. As expected from the phenomenological discussion in Sect. 4.1, the first maximum of the average magnetic field in center regions shows a decline with increasing impact parameter (Fig. 6, upper left). This indicates that the magnetic field in the inner disk is strongly amplified only for nearly central collisions. In contrast, differences between the average magnetic field in off-center regions are less pronounced because the outer disks are tidally disrupted even for relatively large impact parameters (Fig. 6, lower left). The largest field strength in an off-center region is observed for αb = 15° (Sim. No. 2). This suggests that the mixture of shock fronts and turbulence is particularly efficient in driving the amplification of the field if the collision is slightly off-center. During the second encounter, particularly strong fields are produced for αb in the range from 15° and αb = 25° (Sim. No. 2 to 4). For αb ≥ 25°, the maximum field strength in center regions is reached in the third peak. These trends show that in the case of gracing interactions several approaches are required to affect the magnetic field throughout the disks.

thumbnail Fig. 6.

Evolution of the magnetic field strength averaged over center (top) and off-center regions (bottom). Left column: Simulations of edge-on collisions (i1 = i2 = 90°). Right column: Simulations with different disk inclinations sharing the impact parameter αb = 15°. The simulation numbers indicated in the legends refer to Table 2.

Comparing simulations from different groups at a fixed impact parameter (αb = 15°), we see that peaks approximately coincide, but the peak height varies with the disk orientation (right plots in Fig. 6). Peaks are particularly pronounced if the primary’s motion is parallel to its disk plane (i1 = 90°, groups 0, I, and II). For the face-on collision (Sim. No. 16 in group III), magnetic field amplification is significantly reduced both in center and off-center regions. The average magnetic field just prior to the first encounter is also smaller in this case. This is a consequence of stronger gas stripping if the disk moves face-on relative to the ambient medium, which removes magnetic energy from the cylindrical volume in which averages are computed. In groups IV and V, both disks are inclined with respect to their relative motion and their relative inclinations are |i2 − i1| = 45° and 90°, respectively. In these two cases, the peaks are comparable to the face-on collision, with the exception of the center region in Sim. No. 20 during the first encounter. Overall, our findings suggest that the most important parameter is the disk inclination (i.e., the angle between the disk velocity and rotation axis); the relative inclination |i2 − i1| also plays a role, but it is subdominant.

For a systematic comparison of peak heights, we created categorical heat maps of the maximal average magnetic field during first and second encounters (see Fig. 7). In horizontal direction, one can see that the field strength tends to decrease with increasing impact parameter during the first encounters of the disks, while second encounters produce the high peaks also in off-center collisions. This is expected because the separation between the disks in an off-center collision becomes sufficiently small only during the second or third encounter (see also Fig. 2). In contrast, the second peak of the central (αb = 0°) edge-on (group 0) collision is significantly lower than first peak in the central region. Among the central collisions in all groups (vertical direction in the heat maps), the first interactions in group I (relative inclination of 90°) stands out both in center and off-center regions. The central collision in group II (relative inclination of 90°) exhibits strong amplification in the center regions during the second encounter (Fig. 7 upper right). In these cases, the proximity of tilted disks creates particularly strong tidal forces. If the disks are aligned, on the other hand, they will largely overlap at the closest approach, resulting in smaller differential forces. For αb ≥ 15° (second columns), one can discern relatively low average field strengths in groups III to V compared to groups 0 to II (i1 = 90°). This highlights the results shown in Fig. 6. Ignoring outliers for specific disk regions, the overall field amplification appears to be reduced if both disks are inclined with respect to their relative motion. In the case of face-on collisions (group III), the phenomenological discussion in Sect. 4.1 suggests that perturbations in the disk are less efficient in amplifying the magnetic field. Groups IV and V are intermediate cases between edge-on and face-on collisions. In these groups, relatively high peaks are found in collisions with large impact parameters (αb = 30° or 45°) during the second encounter.

thumbnail Fig. 7.

Categorical heat maps displaying the magnetic field strength at the first interaction, ordered by impact parameter angle (αb) and simulation group. Simulations are grouped by common initial inclinations of galaxies, i.e., analogous to Table 2 Group 0: 1–6, Group I: No. 7–10, Group II: 11–14, Group III: 15–18, Group IV: 19–22, Group V: 23–26. Shown are results for center (top) and off-center (bottom) regions at the first (left) and second (right) encounters. Parameter combinations that are not covered by our simulation suite or simulations that were not evolved over a sufficiently long time to reach the second interaction are left white.

In virtually all cases, the peaks in the time evolution of the magnetic field strength are narrow (see Fig. 6), indicating that magnetic field amplification induced by interactions is transient. This behavior can be explained by competing processes. As can be seen in Fig. 1, the magnetic field spreads over surrounding regions through tidal interactions, implying a gradual loss of magnetic field energy inside the disks. Field amplification by shock compression and turbulence acts against these losses. However, these processes are relatively short-lived.

4.3. Mean-field analysis

To analyze dynamo activity in our simulations, we decomposed the numerically computed velocity and magnetic field into mean-field and turbulent components:

B = B ¯ + B turb , $$ \begin{aligned} \boldsymbol{B}&=\overline{\boldsymbol{B}} + \boldsymbol{B}_{\mathrm{turb} },\end{aligned} $$(15)

v = v ¯ + v turb . $$ \begin{aligned} \boldsymbol{v}&=\overline{\boldsymbol{v}} + \boldsymbol{v}_{\mathrm{turb} }. \end{aligned} $$(16)

We followed the approach of Ntormousi et al. (2020) and applied a mass-averaged filter with smoothing scale Δ ¯ $ \overline{\Delta} $ in every spatial dimension to compute the mean-field components. We chose Δ ¯ = 300 $ \overline{\Delta} = 300 $ pc, which corresponds to an average over 53 grid cells (five-point stencil in one dimension). The resulting fluctuations, which are given by the difference between the local and mean fields, are associated with the smallest numerically resolved scales. While the smoothing scale applied by Ntormousi et al. (2020) is comparable to the scale height of a thin cooling disk and supernova-driven bubbles in the disk, our choice can be regarded as a compromise between constraining fluctuations to scales comparable to the grid resolution scale and the strong damping of numerical dissipation in this range of scales (see also Grete et al. 2017).

As shown in Brandenburg & Subramanian (2005), Sect. 6.2, the resulting mean-field induction equation reads

B ¯ t = × ( v ¯ × B ¯ + E ) $$ \begin{aligned} \frac{\partial \overline{\boldsymbol{B}}}{\partial t} = \boldsymbol{\nabla }\times \left(\overline{\boldsymbol{v}}\times \overline{\boldsymbol{B}} + \mathcal{E} \right) \end{aligned} $$(17)

if the magnetic diffusivity η is neglected (ideal MHD approximation). The small-scale electromotive force (EMF) is defined by

E = v turb × B turb ¯ . $$ \begin{aligned} \mathcal{E} = \overline{\boldsymbol{v}_{\mathrm{turb} }\times \boldsymbol{B}_{\mathrm{turb} }}. \end{aligned} $$(18)

In mean-field theories, the small-scale EMF is expanded in powers of the mean magnetic field and current density and closures are applied for the transport coefficients (see Brandenburg & Subramanian 2005, Sect. 6.3). Here, we explicitly evaluate the right-hand side of Eq. (18) from the differences between filtered and local fields. This implies that the value of the EMF depends on the numerical resolution. However, for a fixed resolution, we can compare the evolution of the EMF among different interaction scenarios.

Figure 8 shows the evolution of the small-scale EMF averaged over center regions for the same simulations as in Fig. 6. In each case, we normalized ℰ by ℰ500 (i.e., the average vale at t = 500 Myr). At that point, the galaxies are still in the phase of the first approach, but perturbations in the initial disks have largely relaxed. Thereby, changes in ℰ with respect to the state of the galaxies before the first encounter when ℰ/ℰ500 ≈ 1 can be readily seen. Moreover, we concentrate the discussion on center regions. There is a steep rise of ℰ in at the first encounter in nearly central collisions. In this phase, ℰ/ℰ500 increases by up to two orders of magnitude. Only for larger impact parameters, the increase is smaller or even delayed until the second encounter (Sim. No. 5 and 6). In most simulations, the maximum small-scale EMF is reached around the second encounters. An exception is the central edge-on collisions (Sim. No. 1), where by far the strongest amplification occurs during the initial interaction. Apart from that, small-scale perturbations and, thus, turbulence tend to grow in subsequent interactions. However, this evolution is not sustained. In every simulation that was advanced for a sufficiently long period of time, we see that ℰ gradually decreases in the final phase, when the galaxies have merged and small-scale flow is decaying. Although the averaged magnetic field is produced by both the mean-field and the small-scale EMF, its time evolution (Fig. 6, upper left and right) roughly reflects the evolution of ℰ. This suggests that ℰ contributes significantly to the amplification of the magnetic field, particularly during close encounters. In the intermediate phases between encounters, ℰ/ℰ500 ∼ 10 in most cases (see Fig. 8), while the average magnetic field decreases.

thumbnail Fig. 8.

Small-scale EMF, ℰ, averaged over center regions for simulations of edge-on collisions (left) and simulations with the impact parameter αb = 15° (right). Top row: EMF normalized with respect to its value at t = 500 Myr. Bottom row: ℰ relative to v ¯ × B ¯ $ \overline{\mathbf{v}}\times\overline{\mathbf{B}} $ (see Eq. (17)) averaged over center regions for the same simulations as in the top row. The corresponding magnetic field evolution is shown in Fig. 6.

To investigate the contribution of the small-scale part ℰ relative to the mean-field part v ¯ × B ¯ $ \overline{\boldsymbol{v}}\times\overline{\boldsymbol{B}} $ of the EMF in Eq. (17), we averaged the ratio of these two contributions over the disk regions. Results for center regions are shown in the bottom plots in Fig. 8. Generally, first encounters are associated with downward spikes, indicating that magnetic field amplification is mostly driven by large-scale processes, such as compression and tidal stresses. However, the ratio E / ( v ¯ × B ¯ ) $ \mathcal{E}/(\overline{\boldsymbol{v}}\times\overline{\boldsymbol{B}}) $ rises in the aftermath of the initial interaction to levels of about 5–10%. This can be understood as a consequence of turbulence production, as illustrated in Fig. 4. For fixed impact parameter (Figure 8, lower right), one can clearly see that the ratio also drops at the second encounter, but the small-scale EMF continues to grow until coalescence. In the post-merger phase, both contributions gradually decay. In the group 0 (Fig. 8, lower left), the pre-interaction ratio of about 0.03 is reached after about 2.5 Gyr for nearly central collisions. For larger impact parameters, the small-scale EMF tends to persists longer (Sim. No. 4–6). This indicates that turbulence is mainly produced during late stages, when the mean-field contribution is already relatively small.

One of the simplest dynamo models, the so-called α effect, assumes that ℰ is related to the turbulent kinetic helicity (see Brandenburg & Subramanian 2005, Sect. 6.3). The turbulent kinetic helicity is defined by the integral

K turb = v turb · ( × v turb ) d V , $$ \begin{aligned} K_{\mathrm{turb} } = \int {\boldsymbol{v}_\mathrm{turb} \cdot (\boldsymbol{\nabla }\times \boldsymbol{v}_\mathrm{turb} )\,\mathrm{d} V}, \end{aligned} $$(19)

which we evaluated within two adjacent cylindrical volumes that are immediately below ( K turb < $ K_{\mathrm{turb}}^< $) and above ( K turb > $ K_{\mathrm{turb}}^> $) the midplane, respectively (see also Ntormousi et al. 2020). With a thickness of 4 kpc and a radial extent of 20 kpc, the combined volume above and below the midplane encompasses the center and off-center regions for computing averages of the magnetic field strength. Since helicity is a pseudoscalar, reflection across the midplane reverses its sign (Yokoi 2013). Prior to interaction, K turb > $ K_{\mathrm{turb}}^> $ and K turb < $ K_{\mathrm{turb}}^< $ are comparable and have opposite signs. The time evolution is shown for selected simulations in Fig. 9. In the case of a central edge-on collision (Sim. No. 0), K turb > $ K_{\mathrm{turb}}^> $ has a pronounced negative peak at the first encounter, which is roughly balanced by a positive peak of K turb < $ K_{\mathrm{turb}}^< $. This phase is dominated by strong helical flows inside the tidal arms and the tidal bridge connecting the retracting collision partners. Assuming a linear relation between ℰ and Kturb and an initially coherent toroidal mean field, the α effect preferentially amplifies the magnetic field in opposite directions above and below the midplane, which results in a quadrupolar field component (Ntormousi et al. 2020). After the first encounter, both K turb > $ K_{\mathrm{turb}}^> $ and K turb < $ K_{\mathrm{turb}}^< $ become small and irregular. For larger impact parameters (Fig. 9, Sim. No. 3 and 6), K turb > $ K_{\mathrm{turb}}^> $ and K turb < $ K_{\mathrm{turb}}^< $ remain antisymmetric in the post-interaction phase and throughout the second encounter. In the case of Sim. No. 3 (αb = 20°), antisymmetry is eventually broken in the coalescence phase. A similar evolution is seen in Sim. No. 6 (αb = 45°), but instead of the pronounced peaks at the second encounter, the average kinetic helicity decreases gradually at late time. Comparing to Fig. 8, pronounced antisymmetric peaks of K turb > $ K_{\mathrm{turb}}^> $ and K turb < $ K_{\mathrm{turb}}^< $ are associated with a decrease in E / ( v ¯ × B ¯ ) $ \mathcal{E}/(\overline{\boldsymbol{v}}\times\overline{\boldsymbol{B}}) $. If the magnetic field is still more or less toroidal at that point, the resulting small-scale EMF tends to be antisymmetric and cancels out over the whole cylindrical region. The increasing EMF in post-interaction phases, on the other hand, indicates a growing influence of non-helical isotropic turbulence.

thumbnail Fig. 9.

Turbulent kinetic helicity defined by Eq. (17) averaged over centered cylinders (r ≤ 5 kpc, |z|≤4 kpc) below and above the galaxy disk midplane. Left column: Simulations of edge-on collisions (group 0 with αb = 0°, 20°, and 45°). Right column: Simulations with the impact parameter αb = 15° (groups II, III, and IV).

For scenarios where the disk inclination differs from 90°, we see a markedly different evolution (Fig. 9, right column). As an example for face-on collisions (group III), Sim. No. 16 is shown. In this case, K turb > $ K_{\mathrm{turb}}^> $ and K turb < $ K_{\mathrm{turb}}^< $ have the same sign after the initial approach, meaning that the contributions from above and below the midplane are symmetric rather than antisymmetric. Since the two disks have parallel angular momentum, one would expect that the mean kinetic helicity of the lower half of one disk cancels the upper half of the other disk when they overlap. However, it appears that tidal interactions between the two disks fundamentally change their structure, resulting in a positive (first encounter) or negative (second encounter) net helicity. To a lesser degree this can also be seen in Sim. No. 20 (group IV, both disks are inclined), while Sim. No. 12 (group II, only secondary disk has inclination 0°) is a mixed case with alternating antisymmetric and symmetric episodes. These qualitative differences in the kinetic helicity evolution correspond well to differences in the small-scale EMF and the relatively weak magnetic field amplification in groups III to V, where both disks are inclined.

5. Conclusions

We explored the evolution of galactic magnetic fields during interactions and mergers of disk galaxies by means of numerical simulations. Gaseous disks in hydrostatic equilibrium were initialized using the potential method from Wang et al. (2010). Moreover, we adopted the method from Drakos et al. (2017) to compute particle initial conditions for live dark-matter halos. To compute the evolution of disks and halos, we employed the MHD AMR code Enzo. This allowed us to analyze magnetic field amplification induced by the tidal interactions between the galaxies. In an elaborate parameter study, we varied the impact parameter (i.e., the normal distance of the trajectory of the inbound galaxy to the target galaxy) and the relative orientations of disks.

During interaction phases, we find pronounced peaks in the magnetic field strength averaged over the center and outer regions of the disks. These peaks are a distinctive feature that was also identified by Drzazga et al. (2011) in radio observations of interacting galaxies. However, there are important differences. For all the galaxies in their sample, Drzazga et al. (2011) found a relatively small increase during encounters, and the strongest fields they found are associated with the final coalescence phase. In our suite of simulations, this is only the case for impact parameters αb (the angle defined by the ratio of the linear impact parameter and the initial separation) greater than about 20° or, to a lesser degree, for collisions of disks that are inclined with respect to their orbital motion and whose axes are not aligned (Fig. 6). Nearly central collisions typically exhibit the highest peaks during the first encounter. If these trends are correct, this implies that nearly central collisions are rare and were not present in the small galaxy sample investigated by Drzazga et al. (2011), who ordered the galaxies by interaction stages. However, our results indicate that, even at the same interaction stage, variations in the field strengths may result from different interaction parameters. Since three-dimensional data are required, it will be challenging to lift this degeneracy by classifying observed systems accordingly.

Moreover, the mean field strength in center regions exceeds 10 μG in about half of the scenarios we investigated, while the maximum mean field strength is typically in the range from 5 to 10 μG in off-center regions (Fig. 7). These values are systematically lower than the field strengths inferred by Drzazga et al. (2011) from nonthermal radio emission, which reach 50 μG and 20 μG in center and off-center regions, respectively. This discrepancy might be a consequence of insufficient numerical resolution resulting in a less effective turbulent dynamo. However, comparable field strengths were obtained by Rodenbeck & Schleicher (2016) despite the significantly lower resolution of their simulations (240 pc compared to 15 pc in our simulations). The most likely explanation is the assumption of adiabatic gas dynamics. Since the gas cannot cool, adiabatic disks do not contract into thin disks in our simulations. Thus, they have larger scale heights and the initial magnetic field undergoes less compression (assuming flux conservation), resulting in an overall lower field strength. This is supported by relative amplification factors of about 2–3 (Fig. 6), which compares well to the Drzazga et al. (2011) results. In addition, we neglected magnetic field amplification due to turbulence driven by feedback processes. Whittingham et al. (2021) do reproduce observed galactic field strengths in cosmological zoom-in simulations with nonadiabatic physics and stellar feedback.

As a quantitative indicator, Whittingham et al. (2021) show that the power spectra of the magnetic field energy in the central regions of interacting galaxies approximately agree with Kolmogorov scaling toward small length scales. This suggests that the field is turbulent (Beresnyak 2019). However, power spectra are of limited validity for complex structures such as interacting galaxies, where the intensity of turbulence undergoes substantial spatial variations and rapid changes in time. Moreover, there are localized compression effects, for example shock compression of gas in collision zones. For this reason, we applied a box filter to decompose the velocity and magnetic field into a mean-field part and small-scale fluctuations.

The amplification of the magnetic field during encounters and coalescence is reflected by an enhanced EMF, which subsides in the aftermath of the interaction (Fig. 8, top row). As a result, the amplification is transient. Basically, this agrees with the observational data from Drzazga et al. (2011) and the numerical studies (Pakmor et al. 2014; Rodenbeck & Schleicher 2016; Whittingham et al. 2021). By comparing small-scale and mean-field components, we find that the turbulent, small-scale EMF grows after the first interaction to a fraction on the order of 0.1 of the mean-field contribution. However, when the galaxies approach each other and get close to the orbital pericenter, the ratio of small-scale to mean-field EMF decreases in many cases (Fig. 8, bottom row). These minima are typically associated with antisymmetric peaks of the turbulent kinetic helicity above and below the midplane (Fig. 9). Antisymmetric kinetic helicity was previously reported by Ntormousi et al. (2020) for isolated disk galaxies. As pointed out by Brandenburg & Subramanian (2005), inhomogeneous and anisotropic flows typically have large helicities and act as large-scale dynamos. Such flows are produced by tidal forces, and, consequently, the mean-field contribution should be relatively large during close encounters. The small-scale EMF becomes more significant in post-interaction phases in which the gas becomes increasingly turbulent and the galaxies finally merge. The turbulent kinetic helicity switches from antisymmetry to an irregular regime and decays in the remnant. If the rotation axis of both galaxies are inclined, such as in face-on collisions, a significant net helicity develops during the first encounter. In these scenarios, the maximal amplification of the magnetic field is relatively low (Fig. 7).

Overall, this led to our interpretation that the peaks in the average magnetic field are mainly due to tidally induced helical flows. This process is more efficient if the disk plane is oriented parallel to relative motion. The production of turbulence enhances the magnetic field, especially in post-interaction phases. However, the small-scale dynamo is not sustained for long after coalescence. A lasting effect is the unordered, random structure of the magnetic field in the remnants, which are always ellipsoids. However, in some cases this shape might be an artifact of adiabatic gas dynamics. It is conceivable that radiative cooling allows for disk-like remnants, where a coherent field could be rejuvenated by the α-Ω dynamo. Apart from the highly dynamic evolution of magnetic fields, star formation and, in turn, feedback are triggered by galaxy collisions (Moreno et al. 2021; Whittingham et al. 2021; Linden & Mihos 2022; Renaud et al. 2022). Moreover, the dependence of the small-scale dynamo on numerical resolution in galaxy simulations is a pending problem (Grete et al. 2019; Körtgen et al. 2019; Ntormousi et al. 2020; Whittingham et al. 2021; Liu et al. 2022; Wissing & Shen 2023). As shown by Schober et al. (2013), the small-scale dynamo first amplifies the magnetic field to saturation on the viscous scale. However, this scale is in the sub-parsec range, far below scales that can be resolved in galaxy simulations. Subsequently, the magnetic energy is shifted to larger scales until it saturates on the forcing scale. In a numerical simulation, the smallest scale is set by numerical diffusion. Consequently, the turbulent EMF computed from our simulations initially operates on scales that are orders of magnitude larger than the physical amplification scales. Numerical diffusion also determines the effective resistive scale on which magnetic energy is dissipated. Consequently, growth and decay rates may very well change with resolution. Since we have demonstrated that the small-scale EMF plays a significant role in the turbulent environment of interacting galaxies, improved modeling of the underlying dynamo and dissipation processes will be an important next step.


1

See the website https://enzo-project.org for further details. For the simulation suite presented in this article, the public repository https://github.com/SimonCSelg/enzo-dev_mhdsgs-diskgalaxy, commit fc68caa, was used. In order to maintain consistency between the simulations, this repository was frozen. Subsequent modifications are publicly available at https://github.com/wolfram-schmidt/enzo-dev/tree/diskgalaxy

2

The Python modules for halos and equilibrium disks are publicly available at https://github.com/wolfram-schmidt/isolated-galaxy

3

Using about 1200 MPI tasks, typically 1 million core-h were needed to complete one simulation. Owing to inherent limitations of the parallelization scheme used in Enzo, it is not feasible to run the code with more than a few thousand tasks.

4

Like iron filings surrounding a bar magnet, magnetic field lines are visualized by filtering a texture consisting of white-noise along local streamlines derived from the magnetic field data.

Acknowledgments

We thank Kai Rodenbeck and Dominik Schleicher for their pioneering work and for providing the code for the equilibrium disk setup in Enzo. The work presented in this article was supported by the German Deutsche Forschungsgemeinschaft, DFG project number SCHM 2135/6-1. The simulations were carried out using computational resources granted by HLRN project hhp00049.

References

  1. Arshakian, T. G., Beck, R., Krause, M., & Sokoloff, D. 2009, A&A, 494, 21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Beresnyak, A. 2019, Liv. Rev. Comput. Astrophys., 5, 2 [CrossRef] [Google Scholar]
  3. Bhat, P., & Subramanian, K. 2014, ApJ, 791, L34 [NASA ADS] [CrossRef] [Google Scholar]
  4. Binney, J., & Tremaine, S. 2008, Galactic Dynamics: Second Edition, 2nd edn. (Princton University Press) [Google Scholar]
  5. Blumenthal, K. A., & Barnes, J. E. 2018, MNRAS, 479, 3952 [Google Scholar]
  6. Brandenburg, A., & Subramanian, K. 2005, Phys. Rep., 417, 1 [NASA ADS] [CrossRef] [Google Scholar]
  7. Brummel-Smith, C., Bryan, G., Butsky, I., et al. 2019, J. Open Source Softw., 4, 1636 [NASA ADS] [CrossRef] [Google Scholar]
  8. Bryan, G. L., Norman, M. L., O’Shea, B. W., et al. 2014, ApJS, 211, 19 [Google Scholar]
  9. Chyży, K. T., & Beck, R. 2004, A&A, 417, 541 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Condon, J. J., Helou, G., & Jarrett, T. H. 2002, ApJ, 123, 1881 [CrossRef] [Google Scholar]
  11. Di Matteo, P., Bournaud, F., Martig, M., et al. 2008, A&A, 492, 31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Dolag, K., & Stasyszyn, F. 2009, MNRAS, 398, 1678 [NASA ADS] [CrossRef] [Google Scholar]
  13. Drakos, N. E., Taylor, J. E., & Benson, A. J. 2017, MNRAS, 468, 2345 [NASA ADS] [CrossRef] [Google Scholar]
  14. Drzazga, R. T., Chyży, K. T., Jurusik, W., & Wiórkiewicz, K. 2011, A&A, 533, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Geng, A., Kotarba, H., Bürzle, F., et al. 2012, MNRAS, 419, 3571 [NASA ADS] [CrossRef] [Google Scholar]
  16. Georgakakis, A., Forbes, D. A., & Norris, R. P. 2000, MNRAS, 318, 124 [NASA ADS] [CrossRef] [Google Scholar]
  17. Grete, P., Vlaykov, D. G., Schmidt, W., & Schleicher, D. R. G. 2017, Phys. Rev. E, 95, 033206 [NASA ADS] [CrossRef] [Google Scholar]
  18. Grete, P., Latif, M. A., Schleicher, D. R. G., & Schmidt, W. 2019, MNRAS, 487, 4525 [NASA ADS] [CrossRef] [Google Scholar]
  19. Hayward, C. C., Torrey, P., Springel, V., Hernquist, L., & Vogelsberger, M. 2014, MNRAS, 442, 1992 [NASA ADS] [CrossRef] [Google Scholar]
  20. Hernquist, L. 1990, ApJ, 356, 359 [Google Scholar]
  21. Hopkins, P. F., Kereš, D., & Murray, N. 2013, MNRAS, 432, 2639 [NASA ADS] [CrossRef] [Google Scholar]
  22. Jiménez, E., Lagos, C. D. P., Ludlow, A. D., & Wisnioski, E., 2023, MNRAS, 524, 4346 [CrossRef] [Google Scholar]
  23. Kennicutt, R. C., Jr 1998, ARA&A, 36, 189 [NASA ADS] [CrossRef] [Google Scholar]
  24. Kotarba, H., Karl, S. J., Naab, T., et al. 2010, ApJ, 716, 1438 [NASA ADS] [CrossRef] [Google Scholar]
  25. Krumholz, M. R., & Burkhart, B. 2016, MNRAS, 458, 1671 [NASA ADS] [CrossRef] [Google Scholar]
  26. Körtgen, B., Banerjee, R., Pudritz, R. E., & Schmidt, W. 2019, MNRAS, 489, 5004 [CrossRef] [Google Scholar]
  27. Linden, S. T., & Mihos, J. C. 2022, ApJ, 933, L33 [NASA ADS] [CrossRef] [Google Scholar]
  28. Liu, Y., Kretschmer, M., & Teyssier, R. 2022, MNRAS, 513, 6028 [NASA ADS] [Google Scholar]
  29. Martin, G., Kaviraj, S., Devriendt, J. E. G., Dubois, Y., & Pichon, C. 2018, MNRAS, 480, 2266 [Google Scholar]
  30. Mihos, J. C., & Hernquist, L. 1996, ApJ, 464, 641 [Google Scholar]
  31. Moreno, J., Torrey, P., Ellison, S. L., et al. 2015, MNRAS, 448, 1107 [NASA ADS] [CrossRef] [Google Scholar]
  32. Moreno, J., Torrey, P., Ellison, S. L., et al. 2019, MNRAS, 485, 1320 [NASA ADS] [CrossRef] [Google Scholar]
  33. Moreno, J., Torrey, P., Ellison, S. L., et al. 2021, MNRAS, 503, 3113 [NASA ADS] [CrossRef] [Google Scholar]
  34. Moss, D., Sokoloff, D., Beck, R., & Krause, M. 2014, A&A, 566, A40 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. Ntormousi, E., Tassis, K., Del Sordo, F., Fragkoudi, F., & Pakmor, R. 2020, A&A, 641, A165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Pakmor, R., Marinacci, F., & Springel, V. 2014, ApJ, 783, L20 [CrossRef] [Google Scholar]
  37. Park, K., Blackman, E. G., & Subramanian, K. 2013, Phys. Rev. E, 87, 053110 [NASA ADS] [CrossRef] [Google Scholar]
  38. Renaud, F., Bournaud, F., Kraljic, K., & Duc, P.-A. 2014, MNRAS, 442, L33 [NASA ADS] [CrossRef] [Google Scholar]
  39. Renaud, F., Bournaud, F., & Duc, P.-A. 2015, MNRAS, 446, 2038 [Google Scholar]
  40. Renaud, F., Segovia Otero, Á., & Agertz, O. 2022, MNRAS, 516, 4922 [NASA ADS] [CrossRef] [Google Scholar]
  41. Rodenbeck, K., & Schleicher, D. R. G. 2016, A&A, 593, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Rossa, J., Laine, S., van der Marel, R. P., et al. 2007, ApJ, 134, 2124 [Google Scholar]
  43. Schober, J., Schleicher, D. R. G., & Klessen, R. S. 2013, A&A, 560, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Smith, B. J., Giroux, M. L., Struck, C., & Hancock, M. 2010, ApJ, 139, 1212 [CrossRef] [Google Scholar]
  45. Springel, V., Di Matteo, T., & Hernquist, L. 2005, MNRAS, 361, 776 [Google Scholar]
  46. Steinwandel, U. P., Beck, M. C., Arth, A., et al. 2019, MNRAS, 483, 1008 [NASA ADS] [CrossRef] [Google Scholar]
  47. Steinwandel, U. P., Dolag, K., Lesch, H., et al. 2020, MNRAS, 494, 4393 [NASA ADS] [CrossRef] [Google Scholar]
  48. Struck, C. 1999, Phys. Rep., 321, 1 [NASA ADS] [CrossRef] [Google Scholar]
  49. Teyssier, R., Chapon, D., & Bournaud, F. 2010, ApJ, 720, L149 [NASA ADS] [CrossRef] [Google Scholar]
  50. Toomre, A., & Toomre, J. 1972, ApJ, 178, 623 [Google Scholar]
  51. Wang, H.-H., Klessen, R. S., Dullemond, C. P., van den Bosch, F. C., & Fuchs, B. 2010, MNRAS, 407, 705 [CrossRef] [Google Scholar]
  52. Whittingham, J., Sparre, M., Pfrommer, C., & Pakmor, R. 2021, MNRAS, 506, 229 [NASA ADS] [CrossRef] [Google Scholar]
  53. Wissing, R., & Shen, S. 2023, A&A, 673, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Yokoi, N. 2013, Geophys. Astrophys. Fluid Dyn., 107, 114 [Google Scholar]

All Tables

Table 1.

Overview of the initial disk and halo model parameters.

Table 2.

Overview of varied simulation parameters and minimum distances between galaxy centers.

All Figures

thumbnail Fig. 1.

Slices of the gas density in combination with a line integral convolution of the magnetic field in an edge-on collision of two disk galaxies (Sim. No. 3 in Table 2). The first three panels show stages around the first encounter (t = 762.5 Myr). At the end of the simulation (t = 2.5 Gyr), the two galaxies have merged.

In the text
thumbnail Fig. 2.

Time evolution of the separation, d, between the centers of interacting galaxies for five cases. Sim. No. 1 and 15 are central collisions (αb = 0°) at different initial inclinations (see Table 2). For Sim. No. 3 (αb = 20°), different interaction stages are visualized in Fig. 1. Sim. No. 6 and 18 are examples for grazing interactions at αb = 45°. The dashed gray line marks a separation of 10 kpc, corresponding to the diameter of the center regions of the initial disks.

In the text
thumbnail Fig. 3.

Evolution of galaxy shape parameters s = b/a (left) and q = c/a (right) of Sim. No. 1, 3, and 6, which represent a sequence of simulations with increasing impact parameter for fixed initial inclinations.

In the text
thumbnail Fig. 4.

Slices of the gas density in combination with a line integral convolution of the magnetic field showing post-interaction stages for four different initial orientations. From left to right and top to bottom: Sim. No. 1, 3, 6, and 18 (see Table 2). While the disk inclinations are the same in the first three simulations, Sim. No. 6 and 18 use the impact parameter αb = 45°.

In the text
thumbnail Fig. 5.

Similar stages as in Fig. 1 for a simulation with impact parameter αb = 0° and inclination angles i1 = i2 = 0° (Sim. No. 15).

In the text
thumbnail Fig. 6.

Evolution of the magnetic field strength averaged over center (top) and off-center regions (bottom). Left column: Simulations of edge-on collisions (i1 = i2 = 90°). Right column: Simulations with different disk inclinations sharing the impact parameter αb = 15°. The simulation numbers indicated in the legends refer to Table 2.

In the text
thumbnail Fig. 7.

Categorical heat maps displaying the magnetic field strength at the first interaction, ordered by impact parameter angle (αb) and simulation group. Simulations are grouped by common initial inclinations of galaxies, i.e., analogous to Table 2 Group 0: 1–6, Group I: No. 7–10, Group II: 11–14, Group III: 15–18, Group IV: 19–22, Group V: 23–26. Shown are results for center (top) and off-center (bottom) regions at the first (left) and second (right) encounters. Parameter combinations that are not covered by our simulation suite or simulations that were not evolved over a sufficiently long time to reach the second interaction are left white.

In the text
thumbnail Fig. 8.

Small-scale EMF, ℰ, averaged over center regions for simulations of edge-on collisions (left) and simulations with the impact parameter αb = 15° (right). Top row: EMF normalized with respect to its value at t = 500 Myr. Bottom row: ℰ relative to v ¯ × B ¯ $ \overline{\mathbf{v}}\times\overline{\mathbf{B}} $ (see Eq. (17)) averaged over center regions for the same simulations as in the top row. The corresponding magnetic field evolution is shown in Fig. 6.

In the text
thumbnail Fig. 9.

Turbulent kinetic helicity defined by Eq. (17) averaged over centered cylinders (r ≤ 5 kpc, |z|≤4 kpc) below and above the galaxy disk midplane. Left column: Simulations of edge-on collisions (group 0 with αb = 0°, 20°, and 45°). Right column: Simulations with the impact parameter αb = 15° (groups II, III, and IV).

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.