# Latitude dependence of convection and magnetic field generation in the cube

## Abstract

The 3D thermal convection in the Boussinesq approximation with heating from below and dynamo in the cube are considered. We study dependence of the convection intensity and magnetic field generation on the latitude in -plane approximation. It is shown that kinetic energy gradually increases from the poles to the equator more than order of magnitude. The model predicts the strong azimuthal thermal wind, which direction depends on the sign of the thermal convective fluctuations. The spatial scale of the arising flow is comparable to the scale of the physical domain. The magnetic energy increases as well, however dynamo efficiency, i.e., the ratio of the magnetic energy to the kinetic one decreases to the equator. This effect can explain predominance of the dipole configuration of the magnetic field observed in the planets and stars. The approach is useful for modeling of the magnetohydrodynamic turbulence in planetary cores and stellar convective zones.

LATITUDE DEPENDENCE OF CONVECTION

AND MAGNETIC FIELD GENERATION IN THE CUBE

M. Yu. Reshetnyak

Institute of Physics of the Earth of RAS, Moscow, Russia,

Institute of Terrestrial Magnetism, Ionosphere and Radio Wave Propagation of RAS, Moscow, Russia, m.reshetnyak@gmail.com

## 1 Introduction

3D modeling of convection and dynamo processes in the planetary cores and convective zones of the stars is a modern branch of the magnetohydrodynamic (MHD) simulations [16]. There are two approaches which are usually used. For the large-scale fields modeling the MHD equations are solved in the spherical geometry, see, e.g., [11], Inclusion of the realistic boundary conditions and distributions of the energy sources in the model let to simulate the various observable features of the magnetic fields, like the spatial and temporal spectra of the magnetic fields, the butterfly diagrams in the case of the solar dynamo, reversals of the geomagnetic field.

However similarity of observations and simulations sometimes can be misleading because the parameters used in 3D models are still quite far from the desired ones. Briefly, the main problem is a turbulence which simulation requires resolution of the small scales. The other difficulty is the anisotropy of the flow, concerned with the rapid daily rotation, as in the case of the planetary cores. The spherical models which have dense distribution of the mesh grid points near the poles, require the small time step (because of the Courant condition), and therefore do not suit to the turbulence modeling. Meanwhile, simulations in the Cartesian geometry with the homogeneous grids could be very helpful. These simulations reproduce the cascade processes between the scales, are helpful for estimates of the turbulent coefficients, and demonstrate various remarkable properties of the MHD turbulence [4].

In spite of the fact that these two approaches were used for years, the flat layer dynamo simulations did not take into account the latitude dependence to the moment. In other words the angle between the angular rotation axis and gravity was set to zero [12], [5], [9]. This choice of parameters corresponds to the geographic poles. Having in mind that the most common configuration of the magnetic field in planets and stars is the dipole, and 3D spherical dynamo models predict existence of the toroidal counterpart, concentrated at the mean latitudes, it looks tempting to include the latitude dependence in the flat layer models immediately. The lack of the papers, devoted to this problem in the geodynamo and stellar applications, looks quite surprising because in the meteorology the latitude dependence, known as the -plane approximation, was already used in the Cartesian models for a long time [13].

The other motivation of this paper is to distinguish the physical effects related to the angle between and definitely, leaving apart influence of the inner core and spherical boundaries, which can change hydrodynamics^{1}

Below we consider the standard 3D MHD model in the rapidly rotating cube. The model includes the thermal convection equations in the Boussinesq approximation with the heating from below. The liquid is conductive and at the large magnetic Reynolds numbers the fluid motions can generate the magnetic field. To solve these equations we use MPI C++ code, based on the predictor-corrector method [8]. We check how these processes depend on the angle between and . Some analytical predictions, based on the analogy with the motion of the charged particle in the electromagnetic field are considered as well.

## 2 Dynamo in the cube and numerical methods

The dimensionless dynamo equations for an incompressible fluid () in the cube of the height , rotating with the angular velocity , in the Cartesian system of coordinates have the form:

(1) |

Velocity , magnetic field (derived from the vector potential ), pressure and the typical diffusion time are measured in units of , , and respectively, where is the thermal diffusivity, is the density, the permeability, is the Prandtl number, is the Ekman number, is the kinematic viscosity, is the magnetic diffusivity, and is the Roberts number. is the modified Rayleigh number, is the coefficient of the volume expansion, is the unit of the temperature fluctuations , is the gravitational acceleration, and is the temperature profile, corresponding to the heating from below.

The unit vector defines direction of the angular velocity . The angle between and gravity , directed along the -axis, is equal to the colatitude , which is related to the latitude as .

The problem is closed with the periodical boundary conditions in the -plane. In -direction the following simplified boundary conditions

(2) |

at are used. Conditions for are the so-called pseudo-vacuum boundary conditions, correspond to the following conditions for the magnetic field: . Then the normal component to the boundary of the electric current is zero.

The system 1 was solved using the finite differences of the second-order in space and time. For approximation of the time derivative the three-layer time scheme was used:

(3) |

where denotes the time step .

For and the corresponding equations from 1 were written in the form:

(4) |

where

(5) |

and denotes the corresponding convective term. Eqs(4,5) with respect to were solved using the Gauss-Seidel method.

While using the vector potential provides divergence free of the magnetic field , incompressibility of the velocity field should be provided by some special technique. Here we use the predictor-corrector method [8], [3], introducing the intermediate velocity field by equation:

(6) |

where

(7) |

Then pressure was derived from the continuity equation by solving the another Poisson problem:

(8) |

with the Neumann boundary condition for -coordinate:

(9) |

The last step provides incompressibility of the velocity field:

(10) |

The second order up-wind scheme was used for approximations of the convective terms in the heat and the Navier-Stokes equations:

(11) |

where denotes the index of the grid step . The scheme 11 was used for -,-directions in the same way.

This approach was realised in C++ code with MPI. The whole domain was divided into , subdomains in coordinates, where the MHD equations 1 were solved. The subdomains exchanged by its boundaries at the each time step . In this paper the mesh grid in coordinates, and were used. The simulations were done at two linux workstations Intel(R) Xeon(R) CPU E5-2640 with 40 cores and 128Hb common memory at the station. Each run required 2-4 days, depending on the time step , which was in the range .

## 3 Pure convection

Before to consider the full dynamo the pure turbulent convective regime with , , was studied. Ten runs with step in the latitude were performed. After some intermediate stage solutions reached the quasi-stationary states. To measure the intensity of convection we estimated the mean over the volume kinetic energy, see Figure 1, defined as

(12) |

The considered regimes correspond to the developed turbulence with the Reynolds number . All the energies increase from the poles to the equator. Only at the equator has sharp minimum. This behaviour is expectable because at the equator, , -components of the Coriolis and Archimedean forces are zero, and the non-zero value of is provided by the non-linear term in the Navier-Stokes equation only. The ratio of the maximum and minimum of is . This effect is quite strong and should be explained in some way.

Moreover, analysis of the spatial structure of the flows, see Figure 2, reveals that the small-scale cyclonic convection, existed at the poles, changed to the large-scale convection at the equator. Note, that the Reynolds number in the latter case is larger, and the flow is “more” turbulent, but in the same time it is large-scale. It should be noted that -component at the equator is perpendicular to and it should be twisted at the small scale in the similar way, as it was at the poles.

To find the origin of the large-scale convection the analogy with the motion of the charged particle in the constant electromagnetic field, considered in the next section, is instructive.

## 4 Analogy with the charged particle moving in the electromagnetic field

Similarly to the motion of the charged particle in the constant in time and homogeneous in space electromagnetic field, see [2], we consider two limiting cases, where Archimedean force and angular rotation velocity of the system are directed along the same axis, and the other, where these forces are perpendicular.

In the first, the most studied case, which corresponds to the geographic poles, the flow is accelerated by the Archimedean force . Due to the Coriolice force any motion in the orthogonal plane to the gravity and is twisted with radius , defined by relation , i.e. , where is the velocity orthogonal to . The Rossby number for the Earth’s core . Even with , estimated from the large-scale velocity, based on the west drift of the geomagnetic field, one has in units of the liquid core’s scale. Taking into account decrease of the kinetic energy spectrum will only decrease estimate of . This estimate is valid for the large velocities with negligible viscous dissipation.

The linear analysis at the threshold of convection generation, where viscous diffusion is important, also predicts existence of the small scale in the perpendicular plane: [1], [6].

The both estimates demonstrate that the small-scale convection at the poles is a quite natural phenomenon, and it appears in the rapidly rotating objects with even at the critical Rayleigh numbers.

For the other case, which corresponds to the equator plane, let Archimedean force is still directed along the -axis and along , and the initial velocity is at the -plane. Then the trajectory of the particle remains in the same plane, and its motion is described by equations:

(13) |

where dot is for the time derivative.

Eqs(13), written in the reference system moving with the velocity along the -axis, after substitution , have the form:

(14) |

Choosing , one has circular motion in -plane with frequency . The trajectory in the original -plane is a trochoid, i.e. the superposition of the circular motion and the drift with velocity . In terms of the spherical geometry corresponds to the azimuthal direction. Addition of the initial velocity in -direction, which remains constant because of absence of forces in this direction, does not change situation.

In (14) is a fluctuation of the temperature relative to the non-convective distribution, and it can change the sign. By analogy with the motion of the charge in the magnetic field one has thermal separator, which divides the hot (, ) and cold (, ) flows.

Note that estimates of the radius of rotation in planes perpendicular to the axis of rotation coincide in the both cases. The difference is existence of a thermal wind with a velocity in the azimuthal direction in the latter case. It is this wind has large spatial scale, already detected in Figure 2 at the equator. Due to the continuity equation the other components of the velocity can also posses the large-scale counterpart, as we observe it in -flow.

## 5 Dynamo

Starting from the pure convection velocity and temperature distributions, obtained in Section 3, and the small magnetic field initial seed, the full dynamo system 1 were integrated in time up to the state where all the physical fields stabilised at the quasi-stationary regime. The corresponding estimates of the magnetic energy, defined as

(15) |

are ploted versus lattitude in Figure 3.

The behaviour of the magnetic energy is similar to the kinetic energy up to some details near the equator plane. The largest increase of the magnetic energy with the latitude demonstrates -component, which is stretched by the strong large-scale thermal wind . In contrast to the total kinetic energy, , the magnetic energy slightly decreases at the equator. In some sense it resembles behaviour of the magnetic field in the spherical models with the odd configurations of the magnetic field with respect to the equator, e.g., dipole. Increase of the magnetic energy relative to the poles can reach factor 4.

One can expect that the large-scale flow, based on the thermal wind, and the small-scale cyclonic convection have different efficiencies of the magnetic field generation. To test this hypothesis the ratio was plotted as a function of in Figure 4. It appears that in the range of is approximately constant and then decreases to the equator in one order of the magnitude. We conclude that the large-scale flow with the small gradients is less effective than the cyclonic convection with the non-zero net helicity [14].

## 6 Conclusions

Our simulations clearly demonstrate that the angle between the axis of rotation and gravity changes not only the magnitude of the mean parameters in the flat model, like energies, but the structure of the flow, its spectra, as well. These effects are already notable at the mean latitudes, where the magnetic energy maximum is localised in the spherical models, and should be taken into account in future. Of course, due to variety of the physical effects, application of these results to the spherical shells, should be done carefully. Thus existence of the well-known in geodynamo since [10], see also [15], the large-scaled vortexes in the Taylor cylinder at the large Rayleigh numbers, introduce the new complexity in the model.

The other issue is improvement of the numerical methods for the MHD turbulence modling. The finite difference methods are the most promising approach for the multi-core simulations. The modern higher order approximations of MHD differential equations are already comparable by accuracy to the spectral methods [4], which for years were the best choice for such problems. The main advantage of the finite differences is their scalability at the supercomputers. We hope that this paper can be useful for the further development of the realistic MHD turbulence models with rotation. In spite of the incompressible form of convection’s equations, considered above, the substitute of variables can be used for transition to the anelastic approximation, suitable to the dynamo in the compressible medium.

### Footnotes

- The sign of the boundaries curvature defines the direction of the Rossby waves propagation [7].

### References

- Roberts, P.-H. (1968), On the thermal instability of a rotating-fluid sphere containing heat sources, Phil. Trans. R. Soc, A263, 93–117.
- Arzimovich, L., A., S.Yu. Lukianov (1972), Motion of the charged particles in the electromagnetic fields. In Russian, Nauka.
- Bandaru, V., T. Boeck, D. Krasnov, J. Schumacher (2016), A hybrid finite differenceâboundary element procedure for the simulation of turbulent MHD duct flow at finite magnetic Reynolds number, J. Comp. Phys., 304, 320–339.
- Brandenburg, A., K. Subramanian (2005), Astrophysical magnetic fields and nonlinear dynamo, Phys. Rep., 417, 1–209.
- Buffett, B. (2003), A comparison of subgrid-scale models for large-eddy simulations of convection in the Earth’s core, Geophys. J. Int., 153, 3, 753–765.
- Busse, F.-H. (1970), Thermal instabilities in rapidly rotating systems, Fluid Mech., 44, 441–460.
- Busse, F.H. (2002), Convective flows in rapidly rotating spheres and their dynamo action, Phys. Fluids, 14, 1301–1314.
- Canuto, C., M. Y. Hussini, Q.A. Zang (1988), Spectral Methods in Fluids Dynamics, Springer-Verlag.
- Cattaneo, F., T. Emonet, N. Weiss (2003), On the interaction between convection and magnetic fields, Astrophys. J., 588, 2, 1183–1198.
- Glatzmaier, G., P.H. Roberts (1995), A three-dimensional convective dynamo solution with rotating and finitely conducting inner core and mantle, Phys. Earth Planet. Int., 91, 63–75.
- Jones, C. A. (2000), Convection-driven geodynamo models, Phil. Trans. R. Soc. London, A358, 873–897.
- Jones, C. A., P. H. Roberts (2000), Convection-driven dynamos in a rotating plane layer, J. Fluid Mech., 404, 311–343.
- Pedlosky, J. (2012), Geophysical fluid dynamics, Springer-Verlag.
- Reshetnyak, M. (2017), The anisotropy of hydrodynamical and current helicity, Astronomy Reports, 61, 9, 783–790.
- Reshetnyak, M., V. Pavlov (2016), Evolution of the Dipole Geomagnetic Field. Observations and Models, Geomagnetism and Aeronomy, 56, 1, 110–124.
- Rüdiger, G., R. Hollerbach, L.L. Kitchatinov (2013), Magnetic Processes in Astrophysics: Theory, Simulations, Experiments, Wiley-VCH.