© The Authors 2025

1 Introduction

Disc galaxies, characterised by their flattened stellar distributions and rotationally supported dynamics, represent more than half of the total local galactic population (e.g. Lintott et al. 2008; Willett et al. 2013). They consist of stars, gas, dust, and dark matter, and often display non-axisymmetric structures such as spiral arms and bars (e.g. de Vaucouleurs & Freeman 1972; Bamford et al. 2009). It has long been known that galaxies have likely merged together during the lifetime of our Universe (Tremaine 1981), further compounding the natural complexity of these large-scale systems. Therefore, to gain insight into galactic dynamics, the study of galactic interactions is essential.

Because of their large spatial extent, galaxies experience non-uniform gravitational forces during interactions, which drive both their relative motion and their internal reshaping. The effect induced by the gradient of the gravitational force is called a tidal interaction and, as presented in the seminal article by Toomre & Toomre (1972), can significantly transform galaxies. Further research demonstrated the ability of tidal interactions to induce star formation bursts in galaxies (e.g. Lotz et al. 2008), drastically change their morphology (e.g. White 1978), heat and thicken their disc, strongly weaken radial metallicity gradients (e.g. Bustamante et al. 2018), displace their centre of mass (creating a ‘reflex motion’; e.g. Erkal et al. 2019a; Petersen & Peñarrubia 2021; Brooks et al. 2025b; Yaaqib et al. 2025), impact the morphology and kinematics of nearby stellar streams (e.g. Erkal et al. 2019b; Koposov et al. 2023; Brooks et al. 2024; Jiménez-Arranz et al. 2025b), and create new elliptical galaxies (e.g. Bournaud et al. 2007), among other effects. Thus, understanding the dynamical response of galaxies to tidal interactions constitutes a major objective in contemporary astrophysics.

Tidal interactions are not the only mechanism able to perturb a galactic disc. The interplay between non-axisymmetric patterns and the broad variety of coexisting stellar orbits allows the angular momentum to be redistributed in disc galaxies. Thus, some stars can undergo significant and durable radial excursions: this mechanism is known as radial migration (Sellwood & Binney 2002). The key property for quantifying the radial migration of a star orbiting within a disc galaxy is its guiding radius, R_g(t): the equivalent radius of a circular orbit of a star with the same angular momentum. In practice, the variation in the guiding radii of stars in the galactic disc is called churning (as introduced in Schönrich & Binney 2009) and directly implies a redistribution of the stellar population in the galaxy. This churning was initially described as being caused by transient spiral arms formation (see Sellwood 2014, for a review on the topic) but is now also linked to bar instabilities and resonances with rotating patterns in disc galaxies (e.g. Minchev et al. 2011; Halle et al. 2015; Nagy et al. 2025; Marques et al. 2025). Regarding examples of its consequences: radial migration can significantly weaken the radial metallicity gradient of galaxies (e.g. Kubryk et al. 2013; Kubryk et al. 2015) and increase the metallicity dispersion in the disc (e.g. Sellwood & Binney 2002; Grand et al. 2015).

In this context, radial migration has been primarily examined in simulations of isolated disc galaxies (e.g. Halle et al. 2015; Vera-Ciro et al. 2016; Halle et al. 2018; Mikkola et al. 2020; Iles et al. 2024) and, to a lesser extent, with satellite galaxies (e.g. Quillen et al. 2009; Bird et al. 2012). More recently, radial migration has been investigated within Milky Way (MW)-Sagittarius simulations (Sgr; e.g. Carr et al. 2022) and cosmological simulations (e.g. Lu et al. 2024; Minchev et al. 2025; Bernaldez et al. 2025), where the influence of tidal interactions on stellar churning is explicitly considered. The cosmological simulations show that tidally induced migration can significantly modify the inferred star formation histories of galaxies, as well as other key observational properties. However, while cosmological simulations have provided insight into radial migration, their limited resolution and complex environments make systematic investigations difficult. Thus, the impact of interactions on radial migration in high-resolution simulations remains insufficiently addressed, forming the basis for our analysis.

On the other hand, we see that our neighbourhood is no stranger to tidal interactions, with the MW’s history including numerous mergers and satellite galaxy disruptions (e.g. Unavane et al. 1996; Antoja et al. 2018; Merrow et al. 2024; Hunt & Vasiliev 2025). Among these merged satellites, notable examples are the Sgr dwarf spheroidal galaxy (see Ibata et al. 1995) and the Gaia Enceladus/Sausage (see Helmi et al. 2018; Belokurov et al. 2018). In addition to this merger history, several satellite galaxies populate the surroundings of our Galaxy. Among these, the Large Magellanic Cloud (LMC) – the closest disc galaxy to the MW, located approximately 50 kpc from the Sun (Graczyk et al. 2013; Pietrzyński et al. 2019) and currently interacting with the Small Magellanic Cloud (SMC) – represents a particularly relevant target for our study.

Because of this proximity, the LMC is a uniquely suited laboratory for diagnosing radial migration, especially under the influence of tidal interactions. It is a warped barred Magellanic spiral (de Vaucouleurs & Freeman 1972) featuring a single prominent spiral arm connected to its slightly off-centred bar (El Youssoufi et al. 2019; Gaia Collaboration 2021; Jiménez-Arranz et al. 2023a; Jiménez-Arranz et al. 2025a; Rathore et al. 2025), which has recently been found to be non-rotating (Jiménez-Arranz et al. 2024a). The LMC is undergoing tidal interactions with the SMC, providing a natural explanation for the LMC’s non-rotating bar (see Jiménez-Arranz & Roca-Fàbrega 2025), and the MW (e.g. Garavito-Camargo et al. 2019; Vasiliev 2023). Owing to its favourable position in the sky and proximity, the LMC has been the target of numerous observational surveys (e.g. Zaritsky et al. 1997; Smith et al. 2005; de Jong et al. 2019; Massana et al. 2022). With Gaia resolving around ~15-20 million of its stars (Gaia Collaboration 2021; Jiménez-Arranz et al. 2023a,b; Perryman 2025), the system is not only exceptionally well characterised observationally but has also been the subject of extensive numerical modelling (e.g. Besla et al. 2012; Garavito-Camargo et al. 2019; Brooks et al. 2025a; Sheng et al. 2025). In particular, in this work we used the KRATOS (Kinematic Reconstruction of the mAgellanic sysTem within the OCRE Scenario) suite of simulations (JA24) to characterise how tidal interactions drive radial migration in LMC-like systems.

The KRATOS suite comprises 28 pure N-body simulations of the LMC-SMC-MW system, which naturally reproduce several observed features of the LMC, including the offset of its central non-rotating bar (Jiménez-Arranz & Roca-Fàbrega 2025), the warped structure of its disc (Jiménez-Arranz et al. 2025a), and its kinematic asymmetries (Schölch et al. 2025). The diverse initial conditions across the simulations yield a variety of stellar bar and spiral arm morphologies in the LMC-like systems, enabling an extension of the analysis to the broader subclass of Magellanic spiral galaxies. Furthermore, KRATOS has demonstrated that tidal interactions in LMC analogues can trigger a wide array of mechanisms that drive radial migration (see the main work of KRATOS: JA24), including connection between a bar and spiral arms (e.g. Marques et al. 2025) and changes in the bar pattern speed during mergers (e.g. Zhang et al. 2025; Jiménez-Arranz & Roca-Fàbrega 2025). Consequently, these simulations are expected to host a rich diversity of churning modes, making them particularly well suited for quantifying the impact of tidal interactions on radial migration in disc galaxies. Building on KRATOS, we aim to provide a novel characterisation of how tidal interactions influence the radial redistribution of both stars and metallicities in disc galaxies.

In Sect. 2, we describe the KRATOS simulations that were the basis of this work. In Sect. 3, we present the methodology and assumptions used to quantify radial migration from the KRATOS output data. In Sect. 4, we present our results, including a tidally induced regime of oscillating churning fluxes in LMC-like galaxies. In Sect. 5, we discuss our results, the limitations of our assumptions, and future prospects in this field, before concluding this paper in Sect. 6.

2 The KRATOS suite of simulations

KRATOS (originally presented in JA24) is a high-resolution suite of 28 simulations aiming to study the LMC-SMC-MW interacting system. These simulations are pure N-body simulations, namely, without hydrodynamics and stellar evolution or formation. This suite of 28 simulations allows the evolution of an isolated LMC-like galaxy to be compared with its counterpart, which is interacting with an SMC-mass galaxy and an MW-mass galaxy. In KRATOS, the explored parameters are the mass of the baryonic matter or dark matter content of the different systems and the stellar Toomre parameter Q of the LMC-like system (see Table 2 of JA24, for a summary). In this work we follow the notation introduced in JA24, where the LMC-like galaxy is denoted as G_LMC, and the SMC- and MW-mass systems as G_SMC and G_MW, respectively.

The simulations feature the stellar mass resolution of 4 × 10³ M_⊙, along with a spatial and temporal resolutions of 10 pc and 5000 yr. In the KRATOS suite, none of the galaxies is modelled as an analytical potential and they are simulated for 4.68 Gyr within a 2.85³ Mpc³ box with periodic boundary conditions, using the adaptive refinement mesh Eulerian code ART presented in Kravtsov et al. (1997). For this work, we focused on a subsample of six simulations of KRATOS (see Table 1). They were re-ran so that each simulation provides more than 2000 snapshots¹ of the evolution of the systems, with a time separation lower than 2.6 Myr. The present time t=0 is defined 4.0 Gyr after the initial conditions, and was determined to match the number of G_LMC−G_SMC pericentre passages and the morphology of the G_LMC disc (see Sect. 3 of JA24, for more details). For reference, we mention that the time t=0 was re-evaluated in Sect. 4.3 of Schölch et al. (2025) to match the kinematic asymmetries of the disc, but we used only the original definition. The subsample of the KRATOS simulations used in this work may include the following systems:

• G_LMC: an LMC-like (i.e. similar in both mass and shape to the LMC) galaxy with a stellar mass of 5.0 × 10⁹ M_⊙ and a dark matter mass of 1.8 × 10¹¹ M_⊙. It features 1.2 × 10⁶ stellar particles with a stellar mass resolution of 4 × 10³ M_⊙. G_LMC presents an exponential stellar disc profile with a scale height and length of, respectively, 0.20 kpc and 2.85 kpc truncated at 11.5 kpc from the centre, embedded in a dark matter Navarro-Frenk-White (NFW) halo (Navarro et al. 1996) that evolves during the simulations.

• G_SMC: an SMC-mass (i.e. only aiming to reproduce the mass of the SMC) galaxy with a stellar mass of 2.6 × 10⁸ M_⊙ and a dark matter mass of 1.9 × 10¹⁰ M_⊙ with the same stellar mass resolution of 4 × 10³ M_⊙. Both stars and dark matter follow a NFW profile (presented in Sect. 2.2 of Jiménez-Arranz et al. 2024b).

• G_MW: an MW-mass galaxy without stellar particles, with a total mass of 10¹² M_⊙ distributed in a NFW dark matter halo profile. The mass resolution of its dark matter particles is 5 × 10⁵ M_⊙.

The KRATOS suite includes several groups of simulations featuring an isolated and an interacting group of simulations. The isolated group serves as a control group for this work, including only G_LMC and allowing us to see its secular evolution free from tidal interactions. On the other hand, the interacting group embeds G_LMC within a realistic environment, interacting with the two other systems G_SMC and G_MW. The KRATOS suite includes another group of simulations featuring G_LMC and G_SMC without G_MW. We did not use this group for this work and instead focused on the control isolated group and the most realistic interacting group of simulations including all three systems. In the interacting simulations, G_SMC follows an orbit qualitatively matching the interaction history between the SMC and the LMC (see Sect. 3 of JA24, for a detailed description), while G_MW tidally interacts with G_LMC with a longer range.

In both the isolated and interacting simulation groups, G_LMC is the focus of this study, allowing a comparison of its stellar radial migration across the two scenarios. To that end, we focused on the subsample of the KRATOS simulations presented in Table 1. It includes three isolated and three interacting simulations, each with a Toomre stability parameter of Q=1.0, Q=1.2, or Q=1.5 (K4, K1, and K7 for the isolated simulations, and K6_i, K3_i, and K9_i for the interacting simulations², respectively, for the different Q values). In isolated simulations, the least stable discs (K1 and K4, for instance) are more prone to create non-axisymmetric patterns, and the most stable isolated disc (K7) is not expected to create strong non-axisymmetric patterns. Varying Q therefore allows us to study the impact of non-axisymmetric patterns and stability on radial migration. In summary, the only varying parameters in our work are Q and the fact that the simulations are interacting or isolated. Thus, all masses and orbits (in interacting simulations) are identical, allowing a systematic study.

3 Methods

With the aim of quantifying the amount of stars undergoing radial migration, we first present the method for defining a G_LMC in-plane reference frame and infer their guiding radius, R_g(t), in Sect. 3.1. Then, we present in Sect. 3.2 the methods used to map the non-axisymmetric pattern strengths, the radial migration fluxes (or churning fluxes), and their impact on the metallicity distribution in the disc.

3.1 Defining a G_LMC in-plane reference frame and computing the guiding radius, R_g(t)

In order to quantify properly dynamical effects, it is first required to define a frame of reference, especially in the simulations including interactions between G_LMC+G_SMC+G_MW, G_LMC’s galactic plane gets tilted and sometimes slightly warped due to the interactions (also shown in Jiménez-Arranz et al. 2025a). Since warping is a direct consequence of galactic encounters, its effects cannot be disentangled from those of tidal interactions in this work. Instead, they are treated as a degenerate problem and considered together. However, we corrected for tilting, keeping the original shape of G_LMC and its vertical structure.

For each simulated snapshot, we fitted a plane through the stellar particles of G_LMC’s disc, and defined the z-axis normal to this plane. Then, we defined an orthonormal reference frame with the x – and y-axes within the fitted plane, centred on G_LMC’s centre of stellar mass (as in Schölch et al. 2025). This method is slightly affected by the disc’s tangled warp, but stands as a good compromise to study G_LMC without artificially flattening it, therefore impacting its dynamics. In parallel, we checked the consistency of the method with a plane defined perpendicular to the total angular momentum (L) of G_LMC’s stellar particles (as used in Schölch et al. 2025) and found a difference of less than 3^° between the two planes’ inclinations at all times in all simulations. In the end, we chose to keep the plane fitting method because it appeared more stable especially during the interactions.

The next step consists in reliably estimating the guiding radius, R_g(t), of a star particle in the disc. To do so, we leveraged the method presented in Halle et al. (2015), and used in more recent works such as Khoperskov et al. (2020) and Nagy et al. (2025). This method consists in following the evolution of the galactocentric radius of a particle, R(t) (see the dark blue curves for three random particles in the disc in Fig. 1, top panel), and recording its local maxima and minima. Then, interpolate two curves through the extrema, R_min(t) and R_max(t) (light and dark brown curves, respectively), and estimate R_g(t) as the mean between these curves (pink curves). Because of the simplicity of the method, the R_g(t) estimate remains accurate and stable versus most sudden changes in the gravitational potential, galaxy morphology and tidal interactions with another galaxy. The only requirement is a sufficiently sharp time resolution to properly resolve all epicyclic oscillations of the particles. This condition is fulfilled by the KRATOS suite with a time separation of less than 2.6 Myr between all snapshots (see also the oscillations of the particle closest to the centre in Fig. 1, top panel).

However, this method comes with a two main drawbacks. A first unexpected behaviour of the method can occur when R(t) takes an irregular shape (i.e. non-pseudo-oscillating). These orbits occur mostly during strong tidal interactions, but can also punctually be seen in low Q simulations (i.e. K4, K61, and, to a lesser extent, K1 and K3_i) due to their unstable nature. In that case, the particles can reach local extrema at positions that do not match their epicyclic oscillations. In most orbits, these extrema do not stand as problems to estimate R_g(t) because of the stability of the method. For instance, see the orbit of the middle particle of Fig. 1, top panel: the orbit is irregular during the two last pericentre passages of G_SMC (vertical dashed grey lines) but the pink curve still gives a good estimate of the centre of the orbit. However, the R_g(t) estimate of the particle closer to the galactic centre in the same figure presents unrealistic oscillations between the two last pericentre passages. This effect is over-represented in the very centre of G_LMC’s galactic disc. For this reason and because some slight miscentring of the reference frame can occur, we did not study the radial migration in the region within 1 kpc from G_LMC’s centre.

The second main drawback of the method is its failure to reliably estimate R_min(t), R_max(t) and R_g(t) in short time spans close to the time boundaries of the simulations. The interpolation of R_min(t) and R_max(t) has to extrapolate the existing curve when the next (or previous) extrema is located after (or before) the time boundaries, subsequently affecting R_g(t) in these time spans. Because of this, we avoid drawing strong conclusions on radial migration occurring during the last ∼ 200 Myr of the simulations. The first 1.6 Gyr after the initial conditions are also excluded because some galaxies are still in an unrealistic transient state at the very start of the simulations.

Despite these caveats, this method illustrated itself as the best compromise to study radial migration in G_LMC. We therefore used it to map its galactic disc by estimating the R_g(t) of a large sample of star particles of the disc. Because the KRATOS simulations feature 1.2 × 10⁶ stellar particles in the G_LMC evolving in more than 2000 snapshots, we selected a sample of 1% of all the star particles. We checked the consistency of our results by trying to select different stars and sample sizes, and concluded that 1.2 × 10⁴ stellar particles is a sample size large enough to reliably map G_LMC’s disc.

Fig. 1 Top: visualisation of the method described in Halle et al. (2015) to infer the Rg(t) of the particles: the radial evolution of three random particles in the K9i simulation with respect to time. Rg(t) (pink) is the mean of Rmax(t) (dark brown) and Rmin(t) (light brown), and its value at the start of the simulation is Rg(ti) (horizontal dashed black line). The times of the pericentre passages of GSMC are shown in vertical dashed grey lines, and the evolution of the galactocentric radius of the particles, R(t), is the dark blue curve. The black vertical line represents the t=0 snapshot that is displayed on the bottom panel. Bottom: visualisation of the method for inferring the guiding radius, Rg(t), of a particle at t=0 (same legend as above, with the previous orbit of the particle in white). The shown particle is the particle with the highest radius from the top panel. GS MC is shown in the lower right hand corner, but its particles are not tracked in this work. An animated visualisation is available online.

3.2 Maps of G_LMC

Since its first description in Sellwood & Binney (2002), radial migration has been predominantly attributed to non-axisymmetric patterns and recognised as a driver of radial mixing in the disc (e.g. Minchev & Famaey 2010; Zhang et al. 2025). To illustrate this interplay, we present three types of maps that can be used to survey the interior of G_LMC’s disc.

First, we mapped the non-axisymmetric pattern strength of G_LMC. To that end, we calculated the azimuthal m=2 density Fourier mode (A₂) at each time by integrating it for 30 concentric 0.4 kpc wide rings from the centre of G_LMC to R=12 kpc. We normalised this by the m=0 density mode (A₀), which is the total integrated mass for those rings. Therefore, A₂/A₀ is the fraction of the bisymmetric component over the total mass of each bin.

Second, we mapped net churning fluxes within G_LMC in order to quantify radial migration in the disc (see Halle et al. 2015, for a slightly different implementation of this method). A churning flux (also named radial migration flux) at a radius R is the fraction of disc’s stellar mass whose R_g(t) crosses R between two snapshots, divided by the time separation between these snapshots. The net outward flux is counted positive, and inward flux is counted negative. We present churning fluxes with respect to both t and R by calculating them between every snapshot of the simulations. The high time cadency allows the global tendencies of radial migration during the interactions to be fully resolved. On the other hand, the radial mixing in G_LMC is less well described by the net churning fluxes due to the cancelling effect between the outward and inward migration. This effect, however, is better traced by the metallicity evolution of the disc, whose method of study is described hereafter.

Third, we mapped the metallicity distribution within G_LMC’s disc. Radial migration has long been described as a driver of radial metallicity gradient weakening (e.g. Schönrich & Binney 2009; Kubryk et al. 2013; Grand et al. 2015). We aimed to verify that statement and extend the study to interacting simulations by quantifying the impact of radial migration on the metallicity distribution in the KRATOS suite. Because the KRATOS simulations do not originally consider any form of metallicity, star formation or stellar evolution, some assumptions are necessary.

Firstly, we modelled stellar particles with a fixed metallicity throughout the simulations. This is justified by the fact that we only aim to quantify the impact of churning on a disc’s metallicity profile. Furthermore, we studied migration over a total time span of 3.68 Gyr, which is short enough to consider radial migration as the dominating metallicity-impacting mechanism. Secondly, we considered that particles are born around their guiding radius, R_g(t), at t_i=−3.0 Gyr. The R_g(t) as defined in this work is comparable to an average position of a particle, filtering out its epicyclic oscillations. It is therefore a more reliable estimation of the birth position of the star. Newborn stars are also known for having nearly circular orbits (see Eggen et al. 1962; Yu et al. 2023); hence, it is appropriate to use R_g(t) instead of a random R position in the epicyclic oscillations of the particles for their birth radius. The simulations originally start at t=−4.0 Gyr, and the first pericentre passage of G_SMC happens at t ∼−1.6 Gyr. Therefore, even though some simulations (e.g. K4 and K6_i) are not fully relaxed at t_i, this time stand as a good compromise to use as a stable morphology of G_LMC.

We therefore applied a typical disc galaxy’s radial metallicity gradient of −0.04 dex/kpc to all stellar particles at t_i and followed their mixing with respect to time. This radial metallicity profile is obtained in Tissera et al. (2016) from a catalogue of disc galaxies in a hydrodynamical simulation and is checked consistent with observations. It allows us to write the expression of the initial metallicity profile as follows: $$Z(t_i,R_g)=Z(t_i,0)\times {10}^{-0.04\times R_g}.$$(1) Here Z(t_i, 0) is the initial metallicity at the centre of G_LMC, which is unknown in this work. Because of the unit of the metallicity gradient, R_g(t) must be given in kpc in this formula. We compared the evolving metallicity profile with the initial one at each snapshot t and normalised it by Z(t_i, 0). For testing purposes, we checked different metallicity profile shapes (linear, squared...) and values for the gradient (up to 50% stronger and weaker than observational data), and different values for t_i. None of these changes had a significant impact on the tendencies of our results, as long as the two following requirements are fulfilled. Firstly, the metallicity profile has to be decreasing, which is a very weak assumption confirmed by both theory and observations (e.g. Fu et al. 2013). Secondly, since pericentre passages of G_SMC strongly perturb the morphology and dynamics of G_LMC, t_i has to be selected sufficiently far from such events to better represent a more relaxed, though not fully settled, galactic disc. As an output, this method allows [Z(t, R_g)−Z(t_i, R_g)]/Z(t_i, 0) to be mapped with respect to t and R_g(t). These maps provide a quantification of the impact of tidal interactions on G_LMC’s metallicity profile and, additionally, a direct visualisation of the radial mixing in the galactic disc. Opposite migrations yield zero net churning flux (see the previous paragraph) but still drive radial mixing, as outward migrators enrich the outer disc and inward ones dilute the inner regions. The metallicity maps are therefore complementary to the churning flux maps, enabling us to quantify the flattening of the metallicity profile where net fluxes undergo a cancelling effect.

Additionally, we display on each map the radius of corotation, R_CR(t), with the central bisymmetric pattern. Corotation resonance is known to be able to trap stars in orbits close to R_CR(t), but also to act as a potential barrier, crossed only by stars with highly chaotic or high-energy orbits (e.g. Binney & Tremaine 2011; Brunetti et al. 2011). R_CR therefore is a widely studied parameter with a direct link with churning in isolated galaxy. However, other weaker resonances such as inner and outer Lindblad resonances are beyond the scope of this work. We estimated R_CR by solving the following implicit equation: $$\langle \Omega \rangle (R=R_{\mathrm{CR}},t)-{\Omega }_p\left(t\right)=0.$$(2) Here 〈Ω〉 (R, t) is obtained by averaging the rotation frequencies of all the particles within each of the radial bins used to map A₂/A₀ at each snapshot t. Ω_p is the rotation frequency of the central pattern presented in JA24 calculated with the Dehnen method (see Dehnen et al. 2023, widely used in recent works as Semczuk et al. 2024; Zhang et al. 2024; Jiménez-Arranz & Roca-Fàbrega 2025). A galaxy is considered barred if A₂/A₀ > 0.2 in its central region. If this condition is fulfilled, we refer to the corotation radius with the bar as $R_{\mathrm{bar}}^{\mathrm{CR}}$. On the other hand, when G_LMC does not feature a bar pattern (A₂/A₀<0.2 in the central region), it can still host a non-axisymmetric rotating pattern with which stars might resonate. We hereafter refer to that sort of feature as a weak pattern (w.p.) and to its corotation radius as $R_{\text{w.p.}}^{\text{CR}}$ as displayed in Figs. 2, A.1, and A.2 in metallicity maps only for the sake of readability. For reference, we display the corotation radius for all the snapshots with 0.1<A₂/A₀<0.2. The condition A₂/A₀ > 0.1 is the minimum threshold that allows Ω_p and therefore $R_{\text{w.p.}}^{\text{CR}}$ to be accurately computed.

4 Results

Our main results describe the interplay between nonaxisymmetric structures (Sect. 4.1), radial migration (Sect. 4.2), metallicity distribution (Sect. 4.3) and corotation resonance (Sect. 4.4) under the influence of tidal interactions within G_LMC in the KRATOS simulations displayed on Table 1.

We display our results in maps surveying G_LMC’s interior with respect to time t and galactocentric radius R (or guiding radius, R_g(t)); see Figs. 2 for K7 and K9_i, A.1 for K4 and K6_i, and A.2 for K1 and K3_i. These maps were obtained by applying the method described in Sect. 3.2. Each figure includes three maps of an interacting simulation (the three upper panels) and three maps of its isolated counterpart (the three lower panels). For each of the simulations, the first map is the A₂/A₀ map in the galactic disc, displaying the non-axisymmetric patterns strength. The second map displays the churning fluxes (or radial migration fluxes) in the galaxy. The last map shows their impact on the initial metallicity distribution and displays the radius of corotation with the bar as well as with a weak non-axisymmetric pattern, whereas only corotation with the bar is displayed on the other maps. Additionally, all interacting maps feature vertical dashed white lines and white stars matching the pericentre passages of G_SMC and its impact parameters on G_LMC’s disc, respectively.

Fig. 2 Maps of A2/A0 Fourier mode strength (top panel), difference between the churning fluxes outwards and inwards (central panel, counted positive outwards), and normalised metallicity evolution (bottom panel) with respect to time and galactocentric radius, R, or guiding radius, Rg(t), for the K9i (interacting Q=1.5 simulation, top three panels) and K7 simulations (isolated Q=1.5 simulation, bottom three panels). The grey area masks the first 0.6 Gyr displayed, when the disc is not fully relaxed yet. For the K9i maps, the times of the pericentre passages of GSMC are shown in dashed white lines, and its impact parameters are the white stars. The corotation radius is plotted in green and yellow for resonance with a bar or a weak central pattern, respectively.

4.1 Non-axisymmetric patterns

The correlation between the presence of non-axisymmetric structures and the tidal interactions was already established in the KRATOS presentation paper (JA24). However, their results were obtained by calculating the median of A₂/A₀ on the entire bar region. Contrastingly, our work extends the study to the outer disc and visualises A₂/A₀ on a much thinner binning, providing more details. We present our results for both isolated G_LMC and interacting systems.

We first present the isolated simulations: K7, K4 and K1 (see the three bottom panels of Figs. 2, A.1 and A.2, respectively). Being free from tidal forces with foreign systems, G_LMC only depends on its secular evolution to form non-axisymmetric patterns. None of these simulations show strong spiral arms at the relaxed state (see Fig. 3 of JA24, for a visualisation of the galaxies at t=0). Nonetheless, because their difference in the Toomre parameter, Q, they present dissimilar morphological evolutions. The most stable of them (K7, with Q=1.5) is unable to form a strong bar within the total simulated 4.68 Gyr, only a very weak central pattern that is slowly growing from t ∼−1.0 Gyr until the end of the simulation. Additionally, K7 does not form any spiral arm during its evolution, resulting in a very dark A₂/A₀ map (i.e. with a very low amplitude) in Fig. 2. Contrastingly, the most unstable isolated simulation (K4, with Q=1.0) is able to form, maintain and grow a strong bar pattern longer than 2 kpc from t ≳−2 Gyr (also described in JA24). Additionally, it sparsely showcases small and local spiral arms (i.e. not covering the entire disc, some stains can be seen in K4’s A₂/A₀ map of Fig. A.1). In between both simulations, K1(Q=1.2) develops a weak pattern in its centre from t ≳−1.5 Gyr growing into a strong but very variable bar at t ≳ 1.0 Gyr. This variability can be seen in the low R range of Fig. A.2’s A₂/A₀ map, in particular at the end of the simulation where the upper limit of the strong A₂/A₀ region oscillates. Furthermore, this simulation presents spiral patterns similar to K4’s but with a weaker strength.

Secondly, we present the interacting simulations: K9_i, K6_i, and K3_i (see the three top panels of Figs. 2, A.1, and A.2, respectively). Despite sharing the same Toomre parameter (Q) as their isolated counterparts (K7, K4, and K1, respectively), the interacting galaxies exhibit marked morphological differences. This is triggered by tidal interactions and may include both the formation and destruction of bar patterns and the creation of strong global spirals with a large variety of shapes. This has already been explicitly presented in JA24; the novelty of this work is that we can see the evolution of their A₂/A₀’s radial distribution. The spiral pattern evolution of G_LMC is similar in all interacting simulations. Namely, G_LMC’s outer disc takes a global spiral shape after the first pericentre passage of G_SMC. The high A₂/A₀ area corresponding to this spiral shape appears to propagate outwards in all the interacting simulation’s A₂/A₀ maps. Equivalently, the low A₂/A₀ areas follow the same direction. We note that the first interaction is also the longest-range pericentre passage of G_LMC, and the only one with an impact parameter in the outer disc. Furthermore, these spiral arms remain short-lived and variable, being susceptible to be destroyed and created again. They also are much weakened after the second galactic interaction.

Although these previously listed phenomena are similar for all the interacting simulations, G_LMC’s bar pattern evolves differently depending on the Toomre parameter, Q. For instance, in K9_i(Q=1.5), G_LMC’s bar presents a very singular behaviour. It is formed during the first pericentre passage of G_SMC (the first vertical dashed white line of Fig. 2), initially shorter than 1 kpc and circled by a ring. This bar is not well detected by the Dehnen method but keeps its shape until the second galactic encounter. During this encounter, the bar almost stops rotating for around 100 Myr (see Jiménez-Arranz & Roca-Fàbrega 2025), shifts from the centre and gets rid of its ring. After this encounter, it will relax and continue growing until the end of the simulation. On the other hand, K6_i(Q=1; see Fig. A.1), as its isolated counterpart K4, starts forming a strong steady bar at t ∼−2 Gyr (i.e. before the first interaction). This bar grows during the entire simulation, ultimately being stronger and larger than K4’s bar. However, because of the tidal interactions, K6_i’s bar is more variable, impacting its corotation radius (see Sect. 4.4). Finally, K3_i’s bar (Q=1.2; see Fig. A.2) presents a similar behaviour as K6_i’s. The main difference is that the bar creation starts during the first pericentre passage (i.e. after K6_i’s bar). It therefore starts slightly before K1’s bar (being K3_i’s isolated counterpart).

The morphology of G_LMC at a given time and radius is inherently connected to the motion of its stars. Consequently, A₂/A₀ maps provide complementary information to churning fluxes and metallicity distribution maps presented in the next sections.

4.2 Radial migration

The general tendencies of radial migration are described by the churning fluxes in G_LMC’s disc (see Sect. 3.2). These fluxes are represented in the second panels of each group of maps in Figs. 2, A.1, and A.2. Although some simulations can spontaneously show some exceptions (see the end of this subsection), the general tendencies are similar within each group of simulations: isolated and interacting.

The isolated simulations do not display strong global churning fluxes, usually not exceeding ∼ 15−20% of disc mass per gigayear in any direction. There are still non-zero fluxes in the galactic discs, especially in the least stable simulations (K4 and K1; see Figs. A.1 and A.2, respectively). These fluxes appear to be oscillating with time in the centre of the galaxy – especially at the end of the simulations. On the other hand, the churning fluxes are similar to a diffuse noise in the outer disc, with some spontaneous larger patterns – in particular when the galaxy is not fully relaxed. These features are also present in the interacting simulations, hence the relevance of this isolated control group.

The interacting simulations’ churning fluxes showcase a strong global wave pattern reaching a first minimum at a time perfectly matching the pericentre passages of G_SMC. The churning fluxes decrease before the pericentre passages, then increase back with a similar amplitude, peaking around the impact parameter. Their values reach ∼ 50% of the disc mass per gigayear in both directions during the strongest interaction: the second pericentre passage. This radial migration wave continues to oscillate with a shorter wavelength after the pericentre passages. Although it is most notably visible in K9_i’s second interaction because of the high stability of G_LMC in this simulation (see the second panel of Fig. 2), this wave pattern exists and is visible in all the interacting simulations. It covers the entire galactic disc and oscillates faster near the galactic centre than in the outer disc (i.e. the wavelength is shorter close to the centre, hence the appearance of churning fluxes in diagonal bands on the maps). This statement generally remains valid when the impact parameter of G_SMC is in the outer disc (see for instance the first pericentre passage in K9_i in Fig. 2).

However, some exceptions can be seen in the simulations. One of the most prominent can be seen in K7 at t ∼−1.1 Gyr (see the second last panel of Fig. 2), where a weak wave pattern similar to the ones happening during tidal interactions is triggered. We discuss this phenomenon in Sects. 4.4 and 5.1. This behaviour may also appear in the other simulations, though it is less discernible due its very low amplitude versus the fluxes created by their instabilities.

4.3 Metallicities

As mentioned in Sect. 3.2, the study of the radial metallicity distribution allows the radial migration in G_LMC’s disc to be linked with an observable and provides complementary information on the radial mixing induced by churning. We again present the results for each group of simulations.

In the isolated group, churning tends to enhance the metallicity in the outer disc, and reduce it closer to the centre. It therefore directly weakens the metallicity gradient with time. This effect is much amplified for low Q simulations because of their instability (see K4 and K1 in third panels of Figs. A.1 and A.2, respectively). This metallicity difference can reach ∼ 7−8% of the initial maximum metallicity Z(t_i, 0) in K4 and K1, and ∼ 3% of Z(t_i, 0) in K7. Moreover, a lower Toomre parameter (Q) correlates with a sharper demarcation between regions of enhanced and reduced metallicity. We hereafter refer to these regions of unclear demarcation as ‘mixing zones’ (see for instance the 3<R_g<6 kpc region in almost all of the K1 simulation in the last panel of Fig. A.2). Finally, lower values of Q are associated with a more rapid development of these effects.

The aforementioned effects are also applicable to the interacting group. However, they present some differences due to the tidal interactions. The most visible of which is a metallicity drop of ∼ 3−5% of Z(t_i, 0) during the two first pericentre passages of G_SMC – therefore matching with the minimum of the wave presented in Sect. 4.2. This drop is most distinctly visible in K9_i (see the third panel of Fig. 2), but also appears in K6_i and K3_i (see the third panels of Figs. A.1 and A.2, respectively). After this drop, K9_i features a sudden metallicity increase that propagates outwards across the disc, before decreasing again. This effect is discussed further in Sect. 5. Finally, tidal interactions tend to strengthen the radial mixing presented for isolated simulations. For instance, the radial metallicity gradient is weaker at the end of K9_i than it is for its isolated counterpart K7 at the same time.

We emphasise that these changes in metallicity do not describe blurring but solely churning. They depend only on the guiding radius, R_g(t), and are therefore a metric to quantify the radial migration in the disc. To conclude, the effects presented in this subsection can be seen as the typical impacts of radial migration – with or without tidal interactions – on the metallicity distribution within the galaxy.

4.4 Influence of the corotation resonance

As noted in Sect. 3.2, the corotation radius, R_CR, namely, the radius at which the corotation resonance lies, is directly linked with the pattern speed and strength of the central nonaxisymmetric pattern (a bar when A₂/A₀>0.2 in that region; otherwise, there is a weaker pattern). Thus, the results presented in this section are to be linked with those mentioned in Sect. 4.1, where we show that the non-axisymmetric patterns are highly variable in most of the KRATOS simulations. As a consequence, the corotation radius is also very variable, and is largely affected by tidal interactions. It is, however, still a good indicator of a radius acting both like a trap and a barrier, especially for the isolated simulations where G_LMC hosts a bar.

Due to the high variability of R_CR in interacting simulations, we focused on the results obtained in the isolated simulations. Our main results can be observed in K4 and K1, because of the strength and stability of their bars. First, in K4, R_CR is fairly stable and almost always lies at the demarcation between enhanced and reduced metallicities (i.e. in the small mixing zone; see Fig. A.1). Second, in K1, R_CR is more variable, and is slightly higher than the radius of the mixing zone (see Fig. A.2). In both of these simulations, the amplitude of the variations in R_CR matches the radial width of the mixing zone. This correlation between metallicities and R_CR is discussed in Sect. 5. Furthermore, we state in Sect. 4.2 that the churning fluxes behaved differently closer to the centre than in the outer disc; rapidly oscillating with time near the centre and showing larger noisy patterns in the outer disc. We find that in both K4 and K1, R_CR appears to be the demarcation between these two behaviours. These two statements could also be applied to K6_i, even though the computation of R_CR is complexioned by the interactions.

Even though K7 does never form a strong bar, R_CR suddenly increases at t ∼−0.9 Gyr, suggesting that a weak non-axisymmetric pattern is growing. The time of this growth almost matches the triggering of the weak radial migration wave mentioned in Sect. 4.2. For reference, we note that in K9, R_CR ∼ 11 kpc between t ∼ 0 and t ∼ 0.5 Gyr, due to G_LMC’s slowly rotating bar. This high R_CR correlates with the sudden decrease in metallicity in the same time span at R_g ∼ 8.5 kpc (see Fig. 2). The link of causality between a high R_CR and a sudden drop of metallicity below it remains uncertain. We find a similar behaviour, though in a lower extent, in K3_i and K7 (at t ∼−0.7 Gyr in the third panel of Fig. A.2, and at t ∼ 0.5 Gyr in the last panel of Fig. 2, respectively).

5 Discussion

5.1 Isolated simulations

Isolated simulations (K1, K4, and K7) provide a control group to study tidal interactions. Yet, they display intrinsic dynamics that merit analysis and allow comparison with prior works.

A lower Toomre parameter, Q (Toomre 1964) increases instabilities and non-axisymmetric structure in G_LMC, and in all isolated simulations radial migration correlates with these patterns (which is consistent with e.g. Athanassoula & Sellwood 1986; Fujii et al. 2011). This highlights a clear connection observed in all isolated simulations: radial migration strongly correlates with the presence of non-axisymmetric patterns. Thus, we observe that a lower Q – therefore a stronger non-axisymmetric component – leads to increased churning fluxes (see Sect. 4.2) accompanied by more rapid and intense radial mixing (i.e. weakening of the metallicity gradient; see Sect. 4.3). This tendency is at the basis of the concept of radial migration and has always been identified as its main cause in isolated galaxies. Our results therefore are strongly consistent with those found in the existing literature of the field (e.g. Sellwood & Binney 2002; Schönrich & Binney 2009; Iles et al. 2024; Marques et al. 2025). Furthermore, all simulations exhibit noise-like, incoherent churning throughout (including initial relaxation; grey area in Figs. 2, A.1, and A.2). Given our metallicity definition (dependent on R_g(t_i); see Sect. 3.2), this supports radial migration as the primary driver of metallicity smoothing, especially in the least stable galaxies.

As described in Sect. 4.2, the churning fluxes behave differently in the outer disc than close to the centre of G_LMC: the fluxes show almost coherent oscillations near the galactic centre, whereas they display noise-like patterns varying less rapidly in the outer disc. This feature is especially visible in barred systems, partly because the fluxes themselves are stronger due to their less stable nature (see for example the second last maps of Figs. A.1 and A.2). For instance, the period of these oscillations after t ∼ 0.5 Gyr near the centre of G_LMC’s centre in K1 (∼ 4 full oscillations per gigayear) matches double of the bar rotation period in this simulation. Furthermore, these oscillations are clearer when the bar is stronger (i.e. when R_CR is displayed with a green line in Fig. A.2). To summarise, these oscillations occur each time a strong bar undergoes half a rotation in K1, they therefore may be a consequence of the local changes in the gravitational potential created by the bar’s motion. Thus, it strengthens the link between non-axisymmetric patterns and radial migration. This reasoning can be applied to the other simulations but is harder to quantify. K4’s churning fluxes map features numerous local patterns that may be triggered by instabilities due to its critical value of Q, slightly hiding and perturbing these central oscillations. Regarding K7, it only features very low fluxes, and quantifying the period of these oscillations would be challenging. We still mention that small oscillating patterns are visible in K4 and, to a lesser extent, in K7 near the end of the simulation when the weak non-axisymmetric pattern is formed in its centre. Additionally, it appears in both K4 and K1 that the boundary between the central oscillation-dominated zone and the outer disc is R_CR. This is consistent with the fact that R_CR is known to act as a potential barrier (e.g. Binney & Tremaine 2011; Brunetti et al. 2011), it therefore also serves as a barrier between the different regimes of radial migration.

There is a strong correlation between the radial width of the mixing zone (as defined in Sect. 3.2) and the amplitude of R_CR’s rapid pseudo-oscillations. This effect is particularly visible in K4 and K1, which feature a stable R_CR with a very small mixing zone and an oscillating R_CR with a wider mixing zone, respectively. This correlation could be explained by the fact that the corotation radius tends to trap particles in its radial range. Thus, when R_CR oscillates with some amplitude, stars trapped in corotation mix together, creating this mixing zone, whose width matches this amplitude. Although the most reliable occurrences of this correlation happen when R_CR is calculated with a strong bar pattern, we note that a similar effect is visible for a weaker pattern in K1 before t ∼−1 Gyr. In this time span, K1’s R_CR is strongly variable, correlating with a mixing zone covering almost the entire galactic disc with metallicity enhancements near the centre and reductions in the outer disc (see the last panel of Fig. A.2). The mixing zone width tends to reduce when the bar grows stronger and when R_CR is less strongly varying. A similar tendency is visible in K7 after t ∼−0.5 Gyr, but because the mixing zone is less distinguishable and the central pattern is very weak in this simulation, we refrain from drawing strong conclusions on this effect in K7.

As mentioned in Sect. 4.2, K7’s radial migration hosts a major peculiarity at t ∼−1.1 Gyr (see K7’s churning flux map in Fig. 2): A small radial migration wave seems to be triggered at that time (that kind of wave is discussed more extensively in Sect. 5.2), correlating with a weak metallicity drop (<1% of Z(t_i, 0)) in the inner disc. A few tenths of a megayear after the wave is triggered, a weak non-axisymmetric central pattern is formed in G_LMC’s centre. These two events are likely linked because bars are created by instabilities known for being a cause of radial migration (e.g. Di Matteo et al. 2014; Haywood et al. 2024). Moreover, bar formation is often considered as the strongest cause of churning in isolated galaxies (e.g. Halle et al. 2015). The churning flux wave triggered shortly before the appearance of a weak bar-like pattern in K7 is therefore in agreement with the scientific literature preceding this work, and we conclude that it may be triggered by the bar forming instability. It is also possible that K7 could develop a stronger bar over a longer timescale than that covered in the simulation, as we did not observe its weak non-axisymmetric pattern ceasing to grow. Similar churning patterns are slightly visible in K4 and K1 at similar times, they are less visible because of their higher background churning flux intensity. It is therefore not obvious to know if the bar formation is at their origin too, but we do not exclude this possibility. We also note the bar buckling instability, which arises after bar formation, reshapes the vertical structure into an ‘X’ or peanut shape (Raha et al. 1991; Beraldo e Silva et al. 2023). This instability can alter bar strength and pattern speed, shifting resonances and releasing stars previously trapped in resonant orbits, thereby punctually impacting radial migration (e.g. Halle et al. 2018; Khoperskov et al. 2020). It is therefore another instability expected to generate coherent fluxes in K4 and K1.

Overall, our results with isolated simulations reproduce established radial-migration behaviours, providing a baseline for the interacting cases and reinforcing the robustness of our subsequent conclusions.

5.2 Interacting simulations

The interacting group of simulations (K3_i, K6_i, and K9_i) hosts the main results of this work and therefore deserve an extensive discussion. In this section, we focus on the consequences of the tidal interactions. We expect that most of the results discussed for the isolated simulations (see Sect. 5.1) are also applicable to their interacting equivalents, although usually not directly visible because hidden behind the stronger effects we discuss in this section. Thus, the goal of this section is to provide an interpretation of the effects of tidal interactions on galaxy morphology, radial migration, radial metallicity distribution, and the interplay between all of them.

In all interacting simulations, the first tidal interaction creates a spiral structure remaining ubiquitous in G_LMC until the next pericentre passage of G_SMC. Two main reasons are to be linked with the creation of this structure. First, the impact parameter (b_SMC ∼ 8.0 kpc) means that G_SMC is located near the outer disc of G_LMC during this interaction. Therefore, the spiral structure is created where the tidal interaction is the strongest. Second, this is the longest-ranged pericentre passage of G_SMC with a total distance of 11.8 kpc between the galaxies and a vertical distance of z_peri=−8.7 kpc whereas z_peri ∼ 0 for the two later passages. To summarise, the only interaction that can be seen creating a strong large-scale spiral structure is the longest-range interaction and has an impact parameter in the outer disc. We also note that a similar churning flux wave pattern is triggered during both first interactions, but only the first is seen creating global spiral arms. However, during the first interaction, the churning wave is triggered before the formation of spiral arms suggesting that the non-axisymmetric pattern is a consequence of radial migration (see the radial migration fluxes in green and orange around the first vertical dashed white line, both appearing before the orange overdensity of A₂/A₀ in Figs. 2, A.1, and A.2). In other words, the creation of large-scale spiral arms in the KRATOS simulations is not merely a consequence of a change in the stars’ position but also a radial shift in their guiding centre. This shift is consistent with the pre-established link between radial migration and transient spiral arms (e.g. Sellwood & Binney 2002; Solway et al. 2012; Marques et al. 2025).

A radial migration wave triggered by the two first interactions between G_LMC and G_SMC is visible in all interacting simulations’ churning flux maps (see Figs. 2, A.1, and A.2). These waves systematically show inward fluxes, then a minimum matching with the pericentre passage and outward fluxes, eventually showing more oscillations in the second and stronger tidal interaction. As mentioned in Sect. 4.3, all pericentre passages match with a metallicity drop of a few percent of Z(t_i, 0) in almost all the galactic disc. This drop is a direct consequence of the radial migration wave: The initial metallicity profile decreases with R_g(t). Therefore, when the net churning fluxes describe a general inwards migration, regions that were initially more metal-rich become populated by particles that were initially less metal-rich, leading to a drop in metallicity in those regions. Furthermore, during the second interaction in K9, metallicity increases again following this drop and mirrors the shape of the churning flux wave – peaking after the outward fluxes (see the third panel of Fig. 2). This reversal further reinforces the causal link between radial migration and variations in metallicity. We also mention that these metallicity drops and increases are just redistributions of the metallicity in the disc, because star formation and evolution were overlooked (see Sect. 3.2), the total amount of metals in the model remains unchanged.

In all of the interacting simulations, the churning wave starts with inward fluxes and is followed by a bouncing effect that creates a slightly stronger outward migration. We expect the order of events to be due to the nature of tidal interactions. Indeed, in the KRATOS simulations, G_SMC follows an almost polar orbit around G_LMC. This implies that the tidal forces apply in-plane components towards the centre of G_LMC, and normal components tending to thicken it. Thus, the in-plane component is in agreement with the first direction of the churning fluxes and the normal component explains the increase in G_LMC’s disc scale height after the interactions (see Fig. 6 of JA24). We attribute the subsequent outward fluxes to the relaxation of particles that gained energy during the interaction. In fact, the disc height reduces right after the pericentre passages-especially the first-thus meaning that stars lose some of the vertical kinetic energy they gained after the interaction. Because the stars are dragged towards G_LMC’s centre, therefore with a radial component, some of the vertical kinetic energy of the stars can be converted into azimuthal kinetic energy; this increases their R_g(t), driving the outward churning fluxes. The in-depth link between radial migration and vertical motion remains unclear and is out of the scope of this work. While some works show that the strongest migrators are the particles with the least vertical motion (e.g. Vera-Ciro et al. 2014; Vera-Ciro et al. 2016), more recent studies present other regimes where the vertical motion has less or no impact on the radial migration (e.g. Halle et al. 2018; Mikkola et al. 2020). Finally, we mention that the churning fluxes oscillate increasingly faster over time towards the galactic centre. This may be attributed to the shorter dynamical timescales governing orbital evolution near the galactic centre.

The estimation of R_CR in interacting simulations is challenging; we therefore treat our results with caution. The corotation radius exhibits abrupt changes and even discontinuities in the interacting simulations, especially near the pericentre passages of G_SMC. We attribute these irregularities to the unpredictable nature of morphological changes in the galaxy, the sudden changes in the bar pattern speed, the off-centring of the bar, and the asymmetries induced by tidal interactions (see Schölch et al. 2025). All of these parameters may interfere with the measurement of the bar pattern speed and the estimation of G_LMC’s rotation curve, both of which appearing in R_CR’s calculation (see Sect. 3.2). However, we still notice some local correlations between the resonance with a strong bar and local patterns in the metallicity maps. As described in Sect. 4.4, K3_i and K9_i exhibit local metallicity drops in the outer disc (see the blue horizontal plumes in the third panels of Fig. 2 at t ∼ 0.3 Gyr and Fig. A.2 at t ∼−1.0 Gyr). Both of these oddities happen at ∼ t_SMC+500 Myrs (i.e. after a G_LMC−G_SMC pericentre passage) and correlate with the highest values of R_CR with a strong bar in their corresponding host simulations. The absence of both metallicity drop plumes in the outer disc and very high R_CR values in K6_i further supports our interpretation of this correlation. Although the physical mechanism underlying this interplay remains uncertain, we recall that R_CR is able to drag and trap particles in its orbital range. Because both of the described local metallicity drops start during a tidal interaction, itself creating a larger scale metallicity drop, it is plausible that R_CR may contribute to maintaining the particles’ positions in the outer disc as they were during the radial migration wave. Moreover, we do not draw strong conclusions on the effect of resonance with weak patterns in interacting simulations because of their high variability, and therefore only show their corotation radius in the metallicity maps for reference. Ultimately, we considered only the corotation with the bar and not the spiral arms as these lie outside the scope of this work, alongside the bar’s inner and outer Lindblad resonances.

Regarding the long-term outcome of the tidal interactions, we note from the metallicity maps that radial mixing is stronger at the end of the interacting simulations compared to their isolated counterparts. Although this effect is most evident in the K7 and K9_i pair (see Fig. 2), that tendency exists in all pairs of simulations with at least a 1% stronger metallicity increase in the interacting G_LMC’s outer discs. This result is to be expected: tidal interactions have long been known for inducing dynamical disc heating (e.g. Gerhard & Michael Fall 1983; Gnedin 2003). Thus, a disc that has experienced tidal interactions contains increased random motion, leading to enhanced radial mixing. Eventually, it implies a weaker metallicity gradient.

The results we present for the impact of radial migration on the metallicity distribution may extend to other observable quantities linked to the stars and presenting a gradient. Given our approach to estimating the weakening of the metallicity gradient (see Sect. 3.2), our findings can be applied to any quantity that shows a decreasing (or, conversely, increasing – yielding mirrored results) trend with guiding radius R_g(t). It also requires the assumption that it remains constant over time for a given star to hold reasonably well. Thus, because of radial migration, we expect a smoothing of the radial stellar age profile and the colour gradients in the galaxy with similar trends than what we observe for metallicity. Nevertheless, a more realistic quantification of these effects remains only a future prospect of this work and could require more sophisticated simulations.

In the end, the interacting simulations show very similar results during G_SMC’s pericentre passages with a similar wave triggered in K9_i, K6_i, and K3_i during all interactions. Thus, we conclude that close-range tidal interactions can reshape the galaxy and dominate the dynamics of LMC-like galaxies from Q=1.0 to Q=1.5. Currently, models studying galactic dynamics and radial migration usually consider isolated galaxies (e.g Halle et al. 2015; Halle et al. 2018; Nagy et al. 2025; Marques et al. 2025) while galactic encounters are known to be common in our Universe (e.g. Barnes 2001; Unavane et al. 1996; Antoja et al. 2018; Merrow et al. 2024). We present results showing that they can be the dominant radial migration-driving mechanism in galactic discs, being able to generate churning fluxes of ∼ 40% of disc mass per gigayear. Thus, our work suggests that tidal interactions cannot be neglected to study galactic evolution in a realistic environment (see Williamson et al. 2016; Carr et al. 2022, for similar conclusions using other methods and scales).

5.3 Limitations and future prospects

We remain far from fully understanding the principles underlying the impact of tidal interactions on galactic discs. This work overlooks some key mechanisms that would have to be implemented for a more realistic study of radial migration. These are as follows:

1. The most obvious caveat of the KRATOS suite is the absence of hydrodynamics, and hence stellar formation and evolution, sacrificed for a higher resolution. Because stars are usually formed in nearly circular orbits (e.g. Eggen et al. 1962; Yu et al. 2023), they would populate almost nonexistent orbital types in KRATOS and may have an impact on the overall migration. Furthermore, due to the lower rotation rate of the gas component, it would apply some damping to the stars’ motion that was neglected in this work. Most importantly, incorporating star formation and stellar evolution would provide a far more accurate estimation of key observables of the galaxy. It would free this work from many of the assumptions we relied on to derive an estimate of the metallicity profile’s evolution. Thus, the authors of this work are already developing new simulations that will study the effect of adding hydrodynamics, star formation and evolution physics, as a continuation of this project aiming to achieve more sophisticated and comprehensive simulations;

2. As discussed in Sect. 3.2, all of our maps and most of our results are represented with respect to R or R_g(t) and t. However, it is known that both LMC and its simulated counterparts exhibit asymmetries (e.g. Gaia Collaboration 2021; Jiménez-Arranz et al. 2023a; Jiménez-Arranz et al. 2025a) that we overlooked in this work. Thus, a finer description of migration accounting for variations both with the radius and the azimuthal angle across the disc, would be highly beneficial. This further study might be able to show that the impact position in the galactic plane correlates – or not with the highest local churning fluxes. Some work has been done to quantify the asymmetries in the LMC disc (Schölch et al. 2025) and a clear future prospect of this study would be to combine their results and methods with ours, enabling an investigation of the correlations between radial migration and asymmetries;

3. We only presented results of the estimation of the bar or central pattern’s corotation radius, therefore lacking an estimation of the Lindblad resonances and resonances with the spiral arms. Determining the precise locations of the resonances under the influence of tidal interactions is not straightforward, yet it is an aspect that more advanced hypothetical future studies should take into consideration. However, such an endeavour may require simulating a larger galaxy than G_LMC as the corotation resonance already reach the outer disc in some interacting simulations, therefore increasing the computational cost for the same resolution. This issue would be more pronounced when considering the outer Lindblad resonance or corotation with the spiral arms, with the risk of these radii exceeding the size of the galaxy itself;

4. Because the only varying input parameters of this work are the Toomre parameter, Q, and the presence of interactions, this study could be extended to additional simulations exploring the impact of other parameters. For instance, how is radial migration influenced by changes in the mass or shape of the dark matter halo? How would our results be affected by changes in the orbit of G_SMC ? Would the presence of a supermassive black hole at the galactic centre alter radial migration or reduce the impact of tidal interactions? Do tidal interactions produce similar effects in elliptical or irregular galaxies? Such further research could offer insights into regimes less dominated by tidal interactions than those explored in the present study. In particular, the impact of G_SMC’s orbit and G_LMC’s halo could be investigated using other KRATOS simulations as a natural continuation of this work.

By leveraging simulations and observations, future studies could better characterise tidally induced radial migration waves, which drive inward migration before pericentre passages and outward migration afterwards. They may lead to enhanced surface brightness in the inner disc and reduced brightness in the outer disc during pericentre passages, along with the metallicity drop of ∼ 3−5% of the maximum metallicity in the inner disc (see Figs. 2, A.1, and A.2) and similar changes in stellar age distributions, sketching a lasting imprint on the disc. In the case of G_LMC, these effects peak for ∼ 200 Myr around the pericentre passages and last for up to 1 Gyr. If they are systematic and reliably modelled, they would provide a method for estimating the mass and orbit of a satellite galaxy that interacted with the main system, since the strength of the wave depends on these parameters. By analysing the surface brightness, metallicity, and stellar age distributions in simulations of isolated galaxies and comparing them with the observations, we can use the radial migration wave to infer independent constraints on satellite galaxy properties (e.g. its mass and orbit) and demonstrate the role of harassment as a satellite quenching mechanism (e.g. Rodríguez-Cardoso et al. 2025). Additionally, we note longer-lasting impacts of tidally induced migration in K9_i(Q=1.5) with a ∼ 4% metallicity shift relative to K7, lowering it in the outer disc and enhancing it in the inner disc (see Fig. 2).

6 Conclusions

In this work, we leveraged the high-resolution simulation suite KRATOS (JA24) to study radial migration under the influence of tidal interactions in LMC-like galaxies. These pure N-body simulations of the Clouds include a control group of isolated simulations of G_LMC, an LMC-like galaxy and an interacting group where G_LMC interacts with the SMC-mass system G_SMC and the MW-mass system G_SMC. We focused on a subsample of 6 of the 28 KRATOS simulations (see Table 1), chosen to be suited for a systematic comparison of the effects of the tidal interactions (see Sect. 2). In order to quantify the impact of these interactions, we mapped the non-axisymmetric pattern strength, the net churning fluxes, and the impact of radial migration on the radial metallicity distribution in G_LMC in three interacting simulations and their three isolated counterparts. In this context, our key results can be summarised as follows:

1. We present a description of a global radial migration wave pattern triggered by tidal interactions, along with comparative mapping of the effects of galactic encounters in LMC-like systems. Its amplitude can reach up to ∼ 40% of disc mass per gigayear at times close to the G_LMC−G_SMC pericentre passages;

2. These oscillating net churning fluxes come with a metallicity drop in the galactic disc that matches the pericentre passages of the satellite, amounting to ∼ 3−5% of the maximum value in the galaxy’s radial metallicity distribution prior to the interaction;

3. During the interactions, all simulations behave similarly in the studied range of Toomre parameters (1.0≤ Q ≤ 1.5), suggesting that strong tidal interactions tend to govern the radial migration in galactic discs;

4. Radial migration tends to behave differently in the region below corotation – which is dominated by the bar – compared to the above-corotation resonance.

Using both simulations and observations, future studies could better characterise tidally induced radial migration waves, which drive inward migration before pericentre passages and outward migration afterwards. They may lead to a higher surface brightness in the inner disc and reduced brightness in the outer disc during pericentre passages, along with the metallicity drop we see in the inner disc and similar changes in stellar age distributions. If these biases are systematic and reliably modelled, they could provide a method for estimating the mass and distance of a satellite galaxy, since the strength of the wave depends on these parameters. With sufficiently precise observations and accurate models of surface brightness, metallicity, and stellar age distributions for isolated galaxies, the measured biases could yield independent constraints on satellite galaxy properties.

Within KRATOS’s high-resolution framework, we uncovered oscillations of radial migration, estimated their imprint on the radial metallicity distribution in LMC-like galaxies, and anticipated the biases they could induce on observational properties such as surface brightness and stellar age distribution. Yet, our analysis is limited to the migration of immutable stars, overlooking their rich environment, asymmetries, evolutions, births, and deaths. As a result, future efforts must integrate these missing elements to reveal novel insights into the complexity of galactic dynamics, relying on a combination of simulations and observations. The kinematic data in Gaia Data Release 4, expected December 2026, will have improved precision for nearly two billion stars. In parallel, high-resolution hydrodynamical simulations may improve the theoretical modelling of stellar kinematics. The integration of these approaches will advance our understanding of the dynamical processes that govern disc galaxies.

Data availability

The movie associated to Fig. 1 is available at https://www.aanda.org

Acknowledgements

This work is based on the Master’s thesis of DH (2025, Université de Toulouse) conducted at Lund University’s Division of Astrophysics, under the supervision of OJA and SRF. The computer “Brigit” in the “Centro de Datos de la Universidad Complutense de Madrid” was used to re-run KRATOS simulations with an increased time cadency. OJA acknowledges funding from “Swedish National Space Agency 2023-00154 David Hobbs The GaiaNIR Mission” and “Swedish National Space Agency 2023-00137 David Hobbs The Extended Gaia Mission”.

References

Appendix A Additional maps

Fig. A.1 Same as Fig. 2 but for K6i (top three panels, Q=1 interacting simulation) and K4 (bottom three panels, Q=1 isolated simulation).

Fig. A.2 Same as Figs. 2 and A.1 but for K3i (top three panels, Q=1.2 interacting simulation) and K1 (bottom three panels, Q=1.2 isolated simulation).

All Tables

All Figures

Fig. 1 Top: visualisation of the method described in Halle et al. (2015) to infer the Rg(t) of the particles: the radial evolution of three random particles in the K9i simulation with respect to time. Rg(t) (pink) is the mean of Rmax(t) (dark brown) and Rmin(t) (light brown), and its value at the start of the simulation is Rg(ti) (horizontal dashed black line). The times of the pericentre passages of GSMC are shown in vertical dashed grey lines, and the evolution of the galactocentric radius of the particles, R(t), is the dark blue curve. The black vertical line represents the t=0 snapshot that is displayed on the bottom panel. Bottom: visualisation of the method for inferring the guiding radius, Rg(t), of a particle at t=0 (same legend as above, with the previous orbit of the particle in white). The shown particle is the particle with the highest radius from the top panel. GS MC is shown in the lower right hand corner, but its particles are not tracked in this work. An animated visualisation is available online.

In the text

Fig. 2 Maps of A2/A0 Fourier mode strength (top panel), difference between the churning fluxes outwards and inwards (central panel, counted positive outwards), and normalised metallicity evolution (bottom panel) with respect to time and galactocentric radius, R, or guiding radius, Rg(t), for the K9i (interacting Q=1.5 simulation, top three panels) and K7 simulations (isolated Q=1.5 simulation, bottom three panels). The grey area masks the first 0.6 Gyr displayed, when the disc is not fully relaxed yet. For the K9i maps, the times of the pericentre passages of GSMC are shown in dashed white lines, and its impact parameters are the white stars. The corotation radius is plotted in green and yellow for resonance with a bar or a weak central pattern, respectively.

In the text

	Fig. A.1 Same as Fig. 2 but for K6i (top three panels, Q=1 interacting simulation) and K4 (bottom three panels, Q=1 isolated simulation).
In the text

	Fig. A.2 Same as Figs. 2 and A.1 but for K3i (top three panels, Q=1.2 interacting simulation) and K1 (bottom three panels, Q=1.2 isolated simulation).
In the text