The small-scale dynamo in a multiphase supernova-driven medium

Magnetic fields grow quickly, even at early cosmological times, suggesting the action of a small-scale dynamo (SSD) in the interstellar medium (ISM) of galaxies. Many studies have focused on idealized, isotropic, homogeneous, turbulent driving of the SSD. Here we analyze more realistic simulations of supernova-driven turbulence to understand how it drives an SSD. We find that SSD growth rates are intermittently variable as a result of the evolving multiphase ISM structure. Rapid growth in the magnetic field typically occurs in hot gas, with the highest overall growth rates occurring when the fractional volume of hot gas is large. SSD growth rates correlate most strongly with vorticity and fluid Reynolds number, which also both correlate strongly with gas temperature. Rotational energy exceeds irrotational energy in all phases, but particularly in the hot phase while SSD growth is most rapid. Supernova (SN) rate does not significantly affect the ISM average kinetic energy density. Rather, higher temperatures associated with high SN rates tend to increase SSD growth rates. SSD saturates with total magnetic energy density around 5% of equipartition to kinetic energy density, increasing slightly with magnetic Prandtl number. While magnetic energy density in the hot gas can exceed that of the other phases when SSD grows most rapidly, it saturates below 5% of equipartition with kinetic energy in the hot gas, while in the cold gas it attains 100%. Fast, intermittent growth of the magnetic field appears to be a characteristic behavior of SN-driven, multiphase turbulence.


INTRODUCTION
Magnetic fields pervade the Universe at all scales.Many astrophysical systems consist of plasma, in which the highly turbulent motions drive small-scale dynamos (SSD) that rapidly grow magnetic fluctuations at the scales characteristic of the turbulence.Such magnetic fields influence the structure and dynamics of, for example, star forming molecular clouds (Mac Low & Klessen 2004;van Loo et al. 2012;Federrath & Klessen 2012;Sridharan et al. 2014) and spiral arms (Chandrasekhar & Fermi 1953;Beck et al. 1996;Fletcher et al. 2011).Magnetic fields also grow at larger scales relevant to the shape and structure of their host, such as the polar fields in stars aligned to the rotational axis, or galactic fields aligned to the spiral arms or disk (e.g., Harnett et al. 2004;Fletcher et al. 2011;Beck 2016).Such fields are generated by large-scale dynamos (LSDs).To model galactic LSDs self-consistently (e.g., Gressel et al. 2008;Hanasz et al. 2009;Gent et al. 2013b) requires its evolving entanglement with the SSD be included.The vast separation in scale and growth time makes this challenging.
Amplification and decay of magnetic fluctuations in highly compressible fluids can occur in-dependent of the presence of an SSD through a process called tangling, where a large-scale field is pushed around by the turbulent flow, and a fluctuating contribution is generated (Rogachevskii & Kleeorin 2007;Karak & Brandenburg 2016).In regions of SN compression, magnetic field strength scales with density, the exponent of its proportionality depending on the geometry of the compression.However, observations of the diffuse ISM show no correlation between field strength and density (e.g., Crutcher et al. 2003).Observations of galaxies indicate that the turbulent magnetic field strength is typically larger than that of the large-scale fields.The Solar neighbourhood random field, for example, is about 1.3 times larger than the local regular field (Beck et al. 2003).These fields are roughly in equipartition with the estimated turbulent kinetic energy density.
It is somewhat uncertain how LSD, tangling and SSD interact and contribute to this picture.Based on a suite of numerical experiments, Karak & Brandenburg (2016) reported that while tangling, as expected, is positively correlated with the large scale magnetic fields, the SSD shows an anti-correlation when the mean component of the magnetic field becomes strong.Super-equipartition mean fields, which could arise in presence of fluxes of small scale magnetic helicity, tend to suppress SSD.Some properties of tangling produced magnetic fluctuations are discussed in Gent et al. (2021) where it was noted that tangling produces only a linear growth for a given background mean magnetic field.If the background field itself grows exponentially, the associated fluctuations by tangling are also expected to show an exponential growth.However, in the absence of any regular or mean magnetic field, as is the case in the present work, exponentially growing so-lutions for small scale magnetic fields must be attributed to SSD.
The high resolution simulations of Bhat et al. (2016) include SSD and LSD simultaneously, but only for isothermal, helically-driven turbulence at 0.1 < Pm < 10.Galaxy simulations including halo-disk scale flows (e.g., Rieder & Teyssier 2016;Steinwandel et al. 2019) find SSD but capture no LSD.Their multiphase structure is parameterized, not evolved explicitly.While lower resolution models of LSD in SN-driven turbulence with galactocentric differential rotation do not include SSD, Gent et al. (2013b) appeared to do so in a non-isothermal model, as confirmed by Gent et al. (2021).We do not yet investigate the interaction of the LSD and SSD here, but examine only the properties of the SSD.Seta & Federrath (2022) model SSD in a twophase ISM (cold and warm) using large-scale momentum injection to drive the turbulence, rather than point-source thermal injection as primarily applied here.In the case of both solenoidal and compressive forcing they find a large dispersion in the correlation of magnetic field strength to gas density (see their Figure 5).Compressive forcing yields a perceptible trend of B ∝ ρ 0.5 for warm gas, and B ∝ ρ 0.7 for cold, consistent with the relations for compression along field lines and spherically, respectively.SSD in the cold and warm phases has the same growth rate in Seta & Federrath (2022).Their cold gas has typically higher M, so, based on their previous isothermal models (Seta & Federrath 2021), SSD should be slower.Overall growth is slower in the multiphase medium than in an isothermal gas of similar mean M. Separation into cold and warm phases driven by thermal instability (Field et al. 1969;Sánchez-Salcedo et al. 2002;Brandenburg et al. 2007;McCourt et al. 2012).The Seta & Federrath (2021) cooling function has a slightly reduced range of instability than what we model here.
Of greater significance is that without SN thermal energy injection they also have no hot gas in a quasi-stable third phase (McKee & Ostriker 1977).The models we present here suggest this is a crucial difference.
The SSD in SN-driven turbulence has been modelled by Balsara et al. (2004) and Gent et al. (2021).The second of these included a multi-phase ISM with fractional volumes of cold, warm and hot gas somewhat consistent with observations.Their key finding was to confirm that under the conditions prevailing in an ISM heated by SNe, the SSD is easy to excite and amplifies magnetic fluctuations up to sub-equipartition levels within a few tens of megayears.A critical issue that was not resolved, though, was to explain the erratic growth rate of the dynamo in these simulations, which was not observed by Balsara et al. (2004).
We illustrate this intermittent dynamo in Figure 1 where we reproduce data from Gent et al. (2021) Figure 3, but with the inclusion of leastsquares fits of magnetic energy density e B ∝ exp(γt).False convergence (Fryxell et al. 1991) is apparent at low resolution, although the 2 pc resolution model does exhibit a brief spell of growth at around 200 megayears.At high resolution the solutions are convergent at all times, with slightly more diffusivity apparent at 1 pc resolution.Of interest for our study is the sporadic growth and decay apparent at high resolution between 20 and 40 megayears, and to some degree also after 50 megayears.In particular, we seek to understand how the SSD in the multiphase ISM differs from solutions modelled in isothermal plasmas, and to explain how this drives such intermittency in SSD behaviour.In Gent et al. (2021) only the volume averaged growth rates γ were considered.In results presented here we include a breakdown by phase.
In Section 2 we explain the motivation and design of our experiments.In Section 3 we report on the key drivers that change SSD Volume averaged magnetic energy density e B for Models UP0r10sl, HP0r10sl, MP0r10sl and LP0r10sl (see Table 1).The normalization of e B is given by the time and volume averaged kinetic energy density e K .Fits to the exponential growth rate γ, see Eq.( 5), are shown as solid lines of the appropriate color.Time scale is (a) log, to ease comparison at early times for high resolution, and (b) linear, to show the late time exponential decay of the low-resolution models.For δx ≤ 1 pc the fits span 40 < t < 48 Myr, while for δx ≥ 2 pc they span 110 < t < 300 Myr.
growth within and between models.This includes consideration of how magnetic Prandtl number varies in the multiphase environment of the ISM.We conclude with a discussion of the significance of these results in Section 4.

MHD equations
We use the Pencil Code (Brandenburg & Dobler 2002;Brandenburg et al. 2020) to model SN-driven turbulence as described previously in Gent (2012), Gent et al. (2013a), andGent et al. (2020).We solve the system of nonideal, compressible, nonadiabatic, MHD equations with the ideal gas equation of state closing the system, assuming an adiabatic index (ratio of specific heats) c p /c v = 5/3.Most variables take their usual meanings; a list of notations used is given in Table 2. Terms containing ν 6 , χ 6 and η 6 apply sixthorder hyperdiffusion to resolve grid-scale instabilities (see, e.g., Brandenburg & Sarson 2002;Haugen & Brandenburg 2004), with mesh Reynolds number set to be 1 at the scale of the zone size δx.Terms including ζ D , ζ ν and ζ χ resolve shock discontinuities with artificial diffusion of mass, momentum, and thermal energy, respectively.They depend quadratically on the local strength of the shock (see Gent et al. 2020, for details).Equations ( 2) and ( 4) include momentum and energy conserving corrections for the artificial mass diffusion ζ D applying in Equation (1).Following Gent et al.  (2021), resistive shock diffusion is omitted from Eq. ( 3).SNe at rate σ inject E th = 10 51 erg thermal energy.In dense regions up to 5% is instead injected as kinetic energy E kin (see Kim & Ostriker 2015;Gent et al. 2020).Nonadiabatic heating Γ UV and cooling Λ(T ) (as detailed in Gent et al. 2013a) follow Wolfire et al. (1995) and Sarazin & White (1987).

Model parameters
We simulate a nonstratified, nonrotating domain with initial gas number density n = 1 cm −3 covering (256 pc) 3 and periodic across all boundaries.Resolution along each edge spans 0.5 to 4 pc, corresponding to grid sizes of 64 3 to 512 3 zones.Each model SN rate is given relative to the estimated rate in the Solar neighbourhood of the Milky Way, σSn = 50 Myr −1 kpc −3 .SNe occur uniform randomly in space at times following a Poisson process.Models with common σ use the same SN schedule and locations.For ambient gas number density n 1 cm −3 and σ σSn an estimated forcing scale of 60-100 pc (Joung & Mac Low 2006;de Avillez & Breitschwerdt 2007;Hollins et al. 2017) should support at least 2-4 turbulent cells.For higher SN rates the forcing scale reduces, increasing the number of turbulent cells (Joung et al. 2009), but the forcing scale remains unchanged for the lower SN rates applying here, where an individual SN at gas number density 1 cm −3 W traceless rate of strain tensor 6)  6th order rate of strain tensor merges with the local sound speed within a 70 pc radius (Gent et al. 2020) In addition, Model LP5r8sm-S, is vertically stratified without rotation and has open vertical boundaries.Initially n 1 cm −3 at the midplane.With this we examine the effect that vertical expansion of SN remnants or superbubbles has on the SSD.Here Γ UV is amplified to 3.5× that of Wolfire et al. (1995) to support the thickness of the disk in the absence of ionization heating and cosmic ray pressure gradients (see Hill et al. 2018).The domain size of this model is 512 pc along the disk and ±1.534 kpc perpendicular to it.SN are located uniform random horizontally and normal random with scale heights of 90 pc and 325 pc for Type I and Type II SNe, respectively (Ferrière 2001).This model otherwise is as described in Gent et al. (2013a,b).
All models are listed in Table 1 where the model labelling convention is explained.

Averaging conventions
Angular brackets indicate the quantity is averaged over the volume, or also with a subscript T when over the volume of individual temperature phases.An overbar indicates averaging over a domain and an interval of time, intervals which may vary, as explained in the text.
In the case of e K the interval is selected for each model to exclude initial transients and the period after the SSD approaches saturation.Even with saturation around 5% of equipartition energy, some models show damping in kinetic energy.Models with identical SN rates, have near identical kinetic energy density evolution, so direct comparisons are not affected by this choice.The magnetic energy density e B is then normalized by e K to ease model comparison.

Growth rates and Reynolds number
The erratic growth of the volume averaged magnetic energy shown in Figure 1 and dis-cussed in Sect. 3 shows that the kinematic stage of the SSD does not have a well-defined single growth rate (Gent et al. 2021).SSD growth typically follows an exponential of the form e B (t) = e B (t 0 ) exp[γ(t − t 0 )]. (5) In Figure 1 we fit such a function to specified time intervals for each model, but the SSD growth (or decay) at other times differs, with γ varying or no clear exponential behaviour occurring.Indeed within the domain, we can find very different growth patterns within and between thermal phases.
To interpret how the time and volume averaged growth rate γ is influenced by the phases and dynamical properties, we can use Eq. ( 3) to identify instantaneous changes to the magnetic energy.Hyperdiffusion, being purely numerical, is neglected.Taking the curl of Eq. ( 3) and contracting it with Bµ −1 0 , we obtain an equation for the change of the magnetic energy density Negative values represent decay of the magnetic field and positive its amplification.Dividing by e B , we obtain an equation for the relative growth rate exponent, of the form e B ∝ exp(Γt), at each instant in time and location in space Here, Γ(x, t) is a function of position x, distinct from γ defined above, which is the volume and time averaged quantity.Statistical analysis of Γ in relation to various physical properties will help determine how the varying growth depends on the multiphase structure of the ISM.
The growth rate Γ does not indicate the absolute magnitude of the energy change.We also need to identify where the largest changes in magnetic energy density occur, since these need not correlate with the growth rates.We therefore also define the absolute growth rate replacing e B in the denominator with the time dependent volume averaged e B (t). Rescaling by time dependent e B assists comparison between all stages of the SSD.Both Γ and Γ have values that span orders of magnitude of both signs.To take advantage of logarithmic scales in our histograms we omit negative and negligible growth rates.
Similarly, a field of values for the magnetic Reynolds number Rm is calculated directly from the induction equation by taking the ratio of the advection to the diffusion terms, where In a related approach, Evirgen et al. (2019) decompose the terms of Equation ( 2) to identify separately the spatial variation of each force.

Phase fractional volume
There are various approaches to measuring the proportion of the ISM that contains gas in each phase, as discussed in detail within Gent et al. (2013a).Here, we measure the fractional volume of SSD activity for each phase i defined where V i is the volume occupied by the gas in the temperature range defining phase i, as listed in Table 3, and V is the total volume.Both the phase volume and the total volume omit locations with negative or negligible Γ to focus on SSD growth.For each phase, f V,i is computed from the snapshot within the time intervals listed in Table 3.

Sorting time by growth rate
To examine how properties of the ISM in each phase differ between times with rapid growth and slow growth or decay, for each model we compute the average growth rate γ(t) as the rate of change in e B during the specified intervals from Table 3.The intervals are chosen at each resolution to exclude initial transient magnetic energy decay, prior to the onset of any SSD, and most of the subsequent saturated dynamo stage, so that we capture data of most relevance to SSD.The time series is then binned according to where γ(t) is lowest, median and highest. 1istograms for each phase from snapshots at times belonging to each time bin are then accumulated.Fractional volume of SSD activity f V,i for each phase i for a cumulative histogram is the mean of f V,i in snapshots contributing to that histogram.Summary statistics of mean log Γ or log Γ and relevant averaged physical quantities are calculated from cumulative histograms in each bin.

RESULTS
A hint toward explanation for the intermittent growth in the SSD in these models appears on inspection of slices of the simulation data as displayed in Figure 2 for Model HP2r5sl.During this period of fast growth, however, magnetic energy is high in hot diffuse gas annotated by 'hi-hot'.
The expected response of the magnetic energy to compressive flows is evident in regions highlighted by labels 'hi-cold', 'hi-warm' and 'lohot'.What is anomalous in this scenario are the regions highlighted as 'hi-hot'.In these regions the strongest magnetic field is associated with the most diffuse and hottest ISM.This cannot be explained by passive compression of the field and suggests strong SSD activity in these regions.The snapshot shown is at 20 Myr, a period in the simulation when there is a burst of magnetic field growth.Why is SSD present in these and not in other regions of hot gas, and why is field growth strong during this period and not others?For our core analysis of the multiphase structure of the SSD we primarily focus on our models with resolution of 1 pc.We include the highest resolution 0.5 pc model to demonstrate how well our solutions converge.The lower resolution runs support insights into the dependence of the SSD on ISM structure and SN rate.
We hypothesise that the erratic behaviour of the SSD is due to the changing multiphase structure of the ISM.To test this hypothesis, we compute joint histograms by thermal phase of various physical properties in the total domain from snapshots of each model alongside growth rates computed using Equations ( 7) and (8).

Vorticity
In Figure 3 we display the set of cumulative joint histograms of absolute growth Γ and (twice the) enstrophy or the norm squared vorticity |ω| 2 of cold, warm and hot gas from Model HP5r2sl.As is common across all models and variables, there is large variance and histograms overlap between the phases.However, a clear trend appears of increasing vorticity and growth rate with temperature.In the warm and hot phases there is discernible positive correlation between vorticity and growth rate within the phase.The phases and time intervals used for the histograms are listed in Table 3.  3, in this case 15 < t < 55 Myr.
The most striking correlation we find is between relative growth rate Γ and the norm squared vorticity |ω| 2 , which we show in Figure 4(a).We find that the growth rate log Γ is strongly proportional to log |ω| 2 , and vorticity increases as we move from cold to hot gas at all times.At times belonging to the high growth rate γ(t) bin (right panel) the hot gas increases in vorticity, while the vorticity in the cold and warm gas does not differ so much between bins.The fractional volume of the hot gas increases slightly for high resolution within the upper bin.Summary statistics of mean log norm squared vorticity |ω| 2 for all runs compared to mean log magnetic energy growth rates (a) relative Γ (Eq.[7]) and (b) absolute Γ (Eq.[8]).Averages use the cumulative histograms for each phase (denoted by angle brackets with T subscript; see Tab. 3).Histograms accumulate snapshots binned by the (left) lower, (center) median and (right) upper growth rate γ(t), as explained in Sect.2.6.The phase is indicated by color: cold (blue), warm (green), and hot (orange).The resolution is identified by shape as shown in the legend.The fractional volume f V,i is proportional to the symbol area.
So the efficiency of the SSD is linked to high vorticity and high local and average growth rates Γ and γ are correlated with increased vorticity in the hot gas, consistent with results from models using more idealised turbulence (Federrath et al. 2011;Achikanath Chirakkara et al. 2021) To confirm whether this correlation is reflected in the absolute gains in magnetic energy we plot the same summary for Γ in Figure 4(b).At high resolution the same trends are preserved as for Γ, with the changes in Γ most associated with changes in the characteristics of the hot gas.However for low resolution Γ is actually anti-correlated with vorticity.Increased mixing at low resolution dampens vorticity and inhibits the formation of hot gas, damping the SSD.
Energy growth at low resolution is better correlated with gas density, the cold and then warm phases exhibiting higher Γ, consistent with results from isothermal modelling, except at δx = 2 pc for the upper bin of γ(t) (right) when vorticity and fractional volume increases for the hot gas.Mean vorticity in the warm and cold phases does not change markedly between bins.
In Figure 5(a) we display the summary statistics for Mach number M against relative growth rate Γ for each bin of volume-averaged magnetic energy growth rate.Overall the expected trend of Γ reducing as M increases is visible, more so in the lowest growth rate regions (left panel).We observe that hot gas with high sound speed has low Ms and high relative growth rate Γ and cold gas the inverse.What is, perhaps, unexpected is to see that within each phase there is a trend of increasing Γ with increasing M.This can be explained by the counter effect of increased velocities, which drive higher Mach numbers for a given sound speed, but also higher Rm.So increased Γ with M should have complementary trends in Rm.
In this highly compressible system it is challenging to calculate Rm as would convention- ally be defined for weakly compressible simulations.Instead, we use the ratio of advective to diffusive terms given by Eq. ( 9).With strongly fluctuating characteristic velocities and length scales applying between the continuously shifting phases, and the very sporadic and flex-  8), vs log Rm, Eq. ( 9), for Model HP5r2sl.Otherwise as described in Figure 3.
ible scale of the SN forcing, Rm varies as a function of position and time.An illustration of the extent to which Rm varies is shown in the histograms of Figure 6.While the mean Rm is quite modest of order 16 and 45 in the cold and warm phases, there are regions in all phases where Rm > 10 3 , and in the hot phase Rm > 10 4 .For the summary statistics the log mean values are quite informative to identify the trends in SSD dependencies across simulations over time and between phases.
Results for Rm versus Γ are shown in Figure 5(b).Velocities increase with temperature, increasing Rm and supporting faster SSD, an overall trend as would be expected for Γ.Within each phase and within each γ(t) bin, however, Γ is weakly anticorrelated with Rm, and across each bin the patterns are largely unaffected for cold and warm gas, suggesting higher velocities may be associated with higher Mach numbers.Results are similar for Re versus Γ as shown in Figure 5(c) for samples from models that include explicit viscosity ν.The higher resolution models (δx = 1) show a good correlation between Re and growth rate.The lower normalization of that correlation seen in the lower resolution models may well be due to reduced kinetic energy from excess cooling in under-resolved density gradients.The comparison between the 1 pc and 0.5 pc models without viscosity suggests that the 1 pc model is well converged.For the hot gas Rm, Re, and Γ all increase for the upper growth rate bin, suggesting increased velocities are associated with rotational rather than compressive flows.
If this were the case, we would expect to see this reflected in the structure of the flow, which we consider in Figure 7 (see analysis of hydrodynamic flow by Käpylä et al. 2018).We use a Helmholtz decomposition to separate the velocity field u = u 1 + u 2 into purely solenoidal u 1 = ∇ × F and compressive flows u 2 = ∇Φ, in which F and Φ are, respectively, a vector potential and scalar potential pair from which orthogonal flows u 1 and u 2 can be derived.
Both flows appear strongly correlated with Γ, although not as clearly as vorticity in Figure 4(a).It is reasonable to expect that the magnitude of compressive and solenoidal flows would be correlated, and as such some common correlation to Γ would be evident.The clearer correlations for (b) rotational rather than (a) compressive flows indicates that the latter is coincidental.Growth rate trends in the warm and cold gas do not alter significantly between bins of γ(t) even as flow strength varies.In the hot gas Γ clearly increases when |∇ × F | increases between low (left) and high (right) γ(t) bins, especially for high resolution models.In Figure 7(c) the summary statistics are plotted of the compressive ratio (Kritsuk et al. 2007) defined by (11)

Gent et al.
This shows that the ratio of solenoidal energy to total compressive and solenoidal energy (1−R cs ) is typically between 60% (cold) and 90% (hot).At higher resolution the fraction of energy in solenoidal flow is higher.The compressive ratio reduces from cold to hot phases.The growth rate is directly correlated with the fraction of solenoidal energy.

Prandtl number variability
The multiphase structure and intermittent forcing of the ISM causes the fluid Reynolds number Re to vary just as was shown for Rm in Figure 6.This results in a wide variation of magnetic Prandtl number Pm.In Gent et al. (2021, see Figures 5[b1] and [b2]) we examined the role of Pm on the SSD in ISM simulations, by varying both viscosity ν and resistivity η.For a fixed η, the SSD did not depend on Pm and the energy at which the dynamo saturated was consistently around 5% of equipartition.For a fixed ν, on the other hand, SSD growth rate did increase with Pm for Pm ≤ 5 in the parameter space considered.The energy at saturation increased with Pm, remaining steady for Pm 5.This suggested that it is the ratio of the Lorentz force to Rm rather than Pm that is the controlling parameter (see also Oishi & Mac Low 2011, for a similar suggestion in magnetorotational dynamos).
Therefore, here we vary Pm by varying only η.Where the viscosity, ν, and the resistivity, η, are explicitly defined, we can nominally express the magnetic Prandtl number, Pm = ν/η, which is what we shall use here.Alternatively, as a function of position Pm can also vary due to the flows and length scales characteristic at different regions and times altering the effective ratios of Rm to Re.In the Appendix we determine the effective magnetic Prandtl number Pm eff in various models to examine the relevance the explicit diffusivities to the dynamics (see Figure 12).
In Figure 8(a) we show the volume averaged magnetic energy density e B for models with Pm ∈ [1, 5] with resolution δx = 1 pc.Around 20 Myr there is a burst of SSD, common to all three runs, highlighted between the two leftmost vertical dotted lines.Growth rates γ applying for 18 < t < 20 Myr depend on Pm in this range.Sensitivity to Pm reduces as Pm approaches 5, consistent with the idea that growth rates become asymptotic for higher Pm.We also fit γ for Pm = 1 at 41 < t < 44 Myr, where the growth rate exceeds that for Pm = 5 earlier.
Under most circumstances we would expect γ to be lower for the larger η, but during this period the fractional volume of the hot gas f V,i is higher (Figure 8(b)) in all models than during the earlier subinterval.A surge aligned to the rightmost vertical dotted line is also evident in the high Pm models, although it is inhibited in those by the onset of saturation.
The critical Rm above which SSD can be excited and the typical scaling relations with Rm for the rate of magnetic energy growth are not well determined.For Pm = 1 SSD actually decays during two periods, before its final surge.This could be due to Rm varying or changing conditions resulting in a lower critical Rm.
We plot the hot gas fractional volumes f V,h for the models with Pm ∈ [1, 5] in Figure 8(b).At first, the hot gas fractions f V,h for all models coincide.Despite the initial weakness of magnetic effects, the fractional volumes diverge slightly within 20 megayears.The resistive timestep differs between the models; small changes in timestep lead to the chaotic solutions adopting alternate trajectories.For example a delay of a few decades in scheduling a single SN, even at the identical location, alters the specific ambient conditions modifying the remnant evolution and subsequent dynamics, propagating into diverging trajectories between models.
Nevertheless the trends are consistent between models.In the subinterval between the two left-  3.
most vertical dotted lines, f V,h is enhanced, corresponding to a burst in SSD activity.The dip in f V,h between the middle two vertical lines corresponds to decay in the Pm = 1 model and slower growth at higher Pm, followed by another SSD burst aligned with the third vertical line and a peak in f V,h , least so for Pm = 5.The value of f V,h then rises again at 40 megayears, particularly for Pm = 1, with an accompanying boost in SSD.The described behavior supports the hypothesis that higher SSD activity is correlated with increasing fractional volume of the hot gas f V,h .

Growth trends by phase
To examine our hypothesis further, we split the evolution of the magnetic energy between phases e B T and plot these in Figure 8(c) for the models with Pm = 1 and 5.These are normalised by e K .Horizontal dotted lines of matching color for each phase indicate the time-averaged kinetic energy density by phase e KT /e K .These lines, plotted for Pm = 5, are very similar for Pm = 1.
Usually e B T is smallest in the hot gas, but during the SSD activity bursts, its growth is the most rapid and its energy density then exceeds that of the other two phases.The saturation energy of the hot gas is similar for both models, around 0.05e Kh at 40 megayears for Pm = 5, and 55 megayears for Pm = 1.In contrast the saturation energy in the cold and warm phases is affected by Pm.In the nonlinear phase the magnetic energy in the hot phase e B h decays to around 0.01e Kh , while the magnetic energy in the cold phase e B c grows, above e Kc in the case of Pm = 5.
Is the latter consistent with compression in the statistical steady MHD state?Under compression, magnetic field strength relates to density |B| ∝ ρ α , where α ∈ [0, 1] for compression along field lines up to disc or slab-like compression (Tritsis et al. 2015).In Figure 9, for the two models in Figure 8(c) we plot the mean gas number density in the cold gas (left axis) and the root mean square strength of the magnetic field normalised by cold gas density b rms C = B 2 1/2 C ρ −α C (right axis), taking α = 1/2.In the interval 50 Myr < t < 60 Myr the dynamo has clearly saturated for the warm and hot gas (Figure 8).If the continued growth of the magnetic field in the cold gas were due to compression of the field with the gas, then we would expect the ratio of b to remain steady or even decay.Instead, in both models this ratio continues to increase until ∼ 65 Myr, indicating that dynamo remains active in the cold phase.The increasing mean density in the cold gas toward the end of the simulation reflects the increasing fractional volume of hot gas and mean thermal energy density.The cooling in the hot gas is too inefficient to balance the supply of thermal energy from subsequent SNe in the closed system, and thermal runaway will eventually result.Applying 2/3 < α < 1 in the analysis does not qualitatively alter the development of b.

Supernova rate
The strength of the SSD in SN-driven turbulence has been found to be sensitive to the SN rate, σ (Balsara et al. 2004;Gent et al. 2021).Balsara et al. (2004) used an SN rate of 8 ≤ σ ≤ 40 σSn .Because of their confinement within the periodic boundaries, these models were dominated by thermal runaway, as described by Kim & Ostriker (2015), in which the fractional volume of hot gas reaches f V,h > 0.9.Therefore, SSD has a steady exponential growth, because the ISM has effectively become a homogeneous, single phase medium.Indeed, in similar experiments by Gent et al. (2021), thermal runaway was easily excited at high resolution (δx ≤ 1 pc) even with σ = σSn .At low resolution (δx = 4 pc), with consequentially higher numerical diffusion, and thus unphysically strong cooling, thermal runaway was still excited for σ = 10 σSn .Rapid steady acceleration of the SSD occurs, which is excluded here with lower σ.
In Figure 10(a) we compare evolution of e B for σ = 0.2 σSn and 1 σSn at δx = 2 pc and 4 pc, and (b) the fractional volume of the hot gas f V,h .The higher rate σ produces more hot gas and sustains the dynamo, more so for δx = 2 pc (blue) than δx = 4 pc (brown).For lower σ, e B decays, but growth occurs briefly in the subinterval between the first two vertical lines, while f V,h increases, particularly for δx = 2 pc (green).Least-squares fits for γ are shown for this subinterval (green and cyan).Higher f V,h occurs for higher σ with fits for γ (blue and red) indicated in the second subinterval.Another boost in SSD occurs for δx = 4 pc after 270 megayears when f V,h is at an even higher peak.In general for δx = 4 pc there is less hot gas than for δx = 2 pc at each given SN rate, because of higher numerical diffusion and cooling.
Mean kinetic energy density is not strongly affected by the SN rate (Gent et al. 2021, see their Figure 3[d]), so the change in SSD growth rate would appear to be due to the increased heating of the ISM at higher σ. Figure 10 demonstrates that this is indeed the case.High SSD growth rate γ coincides with high f V,h for hot gas within each model and between models, higher γ occurs in models that support higher f V,h .

Stratification
Anticipating future examination of the LSD, in Model LP5r8sm-S we consider how stratification affects the SSD.To exclude the LSD, we omit large-scale rotation and shear.Inhomogeneities due to non-periodic boundaries vertically might still cause an LSD effect, so any residual mean magnetic field volume and horizontal component average is subtracted at each timestep.In contrast to the periodic domains, the open vertical boundary enables gas flow and magnetic helicity flux to and from the disk.
In Figure 11(a) we plot e B for this model with σ = 0.5 σSn alongside models with σ = σSn .Given that γ always increases with σ in otherwise equivalent models, it is evident that the substantially increased growth rate in model LP5r8sm-S is due to its stratification.We also see from Figure 11 From the summary statistics in Figure 4(b), we find for δx = 4 pc that the absolute growth Γ of magnetic energy is fastest in the cold gas, but occurs mostly in the warm gas due to its higher fractional volume.In contrast the stratification supports high γ and high f V,h for the hot gas, depicted in Figures 4-7 by a large orange plus symbol.
As γ increases, the absolute growth rate Γ in the hot gas reduces, as shown in Figure 4(b) (left to right), for δx = 4 pc, While for δx ≥ 2 pc it increases.In the cold gas Γ increases for δx = 4 pc and reduces for higher resolution.At δx = 4 pc, however, γ does not vary over time as much as for higher resolution.At all stages Γ in the warm gas is five to ten times larger in the stratified model than the periodic models.
Vorticity in all phases is much lower (Figure 4) at δx = 4 pc, but of these it is highest in the stratified model.The values of M, Rm, and velocities (Figs. 5 and 7), are also highest in the stratified model.Higher velocities, and vorticity, in all phases are possible with stratification, because the gas can expand into lower pressure diffuse gas away from the midplane.

Differences with isothermal dynamo
Consensus based on studies of isothermal MHD turbulence predicts that Rm crit , the critical Rm above which SSD can be excited, reduces as Pm increases above Pm 0.1 with an asymptotic limit of order Rm crit ∼ 20 for Pm 1 (Schekochihin et al. 2002;Schober et al. 2012;Seta et al. 2020).Less directly relevant to our context of ISM turbulence, for low Pm there is a range 0.1 > Pm > 0.01 where Rm crit has a maximum (Schekochihin et al. 2005;Iskakov et al. 2007;Brandenburg 2011), and then for Pm < 0.01 Rm crit again begins to drop (Warnecke et al. 2022).Once excited, the rate of growth of SSD also increases, up to some asymptote, for higher Pm.In these isothermal systems, it has also been demonstrated that as M increases it becomes more difficult to excite SSD (Rm 60 for Pm ≥ 5 Haugen et al. 2004) and growth rates are correspondingly reduced (Federrath et al. 2011(Federrath et al. , 2014)).Hence, accord-ing to these findings, SSD can be excited in our models, as this threshold value is easily exceeded.
Superficially, these conclusions might appear irrelevant to an SSD in multiphase SN-driven turbulence with respect to the averages across all phases.The response of models with respect to Pm is inconsistent.Higher SN rates with correspondingly higher compressive forcing appear to excite SSD more easily and with larger growth rate.SSD does not appear to follow a dominant eigenmode, establishing a distinct exponential amplification until saturation.
However, detailed analysis reveals that what is inconsistent in this more realistic ISM is the distribution of temperatures, gas densities, and flow properties throughout the thermal phases and over time.Determination of Pm, Rm, Re, and M must be specific to each region of the ISM.It is likely that the variation we see by position in these properties occurs also in isothermal simulations, but that the statistics of these remain quite consistent throughout the domain and duration.In a multiphase system, on the other hand, these statistics vary by orders of magnitude as a function of position and time.
The breakdown of the models by phase and overall SSD growth rate shows that the insights from isothermal SSD experiments do apply, but contrary trends, such as increasing growth in SSD with increasing M, seen in Figure 5(a), are due to the SSD supportive counter trends, such as increasing velocity (with M).The expected effect of M on growth rate is well demonstrated in the anticorrelation going from hot to cold phase.Conversely, in Figure 5(b), Rm seems anticorrelated with SSD growth rate in the hot phase, indicating that SSD is more easily excited as the hot filling fraction and typical temperature of the hot gas increase, even for lower Rm.The expected trend of positive correlation between Rm and SSD growth rate is evident in the differences between phases.

Growth rates and strength of SSD
In their two-phase simulations, Seta & Federrath (2022) determined that the magnetic energy in their models grows more slowly than in isothermal models with the same forcing and equivalent mean Mach numbers M∼ 1 for the warm gas and 5 for the cold gas.Their SSD saturates after around 400 megayears for solenoidal forcing and 800 megayears for compressive forcing, for a grid resolution comparable to Model UP0r10sl.For this and all our 1 pc resolution periodic models, the SSD saturates within 40-60 megayears with mean M∼ 1, 0.6 and 0.25 in the cold, warm and hot gas, respectively.However, in all phases M spans around two orders of magnitude over a large proportion of the volume, and outliers extend even further.The variance is much larger than in the two-phase simulations.The faster SSD in our models may, therefore, be less sensitive to the mean M, as growth occurs predominantly in low Mach number, high vorticity regions, particularly in the hot gas.
The SSD persists longest in the cold phase and next longest in the warm phase (Fig. 8[c]).We speculate that an explanation for the growth rates differing between phases in three-phase and not in two-phase simulations is that high vorticity generated in the hot gas transfers to warm then cold regions, but with a time delay that is significantly longer than from the warm to cold phase.As a result, we find faster growth in all phases.
Estimates for the post inflationary seed magnetic field are in the range 10 −65 -10 −20 µG (Grasso & Rubinstein 2001;Subramanian 2016).With an exponent of γ = 500 Gyr −1 this could be amplified to sub-equipartition saturation of around 1µG within 170 to 590 Myr.Given these results consistently yield saturation of the SSD significantly below equipartition with the kinetic energy density (see also Schober et al. 2015), this does not exclude LSD action being sufficient to amplify even the lowest estimates of post-inflation magnetic field to the values observed in high redshift galaxies.

Interpretation of Rm
It is perhaps surprising, and inconsistent with current understanding of dynamo, to see in Figure 6 locations with Rm ≤ 1, well below isothermal predictions of Rm crit , coincident with growth rate Γ > 10 3 Gyr −1 in cold gas.These are outliers in a highly inhomogeneous multiphase structure, in which there is phase mixing and energy transfer within and between phases.The bulk of the simulation data corresponds to high Rm and high Γ with the positive correlation reflected in the increasing log mean from the cold to hot phase.
To examine whether the high apparent growth rates at low Rm can be explained without SSD, let us approximate the effect of a blast wave following Sedov-Taylor (Sedov 1959;Taylor 1950) travelling at 150 km s −1 , as would occur in an SN remnant with radius near 30 pc in a 10 4 K medium of uniform gas number density 1 cm −3 .From Equation (8) with the normalized quantity B = B µ 0 e B −1/2 we can write Considering the first term with a shock front traveling at 150 km s −1 into warm gas with sound speed around 12 km s −1 resolved in our model over three cells, −∇ • u ≈ 138 km s −1 (3 pc) −1 ≈ 4.6×10 4 Gyr −1 .If the local magnetic field is even twice the volume averaged strength, this would yield Γ 10 5 Gyr −1 .The effect of the other terms and nonadiabatic cooling would adjust this estimate up or down, but nevertheless it remains entirely consistent with the rare high values of Γ and Γ in regions of low Rm in Figure 6.
For a magnetic field perpendicular to the blast wave reaching a maximum at the shock centre, where the first derivative approaches infinity and the second derivative vanishes, but near zero either side of the shock where the inverse is true.Both can coexist with the compression Γ estimated above, explaining these low Rm outliers in the joint histograms.
Conventionally, a single characteristic Rm for a domain is related to γ, associated with the single exponential growth rate of the dominant eigenmode of the SSD.The analysis of Grete et al. (2017) found that the Lorentz forces could induce weak nonlocal energy transfers, but that magnetic and kinetic energy transfers in both directions are predominantly highly localised in k-space, including for isothermal compressible MHD.Hence, the differences in characteristic length scales between the phases would tend to yield Rm supportive of very different dominant dynamo modes.Thus, in the multiphase ISM the SSD growth rate γ represents a superposition of varying dynamo modes, in which the dominant modes can switch from one phase to another over time due to heating, cooling and phase mixing.
The mean Γ of the instantaneous localised growth rates collected while γ(t) is highest have similar magnitude to γ as fitted in each simulation, supporting the concept of a superposition.Sufficiently high mean Γ or Γ correlated with a supercritical mean Rm would dominate the overall growth rates irrespective of the outliers.

SSD dependence on vorticity
We considered the properties in the simulations of velocity, Mach number, fluid and magnetic Reynolds numbers, Prandtl number, kinetic and magnetic helicity, and kinetic energy to identify how they relate to the relative and absolute growth rates of magnetic energy Γ and Γ.By far the clearest correlation to growth rate was in the enstrophy, or norm squared vorticity.
The dependence of growth rate on velocity shows a weaker dependence, which is very similar to that of the Helmholtz decomposed flow displayed in Figure 7, and even less correlation is seen in the kinetic energy density.The latter indicates that it is the magnitude of the flow rather than the kinetic energy density that is important to the SSD efficiency.
We break down the statistics into epochs of various rates of volume averaged growth γ.Periods of higher γ have higher velocities and vorticity.From the stronger correlation between vorticity and growth, we conclude that the vorticity of the flow drives the SSD.Despite the compressive structure of the SN forcing there is, as expected, a high fraction of rotational flow present, as identified from the Helmholtz decomposition.The excess of rotational to compressive flow is indicated by the compressive ratio at less than 40% and even below 10% in the hot gas.The efficiency of vortex generation is primarily driven in the ISM by the strong baroclinicity, the cross of orthogonal gradients of pressure and density, at the SN shock interactions (see Käpylä et al. 2018).Typically the magnitude of the velocity and the proportion of enstrophy increase with temperature.Fitting data from Figure 4(a) and for the high resolution data only from Figure 7(c), we find that the local growth rate Γ ∝ |ω| and Γ ∝ Re .
A concern might exist that vorticity could be artificially generated, as has been identified at the intersection of adaptive grids of mixed mesh refinement (Plewa & Müller 2001), or due to the use of locally varying diffusive schemes (Robertson et al. 2010) at lower resolution.Diffusivity is independent of position in our models, except for its use in resolving shocks.Inspection of the energy spectra (as in Gent et al. 2021, and here in Figure 12 of the Appendix) confirms in all cases that energy peaks are far from the Nyquist frequency and with very low energy near theses scales.Hence, any transfer of energy from the flow to the magnetic field must occur at length scales independent of the shock resolution over a few cells.
The study of vortex generation by Käpylä et al. (2018) was a purely hydrodynamical study, although stratified, with similar levels of vortex generation.It had no hyperdiffusion.Given the relatively low magnetic energy here, even after saturation of the SSD, we can also exclude the generation of vorticity being due to the magnetic field.

SSD dependence on temperature
We have shown that SSD activity increases when the fractional volume of the hot gas increases, and we have shown that the magnetic energy density in the hot gas grows the most rapidly during these epochs.At other times, though, the magnetic field does not grow as quickly in the hot gas as in the warm gas.The vorticity in the hot gas is also particularly high during these rapid bursts of SSD.So, both the volume and the stirring of the hot gas increase, combining to drive strong dynamo action for a limited period.
In Kirchschlager et al. (2022), we use model UP0r10sl to examine dust-processing due to SN blast waves.The progress of individual remnants are followed in detail for isolated explosions in a dense turbulent region or a diffuse turbulent region.The effects on the magnetic field inside each remnant is instructive for the understanding of SSD in hot gas (see their Figure 1).In the case of the explosion in a dense region, the magnetic field is evacuated from within the remnant and is swept along by the the blast wave for up to 1 megayear.In the case of the diffuse region, the surrounding regions of strong magnetic field are pushed away by the blast wave, but inside the remnant the magnetic field grows, with the interior field strength orders of magnitude stronger at 1 megayear than at the Gent et al.
time of the explosion.The SSD is enhanced in the hot gas where the interaction of the blast wave with the inhomogeneous density structure excites strong, interacting reverse shocks that can persist to late times in the diffuse medium.In dense regions, on the other hand, the reverse shocks from the ambient ISM are weak relative to the blast wave and then quickly slow in the relatively high gas density remnant interior.This is consistent with our finding that the SSD is more sensitive to the vorticity than the kinetic energy, as the diffuse remnant contains relatively low kinetic energy.In our stratified model LP5r8sm-S the SSD is more active, with the ISM typically more diffuse.The SSD accelerates when the hot gas develops more solenoidal turbulence through SN remnant interactions and mergers (such as emerging superbubbles in the local and nearby galactic disks Bik et al. 2018;Sofue & Kataoka 2021).Simulated superbubbles with magnetic fields (e.g., Ferrière et al. 1991;Stil et al. 2009) have tended to apply uniform ambient fields and smooth ambient gas, such that amplification of the magnetic fields are concentrated around the outer shell due to compression, as in the case of the dense region.Superbubbles in turbulent MHD have been modelled (e.g., Korpi et al. 1999;de Avillez & Breitschwerdt 2005;Breitschwerdt & de Avillez 2006;Gressel et al. 2008;Gent et al. 2013a,b), but the local structure of the magnetic field not analysed in such detail.

Summary of results
We conclude that the SSD in the SN-driven multiphase turbulence of the ISM: a. has an average growth rate linearly correlated with the average vorticity of the flow, which is efficiently generated by baroclinic effects due to SN blast wave interactions; b. correspondingly has an average growth rate well correlated with the average magnitude of the rotational flow velocity, as well as the compressible and total velocities and the fluid Reynolds number, and anticorrelates with the compressive ratio; c. grows as a superposition of varying dynamo modes in each phase and within different regions of each phase; d. grows intermittently, predominantly during periods when the hot gas, with low sonic Mach number M, has strong enstrophy and large volume filling fraction; e. is thus most efficient in hot gas produced when supernovae explode in diffuse regions, where interacting reverse shocks can have high velocities; f. is more efficient in a stratified disk, where diffuse gas away from the disk is more easily heated and supports high velocities with low M; g. and already appears to approach asymptotic solutions at values of Rm, Re and Pm significantly below the actual values estimated to apply in the ISM, although higher resolution will be required to confirm this point (see Appendix).
Prior to the onset of SSD the magnetic field strength aligns with the gas density, consistent with turbulent mixing, and following saturation of the SSD the field redistributes back into alignment with the gas concentration.The hot gas saturates significantly below equipartition, and even below the the ISM average of around 5% of equipartition with kinetic energy density, while the cold gas saturates at close to the equipartition energy.Pm eff ≈ 1.2.For η = 10 −3 this is sufficient to obtain Pm eff < 1, which should nominally be the case for any η < ν.
For ν = 10 −3 we obtain Pm eff > 1 in all cases, due to the higher effective numerical viscosity than resistivity.Nonetheless, the difference between η = 10 −3 and η ≤ 0.005 is sufficient to explain the differences, between the Pm = 1 model and Pm = 5.However, we should be cautious about concluding whether the SSD is converging above Pm=5 until we increase resolution or reduce the effective numerical diffusivities.
Figure 1.Volume averaged magnetic energy density e B for Models UP0r10sl, HP0r10sl, MP0r10sl and LP0r10sl (see Table1).The normalization of e B is given by the time and volume averaged kinetic energy density e K .Fits to the exponential growth rate γ, see Eq.(5), are shown as solid lines of the appropriate color.Time scale is (a) log, to ease comparison at early times for high resolution, and (b) linear, to show the late time exponential decay of the low-resolution models.For δx ≤ 1 pc the fits span 40 < t < 48 Myr, while for δx ≥ 2 pc they span 110 < t < 300 Myr.

Figure 2 .
Figure2.Slices from Model HP2r5sl of (a) gas number density, (b) temperature and (c) magnetic energy density, during an epoch of rapid magnetic growth.At most times magnetic energy is concentrated in cold and warm phases, somewhat correlated with gas density, annotated as 'hi-cold', 'hi-warm' and 'lo-hot'.During this period of fast growth, however, magnetic energy is high in hot diffuse gas annotated by 'hi-hot'.

Figure 3 .
Figure3.Histograms of log absolute magnetic energy growth rates Γ (Eq.[8]) vs log norm squared vorticity |ω| 2 for Model HP5r2sl for (a) cold gas, (b) warm gas and (c) hot gas.The orange cross identifies the mean log Γ and log |ω| 2 of each distribution.The histograms are cumulative results for all snapshots from the time intervals listed in Table3, in this case 15 < t < 55 Myr.
Figure 4.Summary statistics of mean log norm squared vorticity |ω| 2 for all runs compared to mean log magnetic energy growth rates (a) relative Γ (Eq.[7]) and (b) absolute Γ (Eq.[8]).Averages use the cumulative histograms for each phase (denoted by angle brackets with T subscript; see Tab. 3).Histograms accumulate snapshots binned by the (left) lower, (center) median and (right) upper growth rate γ(t), as explained in Sect.2.6.The phase is indicated by color: cold (blue), warm (green), and hot (orange).The resolution is identified by shape as shown in the legend.The fractional volume f V,i is proportional to the symbol area.

Figure 5 .
Figure 5. Summary statistics of the correlation of mean log Γ to (a) mean log Mach number M, (b) mean log magnetic Reynolds number Rm for all runs and (c) mean log fluid Reynolds number.The statistics are obtained, and symbols and panels have the same interpretation, as in Figure 4.

Figure 7 .
Figure 7. Summary statistics of mean log relative growth rates Γ vs. mean log norm squared (a) compressive and (b) rotational velocities and (c) mean compressive ratio R cs for all runs.The statistics are obtained, and symbols and panels have the same interpretation, as in Figure 4.

Figure 8 .
Figure8.(a) Evolution of magnetic energy density e B for HP1r10sl (brown), HP2r5sl (green) and HP5r2sl (blue), normalized by e K , the hydrodynamic statistically steady average kinetic energy density at 8 < t < 32 Myr, prior to saturation of the SSD.Fits of growth rate γ for each model span a period over which growth can robustly be considered exponential: 18 < t < 20 Myr.A second fit for HP1r10sl spans 41 < t < 44 Myr.The black horizontal dotted line at 5% of e K indicates the energy at which SSD typically saturates in these models.(b) The concurrent evolution of the fractional volume f V,h of hot gas.(c) Evolution of e B T by phase T for HP1r10sl and HP5r2sl.Line color is as (b), and line styles and symbols denote thermal phase, as indicated in the legend.The values of e KT /e K for each phase are indicated by colored horizontal dotted lines, and the black horizontal line here indicates 100% of e K .Vertical dotted lines are included to identify significant epochs, which are subintervals of those intervals listed in Table3.

Figure 9 .
Figure 9.For the cold gas in models of nominal Pm = 1 and 5, as displayed in Figure 8, the mean gas number density n C (left axis) and root mean square mean magnetic field strength b rms C (right axis) normalised by the square root of mean cold gas density.The subscripts indicate the nominal Pm.

Figure 10 .
Figure10.(a) Magnetic energy density with fits for growth rate γ as listed in the legend.The fits for the δx = 4 pc models span the whole subintervals bounded by the vertical dotted lines, while the fits for δx = 2 pc begin later, at t = 90 and 175 Myr for low and high σ, respectively.Panel (b) for the same models, as listed in its legend for both panels, shows the evolution of fractional volume f V,h of the hot gas.
Figure 11.Magnetic energy density (a) with sample fits of growth rate γ as indicated in the legend.Panel (b) shows the evolution of fractional volume f V,h of the hot gas for the same models, as listed in its legend, including the stratified model.
(b) that the stratified Model LP5r8sm-S has the largest fractional volume of hot gas, although none of the low resolution periodic models have large f V,h (see also Figure 10]b]).

F
.A.G. and M.J.K.-L.acknowledge support from the Academy of Finland ReSoLVE Centre of Excellence (grant 307411), the Ministry of Education and Culture Global Programme USA Pilot 9758121 and the ERC under the EU's Horizon 2020 research and innovation programme (Project UniSDyn, grant 818665) and generous computational resources from CSC -IT Center for Science, Finland, under Grand Challenge GDYNS Project 2001062.M.-M.M.L. was partly supported by US NSF grant AST18-15461.We acknowledge the constructive criticism of the anonymous referee, which helped us to improve the paper.Pencil Code 2 (Brandenburg & Dobler 2002; Brandenburg et al. 2020) Gent et al.

Table 2 .
Meanings of variables

Table 3 .
Data ranges for sampling.

Table 4 .
Estimated effective magnetic Prandtl numbers Pm eff for models with various physical viscosity ν and resistivity η.