# Direct numerical simulations of the galactic dynamo in the kinematic growing phase

Christophe Gissinger, Sebastien Fromang, Emmanuel Dormy

► To cite this version:

Christophe Gissinger, Sebastien Fromang, Emmanuel Dormy. Direct numerical simulations of the galactic dynamo in the kinematic growing phase. Monthly Notices of the Royal Astronomical Society, 2009, 394, pp.L84-L88. 10.1111/j.1745-3933.2009.00616.x. insu-03604909

HAL Id: insu-03604909

https://insu.hal.science/insu-03604909

Submitted on 11 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Direct numerical simulations of the galactic dynamo in the kinematic growing phase

Christophe Gissinger,<fnref n="1" />* Sebastien Fromang<fnref n="2" />,<fnref n="3" /> and Emmanuel Dormy<fnref n="4" />

Accepted 2008 December 23. Received 2008 December 23; in original form 2008 November 4

## ABSTRACT

We present kinematic simulations of a galactic dynamo model based on the large-scale differential rotation and the small-scale helical fluctuations due to supernova explosions. We report for the first time direct numerical simulations of the full galactic dynamo using an unparametrized global approach. We argue that the scale of helicity injection is large enough to be directly resolved rather than parametrized. While the actual superbubble characteristics can only be approached, we show that numerical simulations yield magnetic structures which are close to both the observations and the previous parametrized mean field models. In particular, the quadrupolar symmetry and the spiraling properties of the field are reproduced. Moreover, our simulations show that the presence of a vertical inflow plays an essential role to increase the magnetic growth rate. This observation could indicate an important role of the downward flow (possibly linked with galactic fountains) in sustaining galactic magnetic fields.

**Key words:** magnetic field – MHD – ISM: bubbles – galaxies: ISM – galaxies: magnetic fields.

## 1 INTRODUCTION

It is widely accepted that magnetic fields of planets, stars and galaxies are generated by dynamo action, i.e. by the magnetic field amplification due to electromagnetic induction associated with the motion of an electrically conducting fluid (Moffatt 1978). The flow of gas in the interstellar medium appears to convey the essential ingredients for such dynamo action (Rosner & Deluca 1989; Wielebinski 1990). A differential rotation in the galactic disc creates a strong shear along the radial direction. This shear is very efficient at stretching radial magnetic field lines in the azimuthal direction (this is known as the ω-effect). In combination with this large-scale effect, the turbulent motions at small scales provide a cyclonic flow generating poloidal magnetic field (this is the so-called α-effect). Together, both effects suggest the possibility of an α-ω type of dynamo that might be responsible for generating the galactic magnetic field (Parker 1971; Vainshtein & Ruzmaikin 1971). Let us note that alternative models for dynamo action in galaxies have been proposed through the action of cosmic rays (Hanasz et al. 2004) or in a cosmological context (Wang & Abel 2007), which will not be discussed here. The apparent scale separation between the shear and the turbulent motions has often been invoked to introduce a mean field approach for the galactic dynamo (Beck et al. 1996; Ferriere 1998). In such formalism, an equation only for the large-scale magnetic field is

solved, the effect of small scales being parametrized by an α term (Krause & Raedler 1980). Relying on mean field equations has proven to be a very efficient approach to the galactic dynamo problem (Ferriere 1992). It is, for example, an efficient way to achieve moderate simulation time. However, the results of mean field simulations are intrinsically limited by strong assumptions such as scale separation or the statistical properties of turbulence. It is thus interesting to study galactic dynamos with direct simulations of the full problem by properly treating the small-scale flow associated with the turbulence in the interstellar medium and thus solve for the magnetic field at all scales.

It is often assumed that the most important source of turbulence in the interstellar medium comes from supernova explosions (McCray & Snow 1979). The positions of these explosions are not completely random in the disc but they often occur in cluster. This produces giant expanding cavities of gas known as superbubbles. These explosions occurring in a rotating galaxy, the expansion is affected by a Coriolis force. This yields cyclonic motions and thus a strong helicity in the gas flow (Ferriere 1998). In such a framework, however, it is worth noting that the scale separation mentioned above is not dramatic. Superbubbles have typical sizes of the order of a few hundreds parsec (see Oey & Clarke 1997). This is smaller, but not dramatically smaller than the typical vertical scale of the galaxy (~1 kpc). Given modern day computational resources, these numbers suggest that direct numerical simulations (i.e. numerical simulations that do not rely on an ad hoc parameterization of the small scales, e.g., through the α-effect) are within reach. Indeed,

<footnote n="1">¹LPS, CNRS UMR 8550, Ecole Normale Supérieure, 24 Rue Lhomond, 75231 Paris Cedex 05, France</footnote>

<footnote n="2">²DAMTP, University of Cambridge, Centre for Mathematical Sciences, Wilberforce Road, Cambridge CB3 0WA</footnote>

<footnote n="3">³CEA, Irfu, SAp (UMR AIM), Centre de Saclay, F-91191 Gif-sur-Yvette, France</footnote>

<footnote n="4">⁴MAG (IPGP/ENS), CNRS UMR 7154, LRA, Ecole Normale Supérieure, 24 Rue Lhomond, 75231 Paris Cedex 05, France</footnote>

<footnote>*E-mail: christophe.gissinger@lps.ens.fr</footnote>
Gressel et al. (2008) recently presented such simulations. To cope with the large resolution still needed to address this problem, they adopted a local approach based on the shearing box model. Their results indicate a good agreement of the local approach with mean field models. However, the local approach they used precludes any global diagnostics, such as the global structure of the field, to be established.

The purpose of this Letter is to present such global numerical simulations, resolving the magnetic field at all relevant scales in the galaxy (i.e. from 100 pc to 10 kpc). To reduce the computational burden that would be associated with full magnetohydrodynamics simulations, we work in the kinematic regime: we solve the induction equation using a prescribed and time-dependent gas flow. The latter is set by using an analytical velocity field which intends to reproduce the large-scale shear associated with rotation and the effect of superbubble explosions on the interstellar medium.

## 2 NUMERICAL MODEL

The direct numerical simulations presented in this Letter are fully three-dimensional. We solve the induction equation governing the evolution of the solenoidal magnetic field **B** in a cylindrical coordinate system ($r, \phi, z$):

$$ \frac{\partial \mathbf{B}}{\partial t} = \nabla \times (\mathbf{u} \times \mathbf{B}) + \mathrm{Rm}^{-1} \Delta \mathbf{B}, \qquad (1) $$

written in a dimensionless form using the advective time-scale. The magnetic Reynolds number $Rm$ is defined as $Rm = \mu_0 \sigma r_0 U_0$, where $\mu_0$ is the permeability of vacuum, the typical length scale $r_0 = 10$ kpc is the radius of the galactic disc, the typical velocity scale $U_0$ is the velocity of the large-scale flow (i.e. the differential rotation of the galaxy) and $\sigma$ is the conductivity of the plasma.<fnref n="1" /> In our simulations, the vertical extent of the galactic disc is $H = r_0/10$, $z$, thus range from $-0.05r_0$ to $+0.05r_0$. We restrict our attention here to the kinematic problem, ignoring the back reaction of the magnetic field on the flow. The velocity field $\mathbf{u}$ used in equation (1) is analytical and represents the differential rotation of the galaxy and the supernova explosions. This approach also means that we do not explicitly consider important effects such as density stratification in the vertical direction or the induction effect which would be due to interstellar turbulence (Ruzmaikin, Sokolov & Shukurov 1988).

Equation (1) is solved using a finite volume approach. The method is described in details by Teyssier, Fromang & Dormy (2006): it uses the MUSCL–Hancok upwind method. The solenoidal character of the magnetic field **B** is maintained through the constrained transport algorithm (Yee 1966; Evans & Hawley 1988). We rely here on the so-called pseudo-vacuum boundary conditions for the magnetic field. This corresponds to imposing $\mathbf{B} \times \mathbf{n} = \mathbf{0}$ at all boundaries of the computational domain. These boundary conditions are not fully realistic, but they are often used in parametrized models of galactic dynamos and simple to implement. These boundary conditions are known to modify quantitative results (such as the threshold value for dynamo action) but not the global qualitative solution (Gissinger et al. 2008). We now turn to a detailed description of the velocity field being used. It is the sum of two terms: rotation around the vertical axis and modification of the flow by superbubbles. In

our simulations, we use the following prescription for the rotation: $U = U_0 e_\phi$, with a constant $U_0$. This is a good approximation since the angular velocity is observed to be roughly proportional to $1/r$ in galaxies. The effect of supernova explosions is more subtle to implement. We decided to consider the effect of superbubbles only and ignore here isolated supernovae, as the energy input of the former is largely dominant (Ferriere 1998). Considering superbubbles rather than smaller isolated supernovae yields larger scales which directly translate into resolutions affordable with modern days computing resources. Let us consider first the explosion of one superbubble, in a local spherical coordinate system ($r', \theta', \phi'$). Following the work of Ferriere (1998), we work under the simplifying assumption that each explosion remnant has a perfectly spherical shape. We thus use the simple radius evolution law (Weaver et al. 1977):

$$ r'_{\mathrm{sb}} = A r^{\nu}. \qquad (2) $$

During the expansion of each superbubble, the rotation of the galaxy yields a Coriolis force which tends to deflect the initially radial expansion and create cyclonic motions. This is an essential step in classical mean field $\alpha-\omega$ description of the galactic dynamo (Ferriere 1998). This Coriolis effect can be evaluated by solving the equation of gas motion:

$$ \frac{\partial v}{\partial t} = \mathbf{F}_e - 2\Omega \times \mathbf{v}'_{\mathrm{r}}, \qquad (3) $$

where $\mathbf{F}_e$ is a force leading to the radial expansion described by equation (2). Integrating equation (3) in the radial direction leads to the expansion equation (2). The azimuthal velocity is obtained by integrating the equation (3) in the azimuthal direction. In doing so, we made the approximation that the Coriolis force on the superbubble is only due to the radial expansion of the shell. Inside the superbubble, we assume a linear variation of velocity in radius. An important parameter is $r'_c$, the critical size reached by the superbubble, for which the pressure in the cavity becomes comparable to that of the surrounding medium. At this point, we consider that the bubble merges with the interstellar medium. This situation generally occurs when the radial velocity of the shell becomes comparable to the velocity of sound in the medium. In our modelling, this critical velocity numerically determines the end of existence of a superbubble. The velocity field associated with a superbubble therefore vanishes when the radial velocity reaches this critical velocity $v_c$. This radial expansion and the associated Coriolis force totally determine the flow at small scales. In most observed galaxies, the spatial distribution of explosions in the galaxy is rapidly decreasing away from the mid-plane of the disc. For simplicity, we will assume here that all explosions occur in the mid-plane only, but with random position in the disc. In actual galaxies, there is a large observed dispersion of data about superbubbles, yet averaged values for the explosions rate of superbubbles are $f_0 = 4.5 \times 10^{-7} \mathrm{kpc}^{-2} \mathrm{yr}^{-1}$ (Elmegreen & Clemens 1985). Such parameters, however, are still out of reach of present computations (especially because of the high explosion rate which implies large numbers of superbubbles to be handled at the same time). We use here a lower rate of superbubbles, but more powerful explosions, thus leading to a similar helicity input. In the simulations reported here, $f = f_0/50$, $r'_c = 0.4$, $A = 0.35$ and $\nu = 0.6$. This corresponds to about 150 superbubbles expanding in the galactic disc at a given time in the simulations. In some cases, we will also take into account a downward flow. Due to the simplicity of our model, this velocity could be attributed to turbulent diamagnetism (Sokoloff & Shukurov 1990) or to the galactic fountain mechanism (Shapiro & Field 1976; Bregman 1980). As an attempt to describe these effects, we add the following vertical

<footnote n="1">¹ We could have adopted alternative definitions of the Reynolds number, for example $Rm' = \mu_0 \sigma H V_s/2$, where $V_s$ is the sound velocity, which is equal to the terminal velocity of the superbubbles (see later in the text). With this definition, $Rm' = Rm/200$ and the maximum value achieved in this work would be $Rm' = 500$.</footnote>
velocity to the flow:

$$v_z(z) = \frac{-\gamma z}{\sqrt{2\pi 2\beta}} e^{-z^2/2\beta^2}. \quad (4)$$

It is antisymmetric with respect to the mid-plane and vanishes for $z = 0$. Moreover, the infall velocity decreases far away from the mid-plane. The parameter $\beta$ controls the extension of the infall region and we use here $\beta = r'_c/3$ so that the maximum of the infall is near the region where superbubble explosions tend to accumulate the matter. $\gamma$ is a free parameter controlling the amplitude of the vertical velocity. We will use $\gamma = 0.03$ throughout this Letter corresponding to a typical velocity of 6 km s$^{-1}$. Despite the simplifications implied by working in the kinematic regime, large spatial resolutions are still needed in order to correctly describe the evolution of the superbubbles at small scales. In the runs presented here, we used a resolution of $N_r = 200$, $N_\phi = 640$ and $N_z = 36$.

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Growth rate of the magnetic energy as a function of the magnetic Reynolds number. Black circles correspond to simulations without infall ($\gamma = 0$), whereas the filled red square corresponds to a simulation with $\gamma = 0.03$.</figcaption></figure>

# 3 RESULTS

## 3.1 General features

We performed simulations for seven different magnetic Reynolds number ranging from $Rm = 100$ to $10^5$ (this would correspond to $Rm'$ between 0.5 and 500). We choose to stop the simulations after a few resistive times, when the growth rate of the magnetic energy is statistically invariant and the exponential growth is well established.

For all of these simulations, we measure the growth rate of the magnetic energy. It is displayed on Fig. 1 as a function of the magnetic Reynolds number, $Rm$. It is negative when the magnetic Reynolds number, $Rm$, is smaller than $Rm_c \sim 500$. It is positive for larger $Rm$, indicating exponential amplification in that case. For $Rm = 10^5$, the growth rate is $\sigma = 0.6$ Gyr$^{-1}$. Such growth times are comparable to the ones obtained by Gressel et al. (2008), although they seem to be larger in our case.

The result of a typical simulation ($Rm = 10^5$) once the exponentially growing phase is reached is illustrated in Fig. 2, which shows simultaneously the structure of the magnetic field and that of the flow. Many superbubbles (red isosurfaces) are present at a given time in the model. We also show field lines (plotted in blue) of the magnetic field. The observed magnetic structure is the results of the combined effects of the superbubble explosions and the differential rotation of the disc. The coloured slice shows the magnetic energy in the equatorial plane. It is strongly fluctuating due to the complicated nature of the flow. The overall topology of the magnetic field is complex. We now turn to a detailed study of its structure.

## 3.2 Structure of the magnetic field

The structure in the $(r, \phi)$ plane is complicated and varies with the altitude $z$. Fig. 3 shows the magnetic field in the mid-plane of the galaxy (the solid lines represent field lines projected in this plane and the colour code indicates the strength of $B_z$).

In the mid-plane, the field is organized in a spiral structure. The sign of $B_\phi$ is constant along the radial direction. Near the axis of

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Field lines of the magnetic field (blue) are represented for $Rm = 7000$. Both the spiral structure of the magnetic field and its quadrupolar symmetry can be identified. Isosurface corresponding to 1 per cent of the peak kinetic energy is also represented (red).</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Structure of instantaneous magnetic field in (r, φ) plane in the mid-plane of the galaxy at z = 0. Magnetic field lines projected in the (r, φ) plane are represented by black lines and the colour code reflects the strength of Bz.</figcaption></figure>

rotation of the disc, the azimuthal field largely dominates all others components, but rapidly goes to zero at the inner boundary in order to satisfy the boundary conditions. At larger radii, the vertical field is negligible while $B_φ$ and $B_r$ are now comparable. Their relative value is given by the magnetic pitch angle defined as $p_B = atan(B_r/B_φ)$. It is remarkable that, despite the fact that numerical parameters are far from actual values, $p_B$ is very close to the observations: except near the unrealistic boundaries of the domain, the pitch angle is, in general, close to −15°, which is in agreement with the range [−30°, −10°] observed in real galaxies (Shukurov 2007). An average of the pitch angle in radius from r = 0.2 to 0.8 gives $p_B \simeq -15°$. This is also in agreement with Gressel et al. (2008) although smaller, as they report a pitch angle around −10°.

At higher altitudes, the structure of the field is much more complicated. By increasing z, we observe that $B_φ$ can change sign. We always observe opposite sign of $B_φ$ between the mid-plane (z = 0) and the halo (z = +z₀) (see Fig. 4). At intermediate altitudes, $B_φ$ can also reverse sign along the radial direction itself, as it is shown in Fig. 5.

While the structure in the (r, φ) plane is not very sensitive to the resistivity, the magnetic field in the (r, z) plane presents different behaviours depending on the value of the Reynolds number, Rm, as shown in Fig. 6. A quadrupolar structure is ubiquitous in all simulations but the location of the magnetic loops does depend on Rm. Indeed, the effect of superbubbles is located near the mid-plane and produces strong expulsion of magnetic field in the halo of the galaxy. For weak Rm, the magnetic resistivity counteracts this effect through vertical diffusion. For larger Rm, the weak magnetic resistivity cannot balance anymore the strong vertical expulsion of magnetic field due to multiple explosions. As a consequence, the quadrupole becomes unrealistically confined to the halo of the galaxy (Fig. 6c), far from the active region. Although the creation of this external shell does not totally inhibit dynamo action, it clearly decreases the magnetic growth rate.

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Magnetic field just before the top of the domain. Note that the magnetic field is smooth due to the weak effect of the superbubbles at this altitude. The sign of $B_φ$ is reversed compared to the mid-plane.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Magnetic field at z = z₀/2. Note that $B_φ$ changes sign when the radius is increased.</figcaption></figure>

### 3.3 Effect of vertical infall

This behaviour indicates how the diffusion of magnetic field can play two opposite roles: on one hand, it is obviously defavourable to dynamo action by increasing resistivity in the induction equation. On the other hand, diffusion can be favourable by preventing magnetic flux expulsion away from the mid-plane region where the small-scale flow is important. However, for weak resistivity (as is the case in real galaxies), superbubbles expel the magnetic field out of the active region of galactic disc, thus inhibiting dynamo action.
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Magnetic structure of *φ*-averaged magnetic field in the (r, z) plane for: Rm = 1000 (a), Rm = 3000 (b) and Rm = 10⁴ (c and d). Flux expulsion is clearly visible on plot (c). Plot (d) corresponds to a simulation including vertical infall. Flux expulsion can thus be counteracted.</figcaption></figure>

In that case, adding a vertical inflow by using equation (4) proved to be an essential ingredient to dynamo action. This vertical flow indeed pumps the magnetic field from the halo to the mid-plane, which increases considerably the growth rate of magnetic energy. For Rm = 10⁴, for example, the growth rate increases from σ = 0.7 Gyr⁻¹ without inflow to σ = 1.4 Gyr⁻¹ with vertical inflow (filled red square on Fig. 1). As seen on Fig. 6(d), the magnetic field structure is again quadrupolar in that case and spread out over the whole galaxy.

## 4 CONCLUSION

We have shown that according to our simple model, it is possible to perform numerical simulations of the galactic dynamo without the need for a mean field formalism. We thus avoid assumptions in the scale separation and can control more rigorously the origin of the source term in the induction equation. Our simulations yield magnetic field with two main characteristics: a quadrupolar symmetry in the (r, z) plane and a roughly axisymmetric spiral configuration in the (r, φ) plane. Both characteristics are in good agreement with observations and confirm previous studies that used a mean field approach. A detailed study of the magnetic field topology shows a complicated structure, with reversals of B<sub>φ</sub> along the radial or vertical directions. Another interesting feature of this work are the paradoxical role of superbubbles in the limit of very weak magnetic diffusion. Indeed, the turbulent flow due to explosions is, with the differential rotation, an essential ingredient of the α-ω dynamo but

also inhibits dynamo action by confining the magnetic field in the halo of the galaxy. In this context, the vertical inflow of interstellar gas appears as the third main ingredient needed for dynamo action. The downward flow observed in galaxies could thus be an essential mechanism of galactic dynamo theory.

## ACKNOWLEDGMENTS

Computations presented in the article were performed on the IDRIS and CEMAG computing centres. We are grateful to Romain Teyssier for enlightening discussions and Anvar Shukurov for useful comments.

## REFERENCES

Beck R., Brandenburg A., Moss D., Shukurov A., Sokoloff D., 1996, ARA&A, 34, 155

Bregman J. N., 1980, ApJ, 236, 577

Elmegreen B. G., Clemens C., 1985, ApJ, 294, 523

Evans C., Hawley J., 1988, ApJ, 33, 659

Ferriere K., 1992, ApJ, 391, 188

Ferriere K., 1998, A&A, 335, 488

Gissinger C., Iskakov A., Fauve S., Dormy E., 2008, Europhys. Lett., 82, 29001

Gressel O., Elstner D., Ziegler U., Rüdiger G., 2008, A&A, 486, L35

Hanasz M., Kowal G., Otmianowska-Mazur K., Lesch H., 2004, ApJ, 605, L33

Krause F., Raedler K.-H., 1980, Mean-field Magnetohydrodynamics and Dynamo Theory. Pergamon Press, Oxford

McCray R., Snow T. P., Jr, 1979, ARA&A, 17, 213

Moffatt H. K., 1978, Magnetic Field Generation in Electrically Conducting Fluids. Cambridge Univ. Press, Cambridge, England

Oey M. S., Clarke C. J., 1997, MNRAS, 289, 570

Parker E. N., 1971, ApJ, 163, 255

Rosner R., Deluca E., 1989, in Morris M., ed., Proc. IAU Symp. 136, The Center of the Galaxy. Kluwer, Dordrecht, p. 319

Ruzmaikin A. A., Sokolov D. D., Shukurov A. M., 1988, Magnetic Fields of Galaxies. Kluwer, Dordrecht

Shapiro P. R., Field G. B., 1976, ApJ, 205, 762

Shukurov A., 2007, in Dormy E., Soward A. M., eds, Galactic Dynamos. CRC Press, Boca Raton, p. 313

Sokoloff D., Shukurov A., 1990, Nat, 347, 51

Teyssier R., Fromang S., Dormy E., 2006, J. Comp. Phys., 218, 44

Vainshtein S. I., Ruzmaikin A. A., 1971, Astron. Zh., 48, 902

Wang P., Abel T., 2007, preprint (arXiv:0712.0872)

Weaver R., McCray R., Castor J., Shapiro P., Moore R., 1977, ApJ, 218, 377

Wielebinski R., 1990, Geophys. Astrophys. Fluid Dyn., 50, 19

Yee K., 1966, IEEE Trans. Antennas Propag., 14, 302

This paper has been typeset from a TeX/LA TeX file prepared by the author.