Role of magnetic pressure forces in fluctuation dynamo saturation

Using magnetohydrodynamic simulations of fluctuation dynamos in turbulent flows with rms Mach numbers $\mathcal{M}_{\rm rms} = 0.2, 1.1$ and $3$, we show that magnetic pressure forces play a crucial role in dynamo saturation in supersonic flows. First, as expected when pressure forces oppose compression, an increase in anticorrelation between density and magnetic field strengths obtains even in subsonic flows with the anti-correlation arising from the intense but rarer magnetic structures. In supersonic flows, due to stronger compressive motions density and magnetic field strength continue to maintain a positive correlation. However, the degree of positive correlation decreases as the dynamo saturates. Secondly, we find that the unit vectors of $\nabla\rho$ and $\nabla B^{2}$ are preferentially antiparallel to each other in subsonic flows. This is indicative of magnetic pressure opposing compression. This antiparallel alignment persists in transonic and supersonic flows at dynamo saturation. However, compressive motions also lead to the emergence of a parallel alignment in these flows. Finally, we consider the work done against the components of the Lorentz force and the different sources of magnetic energy growth and dissipation. We show that while in subsonic flows, suppression of field line stretching is dominant in saturating the dynamo, the picture is different in supersonic flows. Both field line stretching and compression initially amplifies the field. However, growing magnetic pressure opposes further compression of magnetic flux which tends to reduce the compressive motions. Simultaneously, field line stretching also reduces. But, suppression of compressive amplification dominates the saturation of the dynamo.

Given the ubiquity of the occurrence of FDs in astrophysical objects, we revisit the mechanism of FD saturation in turbulent flows in this work.The standard picture of the non-linear evolution of the dynamo is the following.For magnetic Prandtl numbers Pm > 1 and fluid Reynolds number (Re) large enough such that the velocity spectrum consists of multiple scales, the dynamo is excited initially by random motions on small enough scales ℓ << l f , where Rm(ℓ) > Rmcr (Haugen, Brandenburg & Dobler 2004;Schekochihin et al. 2004;Shukurov & Subramanian 2021), where l f is the driving scale of turbulence.The amplification of the field due to such eddies continue till the local magnetic energy density approaches values comparable to the energy density of fluid motions on similar scales.Further growth of the field by such eddies is suppressed due to non-linear back reaction arising from the Lorentz force.However, turbulent eddies at a neighbouring larger scale can still continue to grow the field as they have a larger energy.Gradually over time, scale-by-scale saturation of the field occurs with the peak of the magnetic energy spectrum approaching a fraction of the outer scale (or the forcing scale of turbulent motions) (Shukurov & Subramanian 2021, and references therein).
Although the above mechanism has been confirmed in numerical simulations by a number of authors mentioned earlier, the exact manner in which the approach to a global non-linear steady state occurs still remains unclear.For this purpose it is instructive to consider a qualitative picture of FD amplification.In an ideal plasma, B/ρ ∝ δr, where B is the magnetic field, ρ the plasma density and δr the infinitesimal separation vector of fluid elements along the field line (Brandenburg & Subramanian 2005).If δr increases in a mean square sense, due to random stretching of the field, and if ρ is assumed strictly constant as would obtain in an incompressible flow, the magnetic energy also increases.Suppression of random stretching of the magnetic field lines by the growing magnetic tension (of the Lorentz force) has been recognised as the key mechanism ushering in the saturation of the dynamo (Schekochihin et al. 2002).
However, there is another possibility to limit growth in case we do not impose strict incompressibility, which obtains due to the effect of magnetic pressure.When a region containing a flux of field lines (or a flux tube) of length 'L' is stretched, magnetic pressure opposes the compression of the flux tube cross-section of area 'A'.Consequently, ρ decreases to conserve mass, and the magnetic field strength B then need not increase.This mechanism could be important for FD saturation in the interstellar medium (ISM) of galaxies, where turbulence is expected to be supersonic at the driving scale.Analysis of the alignment between the magnetic field and the eigenvectors of the rate of strain tensor (Sij) in compressible magnetohydrodynamic (MHD) simulations by Sur, Pan & Scannapieco (2014b) further showed that compressive motions are statistically oriented perpendicular to the direction of stretching.In a similar vein, magnetic pressure can also oppose any purely compressive motions and prevent the increase of ρ hence B. We examine here how important is this mode of saturation for the FD in random compressible flows, even if the compressibility effects are small as in subsonic flows.Note, however, that the role of magnetic pressure in saturating the dynamo cannot be captured in numerical simulations which solves the exactly incompressible MHD equations.
The paper is organised as follows.In Section 2, we present in brief the initial conditions and the set-up of the simulations.Thereafter, results of our investigation is presented under various sections.To guide the readers, we begin by first discussing the evolution of the rms Mach number of the flows, the ratio of magnetic to kinetic energies, and the PDFs of the magnetic energy density in Section 3.With the overarching goal to elucidate the role of magnetic pressure in FD saturation, we then proceed to explore the time evolution of the ρ − B correlation and how regions with varying field strengths (as obtained in the FD) and over densities contribute to this correlation in Section 4. If magnetic pressure plays a role as described earlier, we expect this correlation to decrease as the field becomes dynamically important.Next, in Section 5, we study in detail the nature of the alignment between the magnetic field direction, the gradient of the density, and the components of the Lorentz force, namely the gradient of the magnetic pressure and the magnetic tension and explore how these alignments manifest themselves in the kinematic and non-linear stages of dynamo evolution.
In Section 6, we study the work done against the magnetic tension and gradient of the magnetic pressure and how it underlines the importance of magnetic pressure in dynamo saturation.Next, similar to previous studies (e.g., Seta et al. 2020;Seta & Federrath 2021a), we also explore the relative importance of stretching, advection, compression and dissipation on the growth and saturation of magnetic energy in Section 7. Finally, in Section 8 we conclude with highlights of the key results and emphasise the role of magnetic pressure forces in dynamo saturation.

NUMERICAL METHOD AND INITIAL CONDITIONS
We focus our analysis on three non-ideal MHD simulations of non-helically forced FD simulations in the simplest possible numerical setup.The simulations were performed with an isothermal equation of state and periodic boundaries using the FLASH code1 (Fryxell et al. 2000;Benzi et al. 2008;Eswaran & Pope 1988) (version 4.2), which solves the threedimensional compressible MHD equations using explicit viscosity and resistivity in a conservative form.However, to better elucidate the role of the different terms, we express the MHD equations in a non-conservative form here : Here, ρ, U , p ′ = p + |B| 2 /2, B, denote density, velocity, total pressure (thermal and magnetic) and the magnetic field, respectively.Sij = (1/2)(∂iuj + ∂jui) − (1/3)δij∂ k u k is the traceless part of the rate of strain tensor, ρF is the artificial driving term and ν and η are the constant viscosity and resistivity, respectively.In all the simulations, turbulence is driven solenoidally (i.e., ∇ • F = 0) over a range of wave numbers between 1 ⩽ |k|L/2π ⩽ 3 (such that the average forcing wave number k f L/2π = 2 and 'L' is the length of the box) as a stochastic Orstein-Ulhenbeck (OU) process with a finite time correlation.To cover different regimes of compressibility, we analyse the data from a subsonic, a transonic and one supersonic simulation with steady state rms Mach numbers Mrms = urms/cs = 0.2, 1.1 and 3.0, respectively.Here urms is the turbulent rms velocity and cs is the isothermal sound speed.All simulations have a resolution of 512 3 , Pm = Rm/Re = 1 and are in dimensionless units.Both Re and Rm are defined with respect to the driving scale of turbulence.The initial conditions consist of a box of unit length, the density ρ and cs set to unity with zero initial velocities.The magnetic field is initialised as Binit = B0[0, 0, sin(10πx)], with B0 adjusted to a value such that the initial plasma beta β0 = p th /pmag ≈ 10 6 −10 7 .Thus, the magnetic field is dynamically weak to start with in all the simulations.Further, to maintain ∇ • B to machine precision level, we use the unsplit staggered mesh algorithm in FLASHv4.2 with a constrained transport scheme (Lee & Deane 2009;Lee 2013) and a Harten-Lax-van Leer discontinuities (HLLD) Riemann solver (Miyoshi & Kusano 2005).Table 1 lists the important parameters of these runs.

MAGNETIC ENERGY EVOLUTION AND ITS PROBABILITY DISTRIBUTION FUNCTION
Before delving into the details on the role of magnetic pressure forces in dynamo saturation, it is worthwhile to explore some of the important evolutionary features that we find in our simulations.Specifically, we show in Fig. 1 1.
The bottom panel shows that the evolution of Em/E k has three distinct phases : an initial exponential phase (kinematic) followed by an intermediate stage of slower growth and eventual saturation.It is clearly evident that the kinematic phase lasts for about t/t ed = 8, 13 and 22 in runs A, B and C, respectively.The slope of the curves in this phase are denoted by thin solid lines.Thereafter, an intermediate phase of slower growth ensues.The onset of saturation in runs A, B and C occurs from t/t ed = 15, 22 and 30 respectively.This can also be inferred from the evolution of the magnetic spectra M (k) in Fig. A1 where the spectra tend to bunch together from the aforesaid times.The bottom panel further shows that steady state value of Em/E k decreases as the compressibility of the flow increases.In the interval t/t ed = 15 − 25, the steady state value of Em/E k = 0.23 for run A, while for run B and C the steady state values are 0.22 (for t/t ed = 22 − 31) and 0.07 (for t/t ed = 30 − 47), respectively.
A distinct feature of FD is that it generates both less intense volume filling fields together with much rarer but strong fields.To illustrate our point, we show in Fig. 2, the probability distribution functions (PDFs) of (B/Brms) 3 P (B/Brms) in the kinematic (blue, dotted lines), intermediate (red, dash-dotted) and saturated phases (black, solid lines) for run A (subsonic) in the top panel and run C (supersonic) in the bottom panel.Note that (B/Brms) 2 represents the normalized magnetic energy and multiplying it by an additional factor of B/Brms to obtain (B/Brms) 3 measures the magnetic energy density contributed by fields in a logarithmic interval of B/Brms.Here P (B/Brms) is the PDF of the magnetic field strength.Two features are evident from these plots.First, in the kinematic phase, magnetic fields of a wide range of strengths are present and contribute to the total energy.Especially in the supersonic case, we find fields with strengths ≈ 10Brms contribute comparable energy density to the total as the rms fields.This is not the case with the subsonic run.This feature could be due to the continuous amplification of fields by density compression in addition to amplification by stretching.It is further evident that in both intermediate and saturated phases, the peak of the PDFs in both subsonic and supersonic cases lies at ≈ 2 B/Brms.However, the prob- ability of fields B > 6Brms contributing to the total energy decreases gradually from the kinematic to intermediate and saturated phases in the supersonic case.Strong field regions with B/Brms > 2 also contribute to the magnetic energy with a probability of 10 −2 for B/Brms = 6 in the supersonic case.This implies that together with the general sea of volume filling fields, strong field regions may also play a role in dynamo saturation, particularly when compressibility effects are significant.We also note that in the supersonic case, regions with much higher field strengths (B/Brms > 6) continue to exist albeit with a much lower probability.This could imply that even though the dynamo suppresses the compressive motions when the fields become dynamically important, density compression still manages to amplify fields in the rarer regions.

CORRELATION BETWEEN DENSITY AND MAGNETIC FIELD STRENGTH
In our simulations, density and magnetic field strength are uncorrelated to start with.To explore the temporal evolution of the (anti)correlation between them as the FD evolves, we compute the Pearson correlation coefficient (rp) defined as, z is the magnitude of the field, ρ i,j,k and B i,j,k are the density and the magnetic field strength at a point (i, j, k) in the simulation volume.ρ and B are the mean values of density and B, respectively.
Fig. 3 shows that in the subsonic case (run A, blue dotted line with asterisks) ρ and B evolves to become anticor- related after about one eddy-turnover time and eventually settles at an average value ⟨rp⟩ = −0.3 in the saturated state of dynamo.Thus, when the effects of compressibility are weak, high-density regions correspond to low magnetic field strength regions and vice versa.This is a consequence of the magnetic pressure acting to oppose the compression of the flux tube.
When compressibility effects become important (with the increase in Mrms) the growth of the field is expected to be driven by a combination of amplification by random stretching and compression in converging flows which also enhances the density.The former amplifies fields in regions where the vorticity is strong and the latter in regions where over-densities are large.This is because, the solenoidal nature of turbulent driving leads to the production of vortical motions that drive field amplification by random stretching (Haugen, Brandenburg & Dobler 2004;Federrath et al. 2011a;Porter, Jones & Ryu 2015;Achikanath Chirakkara et al. 2021) in contrast to compressive motions.However, regions of strong vorticity may not correspond to regions with large over-densities.In the kinematic phase, when fields are dynamically weak, amplification due to density compressions seem to dominate the evolution of the correlation coefficient.This leads to an initial positive correlation with the degree of correlation strongly dependent on the rms M number of the flow as shown also by Yoon, Cho & Kim (2016) and Seta & Federrath (2021b).
Indeed, in the transonic case (run B, black line with squares), we observe a positive correlation with ⟨rp⟩ = 0.33 which then gradually declines as the dynamo saturates.After about 25t ed , the density and the magnetic field strength becomes uncorrelated with ⟨rp⟩ = −0.03 in the saturated phase.Thus, with exception to the initial positive correlation, the evolution of rp(ρ, B) in the transonic case bears resemblance to the corresponding evolution in the subsonic case.The decrease in rp can be understood as a consequence of magnetic pressure forces gaining in importance, which acts to further resist the compressive motions resulting in suppression of density enhancements.This results in regions of strong fields amplified by random stretching, that are not necessarily associated with high density regions which arose from compression.On the other hand, in the supersonic case (run C, red dotted lines with triangles) where δρ/ρ ∝ M 2 ∼ 9, the dominance of density compressions amplifying the field pushes the ⟨rp⟩ = 0.6 up to ≈ 20 t ed , which corresponds to the duration of the kinematic phase.Thereafter, the correlation decreases only slightly as the field evolves to the non-linear saturated state with ⟨rp⟩ = 0.43 in the interval t/t ed = (30 − 47).This suggests that when compressible effects are strong, field amplification due to density compression dominates the evolution of rp even when Lorentz forces are strong enough to be dynamically significant.

Evolution of the correlation coefficient for different ranges in B/Brms and over densities
A generic feature of FD is that it generates magnetic fields consisting of rarer, intense field structures embedded in a sea of less intense, volume filling fields.Therefore, it is worthwhile to distinguish the contribution to rp(ρ, B) seen in Fig. 3 from regions with differing magnetic field strengths.Note here that for a given range of B/Brms we focus only on those regions in the simulation volume where the range condition is satisfied.Thus, the data in each of the ranges discussed below correspond to a subset of the full data set and equation ( 2) is solved for the respective data sets.In Fig. 4, we show the evolution of rp for the subsonic (top panel), transonic (middle) and supersonic (bottom panel) for two different ranges of B/Brms while the black solid lines in each panel depicts the variations of rp as in Fig. 3, where no such ranges in B/Brms are considered.For the subsonic case, it is now clearly evident that the strong anticorrelation ⟨rp⟩ = −0.46 in the range t/t ed = 13 − 25 arises from the rarer, intense field regions where B/Brms > 3 (green dashed line with asterisk).In the absence of strong density compressions, such intense magnetic field regions solely arise due to random stretching.Thus, the observed strong anticorrelation implies that magnetic pressure exerted by these strong fields resist further compression of flux tubes.At the same time, the sea of less intense volume filling fields, however, remain uncorrelated with the density at all times.
In the transonic case, initial density enhancements due to compression also enhances the magnetic field.This leads to an initial positive correlation, although the degree of positive correlation is lesser when regions with fields B/Brms > 3 are considered.However, the plots clearly show that once the dynamo builds up the fields, the degree of positive correlation steadily decreases even in regions where B/Brms ⩽ 1, while for the strong field regions, rp(ρ, B) evolves to be anticorrelated.Thus, with exception to the initial positive correlation, the general trend of rp(ρ, B) is very similar to the one observed in the subsonic case.Finally, in the supersonic case ρ and B remain positively correlated throughout the evolution.However, regions with B/Brms > 3 only show a weaker degree of positive correlation compared to regions with B/Brms ⩽ 1.This again reflects the fact that such strong and intense field regions arise due to amplification by field line stretching that obtains in vortical turbulence.However, density enhancements due to strong compressions continue to dominate over the magnetic field enhancements due to random stretching leading to a net albeit weaker positive correlation.
The fact that we consider three simulations with varying degrees of flow compressibility also makes it equally important to probe the degree of (anti)correlation in differ-ent over density ranges.Recall that density inhomegenities δρ/ρ ∝ M 2 , where δρ = (ρ − ρ) is the fluctuation in density.While δρ/ρ is negligible in the subsonic case, significant density fluctuations ∝ M 2 could influence the correlation in transonic and supersonic regimes.To explore this, we compute the evolution of rp(ρ, B) in density ranges segregated by regions : δρ/ρ < M 2 and δρ/ρ ⩾ M 2 .The results are shown in Fig. 5 with panel (a) showing the evolution for the subsonic, panel (b) for transonic, and panel (c) for the supersonic run.
We find that in the subsonic case, the observed anticorrelation arises from regions with δρ/ρ < 0.04 as the magnetic fields in these regions are amplified solely due to vortical motions.Regions with δρ/ρ ⩾ 0.04 remain uncorrelated.In the transonic and supersonic cases, the evolution of rp in regions with δρ/ρ < M 2 follow the pattern seen when no cutoffs are introduced (black solid lines).On the other hand, in both cases, rp starts out with a much weaker positive correlation when over dense regions with δρ/ρ ⩾ M 2 are considered.In the transonic case, magnetic field strength eventually becomes uncorrelated with the density while in the supersonic case, they maintain a weaker positive correlation.This could be due to the fact that such high-density regions are fewer in number in the simulation domain.Furthermore, we find that the values of rp calculated using a mass-weighted version of equation ( 2) are almost identical to the ones shown in Fig. 5 particularly for over densities δρ/ρ > M 2 .

2D slices of density fluctuations and fluctuations in magnetic energy
Further evidence of the anticorrelation between the density and magnetic field strength obtained in run A can be easily discerned from the two-dimensional (2D) slices of the density and magnetic energy.Such slices can also reveal the positive correlation expected in supersonic flows.
To illustrate this, we show in Fig. 6, 2D slices of δρ/ρ (left column) and the corresponding fluctuations in B 2 , δB 2 /B 2 rms = (B 2 − B 2 rms )/B 2 rms (right column).These slices were obtained in the saturated phase of the dynamo in the respective runs.The top row shows the slices in the x − y plane (at z = 0) corresponding to t/t ed = 24.5 for run A, while the bottom row shows the same at t/t ed = 48 from run C.
It appears that fluctuations in both ρ and B 2 resemble that of a cloth which is wrinkled and folded at multiple regions due to action of the turbulent driving.Note that even though both ρ and B 2 are positive definite quantities, there are regions within the simulation volume where ρ < ρ (ρ = 1) and B 2 < B 2 rms .Both the left and right panels of the figure shows a wide variety of high contrast but smoothly varying structures in δρ/ρ and δB 2 /B 2 rms .More importantly, in run A we find that regions of strong density fluctuations correspond to regions of weak fluctuations in magnetic energy and vice versa.The strongest anti-correlation occurs in the rarer, strong field regions with δB 2 /B 2 rms > 7.8.This once again corroborates the results obtained in Section 4.1.On the other hand, in comparison to run A the bottom row shows that fluctuations in both δρ/ρ and δB 2 /B 2 rms are much stronger in run C.It is also clearly evident that density fluctuations are positively correlated with fluctuations in magnetic energy in agreement with the results presented earlier.

ALIGNMENT ANGLES
The evolution of rp(ρ, B) offers a first glimpse into the possible role of magnetic pressure in saturating the dynamo.Additionally, it is also worthwhile to investigate the nature of the alignments between the density gradient, magnetic field, and the components of the Lorentz force.Accordingly, we show in Fig. 7, the PDFs of the cosine of the angles between the magnetic field and the gradient of density (left column), the gradient of density and gradient of the magnetic pressure (middle column), and finally, between the gradient of the density and the magnetic tension (right column).These PDFs are computed by averaging over a number of independent realisations (at different t ed ) in kinematic (blue, dash-dotted), intermediate (yellow, dashed) and saturated phases (red, solid).
The left column of the above figure shows that the distributions of the unit vectors of ∇ρ and B are symmetric about cos(θ) = 0 with the unit vector of the magnetic field (directed along the field line) preferentially perpendicular (90 • or 270 • ) to the unit vector of the density gradient.This is consistent with the simple picture of a flux tube threaded by a field of strength 'B' where the density can be enhanced or depleted in the flux tube.Moreover, the degree of orthogonal alignment appears to be stronger in the saturated and intermediate phases in comparison to the kinematic phase.This reinforces in a statistical sense, the simple qualitative picture described earlier.Due to magnetic pressure gradience opposing compression of the field, density is decreased in strong field regions and density contrasts are enhanced perpendicular to the field.On the other hand, the gradient of the magnetic pressure (∇B 2 ) is directed inwards to the flux tube.The middle column in the figure shows that in the subsonic case (run A), ∇ρ and ∇B 2 are antiparallel (i.e. at 180 • ) to each other in all the three phases.This is a manifestation of the fact that as magnetic pressure gains in importance due to dynamo amplification, the force due to magnetic pressure opposes the compression of the flux tubes.This results in density variations that are anticorrelated with variations in magnetic pressure.Thus, force due to magnetic pressure tends to empty out the flux tubes rather than allowing the cross-sectional area to decrease.Note that in the subsonic case, this antiparallel alignment is caused by the rarer, stronger field regions that are strongly anticorrelated with the density, while the general sea of volume filling fields remain uncorrelated (see Fig. 4).The kinematic phase in the transonic case (run B) shows an antiparallel alignment together with a weak parallel alignment between the two, which evolves to a stronger anti-parallel alignment in the intermediate phase and continuing to the saturated phase.
However, the picture changes considerably in the supersonic case (run C) where compressibility effects are strong enough that in the kinematic phase, density enhancements drive the amplification of the magnetic field.In this case, the strong positive correlation seen in Fig. 3 results in density variations that are also correlated with magnetic pressure variations.Consequently, ∇ρ and ∇B 2 are initially aligned parallel to each other.However, over time random stretch- ing also amplifies magnetic fields in addition to field amplification due to compressions.Similar to the subsonic case, we find that as the dynamo saturates force due to magnetic pressure opposes further compression of the field lines resulting in the emergence of an antiparallel alignment between the two that coexists with the parallel alignment.
The right column shows that magnetic tension does not show any preferred alignment with the density gradient in the subsonic case.However, with the increase in flow compressibility, a weak orthogonal alignment between the two emerges which becomes more prominent at Mrms = 3.Thus, the nature of the alignments discussed above provides further appreciation of the role of magnetic pressure in dynamo saturation.

IMPACT OF LORENTZ FORCES ON MAGNETIC ENERGY EVOLUTION
In this section, we discuss how the energy exchange between the velocity and magnetic fields facilitated by the different components of the Lorentz force elucidates the role of magnetic pressure and magnetic tension in saturating the dynamo.
The volume averaged magnetic energy evolution in terms of the work done by the Lorentz force is where σ is the conductivity and the angular brackets ⟨...⟩ denotes volume averaging.Substituting J ×B = (B •∇)B − ∇(B2 /2) in the above equation and multiplying both sides by t ed /B 2 rms , the magnetic energy evolution in dimensionless form is, The first two terms on the right-hand side correspond to the rate of work done by the magnetic tension LT = LT − ( B • LT) B and the gradient of the magnetic pressure Lp = Lp − ( B • Lp) B. Here, LT = (B • ∇)B and Lp = ∇(B 2 /2) and B is the unit vector of the magnetic field.Thus, LT and Lp as defined above only retain components perpendicular to B. This is because the Lorentz force has no component parallel to B. The last term represents the decrease in magnetic energy due to Joule dissipation. 2The effect of this term can be gauged from Figs 9 and 10 presented in the following section.
Using the instantaneous values of urms and Brms for the normalisation parameter t ed /Brms, we show in Fig. 8 the time evolution of the average values of the first term (blue dotted line with diamonds), and the second term (red, dashed line with '+' symbol).The fluctuations in the velocity are within 10 per cent of the rms value.In all the three cases corresponding to different degrees of flow compressibility, both the terms have positive average values which implies that they contribute to the growth in magnetic energy (see equation 4).In the subsonic and transonic cases (top and middle panel), the work done against the gradient of the magnetic pressure is always small compared to the magnetic tension.The plots also show that both magnetic pressure and magnetic tension are involved in dynamo saturation.We see that the average values of both the terms decrease as the dynamo evolved from the kinematic to the saturated phase.In the subsonic case, both magnetic pressure and tension contributions decrease by identical factors.However, in the transonic case, the suppression of the magnetic pressure term is more compared to the decrease due to the magnetic tension.Interestingly in the supersonic case (bottom panel), the energy gained from the velocity due to ⟨U • LT⟩ and ⟨U • Lp⟩ terms are comparable up to 22t ed .As can be seen from Fig. 1, this corresponds to the time by when the kinematic growth phase of the magnetic energy culminates and the intermediate phase of growth begins.Similar to the transonic case, once the dynamo evolves to the saturated phase, the work done against the gradient of the magnetic pressure is suppressed more in comparison to the magnetic tension.This again reflects the fact that the magnetic pressure gains in importance in opposing compressive motions as the dynamo saturates.
Thus, expressing the equation for magnetic energy evo-lution in terms of the rate of work done against the components of the Lorentz force offers a fresh perspective which highlights the significance of magnetic pressure in the supersonic case.

IMPACT OF LOCAL STRETCHING, ADVECTION, COMPRESSION, AND DISSIPATION ON MAGNETIC ENERGY EVOLUTION
Another approach to understand how the FD saturates is to analyse the effects of local stretching, advection, compression, and dissipation on the growth or decay of magnetic energy (Seta et al. 2020;Seta & Federrath 2021a).
In this context, it is important to clarify if suppression of field line stretching is the only agent responsible for saturating the dynamo.To this end, we take the dot product of the magnetic induction equation with B and then expand the B • [∇ × (U × B)] term to obtain the evolution of the mean magnetic energy.In the following, we present a detailed derivation of the magnetic energy evolution equation and also highlight a subtle point.Assuming η to be constant, the i-th component of the magnetic induction equation is where J = ∇ × B is the current density with µ0 = 1.
Expanding the cross-product in the induction term using Levi-Civita symbols yield Now, substituting ϵ ijk ϵ klm = (δ il δjm − δimδ jl ) and utilising the properties of Kronecker Delta δij's gives, Taking a dot product of equation ( 7) with Bi, and neglecting the term ∝ ∇•B we obtain the evolution equation of the magnetic energy in terms of local stretching, advection, compression, and dissipation terms, Next, we define the rate of the strain tensor Sij = (Ui,j + Uj,i)/2 where Ui,j denote partial derivative of Ui w.r.t.j and vice versa.Now, expressing the stretching term in equation ( 8) in terms of Sij and recasting the third term Bi(Uj∂j)Bi as Uj∂j(B 2 i /2) the volume integrated magnetic energy evolution is given by Note that in deriving equation ( 9), the term ∝ (Ui,j − Uj,i) is anti-symmetric in (i, j) and hence does not contribute to the magnetic energy.
The first term on the right-hand side of equation ( 9) represents the stretching and compression of the magnetic field lines by the turbulent flow which may increase or decrease the magnetic energy.The second term reflects the effects of advection of the field lines by the turbulent flow.The third term representing the effects of compression can also reduce the magnetic energy depending on the divergence of the velocity field and is generally important in transonic and supersonic flows.The last term is dissipative in nature which acts to reduce the magnetic energy.The effects of local stretching, advection, compression, and dissipation on the magnetic energy growth (or decay) can then be analysed from the PDFs of each of the terms on the right-hand side of equation ( 9).
Before proceeding further, we wish to elucidate a subtle point.We note that the advection term in equation ( 9) can be further simplified to which together with the compression term leads to a term ∝ (∇ • U ) and a surface term ∝ (U |B 2 |).While for periodic boundary conditions, the latter term vanishes, it will still contribute when computing the PDFs at each mesh point.Similar reasoning is applicable to the dissipation term in equation ( 9) which can be expressed as a combination of a divergence term ∝ ∇ • (J × B) and a term ∝ |J | 2 .Again, the divergence term vanishes for periodic boundaries but will contribute to the PDFs at each mesh point.Therefore, while computing the PDFs it is essential to retain the contribution from the advection and dissipation terms in the form given in equation ( 9).This is in contrast to Seta & Federrath (2021a), where the aforesaid surface terms were neglected.
Similar to in Section 6 we express equation ( 9) in dimensionless form by multiplying both sides of the equation by t ed /B 2 rms which results in While computing the PDFs of the compression, advection, and dissipation terms are fairly straightforward, to compute the PDFs of the local stretching we first project the local magnetic field along the three eigenvectors of the rate of strain tensor (e1, e2, e3).The corresponding eigenvalues are denoted as λ1 > λ2 > λ3, respectively (Brandenburg 1995; Sur, Pan & Scannapieco 2014b;Seta et al. 2020;Seta & Federrath 2021a).In the incompressible case, λ1 + λ2 + λ3 = 0 while in the supersonic case the sum of the three eigenvalues need not be zero.In general, e1 and e3 correspond to the directions of local stretching and compression while e2 can be either, depending on sign of λ2 (Ashurst et al. 1987).
We then compute the local stretching as, Σi=1,2,3[λit ed (B • ei) 2 ]/B 2 rms , where (B • ei) denotes the component of the local magnetic field along the direction of the eigenvectors, and λi are the eigenvalues corresponding to the eigenvectors ei.
In the kinematic phase, there is a net growth in the magnetic energy, while in the saturated phase the magnetic energy is expected to remain constant.It is expected that this is caused by the mutual adjustments between the local stretching, compression, advection, and dissipation terms.In what follows, we discuss the effects of each of these terms in two regimes of compressibility focussing on run A where Mrms = 0.2 and run C where Mrms = 3.0.Since we normalised equation (10) using the instantaneous values of urms and Brms, we refrain from computing the total PDFs over the range of t ed in the kinematic and saturated phases.Instead, we show in Fig. 9 the time evolution of the mean values of the PDFs of local stretching (µs), advection (µa), compression (µc), and dissipation (µ d ) and the sum total of the mean values (µtot) for the subsonic run with Mrms = 0.2.In the same vein, Fig. 10 shows the evolution for the supersonic run with Mrms = 3.For the sake of completeness, the PDFs of stretching, advection, compression, and dissipation terms at a few representative times in the kinematic and saturated phases in shown in Figs B1 and B2.
In the subsonic case (run A), Fig. 9 shows that field amplification in the kinematic phase is predominantly driven by an increase in local stretching of the field by the turbulent flow.This is accompanied by an increase in dissipation.The mean values of local stretching, advection, and dissipation peaks ≈ 3 − 4t ed beyond which, the growth rate of the magnetic energy changes (see Fig. 1).In fact, Fig. B1 show that dissipation is relatively stronger than both field line stretching and advection initially.Beyond this, the change in growth rate of magnetic energy is manifested by a gradual decrease in µs, µa and µ d .This continues up to t/t ed ≈ 15 beyond which, all three attain steady state values.However, µc remains negligible throughout the evolution.This reconfirms the fact that compression plays a negligible role in dynamo saturation in subsonic flows.The net growth of magnetic energy in the kinematic phase is confirmed from the non-zero value of the total mean.Of course, as the dy- namo approaches the non-linear phase, the net growth consistently decreases.In the saturated phase, field line stretching and dissipation compensate each other (i.e., |µs| ≈ |µ d |) while µa is totally suppressed to negligible values.A quick comparison between the PDFs in Fig. B1 show that both local stretching and dissipation are reduced in the saturated phase.Note that the rate of strain tensor in the stretching term in equation ( 10) is only affected by magnetic tension.Thus, the reduction in stretching observed in the saturated phase suggest that magnetic tension is resisting the further stretching of the field lines.
In contrast to the subsonic case, Fig. 10 shows that amplification of the magnetic field in the kinematic phase in supersonic turbulence (run C) is dominated by compression of magnetic field lines.Growth of magnetic energy due to stretching of the field lines remains subdominant to the effects of compression by almost a factor of two in this phase.The compressible nature of turbulence (Mrms = 3.0) drives local density enhancements which further amplifies the magnetic field.On the other hand, the effects of advection remain marginal compared to compression, stretching, and dissipation.As the magnetic energy transitions to the intermediate phase of growth, both compression and field line stretching reduces marked by a drop in µs and µc after t/t ed ≈ 22.At the same time, advection and dissipation also decreases slightly.In the saturated phase, µs and µc become comparable to each other.The overall evolution of the mean values spanning from kinematic to saturated phase show that field amplification by compressive motions suffer the strongest suppression as µc decreases by a factor ≈ 3.3 (also see Fig. B2) compared to a factor of ≈ 1.5 for the reduction in stretching.This implies that in the compressible regime, in addition to reduction in field line stretching, advection, and dissipation, the suppression of further compression of field lines due to the action of magnetic pressure plays a major role in saturating the dynamo.Interestingly, the aforesaid deduction is consistent with Sur, Pan & Scannapieco (2014b) where it was shown that compressive motions of the strain (perpendicular to the direction of stretching) are suppressed in the saturated phase of the FD.

CONCLUSIONS
Over the course of the last two decades, a great deal of progress has been made to understand the manifold features of FDs and it's application to magnetic fields in the ISM of galaxies and in the intracluster medium of galaxy clusters.A majority of these studies were based on direct numerical simulations that dealt with amplification of dynamically insignificant seed magnetic fields and subsequent saturation in flows encompassing various degrees of compressibility.However, the precise way in which the dynamo evolves to the non-linear saturated state and the factors that play a deciding role in achieving this steady state remained far from understood, more so in compressible flows.While in a broader sense, the back reaction due to the Lorentz force is rightly adjudged to be responsible for saturating the dynamo, it had been argued that the magnetic tension force is the main ingredient that suppresses Lagrangian chaos via reduction of field line stretching which in turns results in the saturation of the dynamo (Cattaneo et al. 1996;Rincon 2019) (also see Eyink (2011) and Eyink et al. (2011)).
Although the role of magnetic tension is appreciated for dynamo saturation in incompressible flows, a better understanding is required in compressible flows which are likely to occur in the galactic ISM where density fluctuations can be significant.In particular, two important questions arisedoes magnetic pressure gradient associated with the Lorentz force also play a role in saturating the dynamo and, is suppression of field line stretching the only mechanism responsible for dynamo saturation when the condition of strict incompressibility is relaxed?
With the aforesaid questions in mind, here we have provided a new perspective on how FDs saturate.To this effect, we performed numerical simulations of FDs at Pm = 1 in turbulent flows driven at half the scale of the box with the amplitude of turbulent driving adjusted such that the steady state Mrms = 0.2, 1.1 and 3.0.We then analysed the data from these simulations using a variety of probes.These include : (i) exploring the correlation between the density and magnetic field strength, (ii) the nature of the alignments between the magnetic field, the gradient of the density, and the components of the Lorentz force, (iii) the work done against different parts of the Lorentz force to grow the magnetic field, and (iv) the impact of local stretching, advection, compression, and dissipation on the evolution of the magnetic energy.Our main conclusions are as follows : (i) The PDFs of the magnetic energy density shown in Fig. 2, demonstrate in a novel manner that together with the general sea of less intense, volume filling fields, strong fields regions are also likely to play a role in dynamo saturation in supersonic turbulence.In fact, as Fig. 2 shows strong fields of strength ≈ 10Brms in the supersonic case contribute comparable energy density as the rms strength fields in the kinematic phase.Even after saturation, fields upto 4Brms contribute significantly to the magnetic energy density.(ii) Density and magnetic fields strength are anticorrelated in subsonic flows.Due to the incompressible nature of the flow, densities cannot be changed by a large factor.Thus, saturation of the dynamo can happen by suppression of field line stretching.This is also confirmed from Fig. 9. (iii) As the effects of compressibility become important (e.g.runs B and C) magnetic fields are amplified by a combi-nation of random stretching and compression.The latter ensures an initial positive value of rp which then decreases as dynamo saturates (see Fig. 3).In transonic flows, ρ and B tend to become uncorrelated in the saturated state while in the supersonic case, rp shows a lesser degree of positive correlation in the saturated phase compared to the value in the kinematic phase.(iv) A quick glance at Fig. 4 reveals that the anti-correlation between ρ and B seen in the subsonic case stems from the strong field regions with B/Brms > 3.Such rarer and intense fields also lead to negligible rp in the transonic case and a much weaker positive correlation in the supersonic case.Interestingly, we find that regions with δρ/ρ < M 2 (see Fig. 5) lead to the anticorrelation in the subsonic case.
In the transonic case, these over dense regions tend to be uncorrelated with the magnetic field while in the supersonic case, they tend to reduce the degree of positive correlation.(v) The fact that magnetic pressure forces play a significant role in dynamo saturation is clearly evident from the PDFs of the alignment presented in Fig. 7.In the subsonic case, we find that n∇ρ is statistically aligned anti-parallel to n ∇B 2 , with the alignment becoming stronger in the saturated phase.In the transonic and more prominently in the supersonic case, apart from the antiparallel alignment, a strong parallel alignment between the two emerges in the kinematic phase.Nevertheless, the antiparallel alignment due to magnetic pressure forces dominate in the saturated phase.In contrast, magnetic tension does not show any preferred alignment with the gradient of the density irrespective of the flow compressibility.(vi) An analysis of the impact of Lorentz forces on the growth of magnetic energy (Fig. 8) reveals that the contribution from the magnetic pressure gradient is always small compared to the magnetic tension up to Mrms = 1.1.However, in the supersonic regime, the two are comparable till the time the non-linear saturated phase ushers in.Thereafter, the work done against the gradient of magnetic pressure (aiding field growth) is suppressed more in comparison to the magnetic tension.This again illustrates the importance of magnetic pressure gradients in dynamo saturation for supersonic turbulence.(vii) The time evolution of the mean values of the PDFs of local stretching, advection, compression, and dissipation terms are shown in Figs 9 and 10.They clearly demonstrates firstly that local stretching of the field lines play a crucial role in amplifying the field in subsonic flows.And the saturation of the dynamo in such flows occurs via a reduction in field line stretching.In contrast, in supersonic flows, amplification of magnetic fields is dominated by compression of the field lines in the kinematic phase with stretching playing a subdominant role.As the magnetic energy evolves to equipartition strengths, the dynamo saturates more due to suppression of compressive amplification than due reduction in field line stretching.
In a nutshell, the results from our work confirm in detail a new route to dynamo saturation based on the simple schematic picture of magnetic flux tubes discussed in Section 1.In incompressible flows, where magnetic tension plays a decisive role in dynamo saturation, magnetic pressure forces nevertheless act simultaneously to oppose the compression of the flux tubes, especially in rare but intense magnetic structures.The pivotal role of magnetic pressure forces in saturating the dynamo comes into prominence in transonic and particularly in supersonic flows.There, growing magnetic pressure forces resist compression of flux tubes and lead to regions of strong magnetic flux being emptied out rather than a decrease of their cross-sectional areas.This fact has not been sufficiently appreciated so far and could be important as many astrophysical systems host compressible turbulent flows.Future work extending these considerations also to cases with Pm ≫ 1, and where turbulence is inhomogeneous as expected in galaxies, would be of interest.
FDs play an important role in amplifying weak seed magnetic fields to strengths of dynamical importance.In astrophysical systems like galaxies, where conditions are favourable for large-scale dynamo action, the fields amplified by FDs provide a seed for the later evolution (Shukurov & Subramanian 2021).Here, turbulence is mainly driven by supernovae explosions.On the other hand, during the formation of the First stars FD action can be triggered due to gravitational collapse which generates turbulent motions (Sur et al. 2010(Sur et al. , 2012;;Federrath et al. 2011b).On even larger scales, the same mechanism may amplify fields during the formation of the first galaxies (Latif et al. 2013;Schober, Schleicher & Klessen 2013;Pakmor et al. 2017) and also in galaxy clusters (Subramanian, Shukurov & Haugen 2006;Cho & Ryu 2009;Vazza et al. 2018).In the latter, such dynamo generated fields crucially impact the properties of polarized synchrotron emission (Basu & Sur 2021;Sur, Basu & Subramanian 2021).Thus, elucidating the role of the magnetic pressure and magnetic tension in dynamo saturation is of lasting importance.

Figure 1 .
Figure 1.Time evolution of Mrms (top panel) and the ratio of Em/E k (bottom panel) for all the runs.The thin horizontal lines in the top panel are the average values of Mrms in the time interval considered.The thin black lines in the bottom panel denote the slopes of Em/E k in the kinematic phase .

Figure 3 .
Figure 3. Evolution of rp(ρ, B) for runs A, B and C. The blue dotted, black solid, and red dashed vertical lines indicate the time of transition to the saturation phase in the respective runs.

Figure 4 .
Figure 4. Evolution of rp(ρ, B) for runs A, B and C for different ranges of B/Brms.The dotted vertical lines denote the transition to the saturated phase in the respective runs.

Figure 5 .
Figure 5. Evolution of rp(ρ, B) for runs A, B and C for over density ranges segregated by regions : δρ/ρ < M 2 and δρ/ρ ⩾ M 2 .As before, the dotted vertical lines in all the panels indicate the transition to the saturated phase in the respective runs.

Figure 6 .
Figure 6.2D slices of δρ/ρ (left column) and δB 2 /B 2 rms (right column) from snapshots in the saturated state at t/t ed = 24.5 from run A (top row) and at t/t ed = 48 from run C (bottom row).The slices are shown in the x − y plane at z = 0.The anticorrelation between the two in subsonic flows in clearly evident with bright regions of δρ/ρ corresponding to dark regions (weaker fluctuations) in δB 2 /B 2 rms

Figure 7 .
Figure 7. PDFs of the cosine of the angles for the subsonic (Mrms = 0.2) in the top row, transonic run (Mrms = 1.1) in the middle, and for the supersonic (Mrms = 3.0) run in the bottom row.The alignments are shown between : n B • n ∇ρ (left column), n ∇ρ • n ∇B 2 (middle column), and n ∇ρ • n (B•∇)B (right column).

Figure 8 .
Figure 8. Evolution of the first and second term in equation 4 for runs A (top panel), B (middle panel), and C (bottom panel).

Figure 9 .
Figure 9. Evolution of the mean values (µ) of the PDFs of local stretching (blue, asterisks), compression (green, filled stars), advection (red, open diamonds) and dissipation term (purple, '+') for the subsonic run with Mrms = 0.2.The sum total of the mean values is depicted by black crosses.

Figure 10 .
Figure10.Same as in Fig.9but now for the supersonic run with Mrms = 3.0.

Figure A1 .
Figure A1.Time evolution of the kinetic K(k), and magnetic M (k) energy spectra in runs A (left panel), B (middle panel), and C (right panel).K(k) is shown by solid, black lines while M (k) is denoted by red, dash-dotted lines.The wave number is normalised to k min = 2π.The thin vertical line in each panel shows the turbulent forcing wave number.

Figure B1 .
Figure B1.Evolution of the PDFs of the stretching term in panel (a), the advection term in panel (b), the compression term in panel (c) and the dissipation term in panel (d) for run A with Mrms = 0.2 (subsonic).The values of µ and σ shown in the legends corresponds to the mean and standard deviation of the distribution in each case.

Table 1 .
Key parameters of simulations used in this study.The resolution in each case is 512 3 .k f L/2π = 2 is the average forcing wave number and Mrms is the average values of the rms Mach number in the steady state.ℓ f = 2π/k f is the forcing scale and Pm and Re are the magnetic Prandtl number and the fluid Reynolds numbers, respectively.