Regular reversals in the six-jet geodynamo model

. A low-mode geodynamo model is developed, controlled by 6-jet convection in the core of the Earth. The model contains only four modes, representing the ﬁelds of temperature, velocity, and two ﬁelds of magnetic induction. The magnetic modes was chosen by combining eight magnetic modes of free decay. There are two noise components in the model. In the model, stable regimes of generation of a magnetic ﬁeld with reversals having a regular character were obtained. These reversals do not cause changes in the convection structure.


Introduction
The process of formation of magnetic fields of planets and stars is successfully explained by the theory of hydromagnetic dynamo [6][7][8]. Dynamo system modeling is an intensively developed section of dynamo theory. Their study is carried out using direct numerical modeling, and on the basis of simplified models. DNS of magnetohydrodynamic equations allows one to reproduce the real magnetic field regimes of planets and stars, but does not answer questions about their causes. In addition, complete equations contain a lot of parameters, the estimates of which often vary by many orders of magnitude, or do not exist at all. Therefore, using simple dynamic systems of small dimension, they are trying to explain the physical cause, the signs, the most important properties of this phenomenon, for example [1][2][3]. We will build a large-scale model in the form of a low-dimensional dynamic system with the help of some modifications of the classical spectral methods [4,9].

Basic equations
Let us consider a spherical layer of an incompressibleuid (liquid core), which rotates about the Oz-axis with an angular velocity Ω. In a spherical coordinate system (r, θ, ϕ) the inner core boundary (ICB) is r = r i and the core-mantle boundary (CMB) is r = r o . The temperature at the ICB and CMB is constant and equal to T i and T o , respectively.
The equations of geodynamo considering the α-effect in the Boussinesq approximation have the form: (1) The governing dimensionless parameters of the model are the Rossby number Ro = u/2ωr o , the Rayleigh number Ra = βg 0 δT r 3 o /νk, the Roberts number q = k/η, the Reinolds number Re = ur o /ν, the magnetic Reinolds number Re m = ur o /η, and the magnitude of α-effect R α . In these expressions, ν denotes the kinematic viscosity, β -the thermal expansion coefficient, g 0 -is the gravitational acceleration at the CMB, δT = T i − T o -is the temperature difference between the ICB and the CMB, k -is the thermal, η -is the magnetic diffusivity end u -characteristic velocity.
This choice of dimensionless parameters implies that r o , r 2 o /η, δT , 2ωρ 0 η and √ 2μμ 0 ρωη, with ρ -being the density and μ -the magnetic permeability, serve as the length scale, time scale, temperature scale, pressure scale and the scale for the magnetic field, respectively. To maintain the correct ratio at this length scale, the radius of the inner coreri r i = 0.35.
In (1) the variable T denotes the temperature deviation from the equilibrium It is assumed that turbulence is isotropic and the α-effect is used for scalar parameterization in the form α = cos θ.
For the velocity v, no-slip boundary conditions are applied. Moreover, it is assumed that the magnetic permeabilities of the outer and inner cores are equal, and the surrounding of the core is electrically non-conducting. Therefore, for the magnetic field there are specified potential boundary conditions at the CMB and finiteness at the center of the Earth. Finally, the temperature deviation is zero at the ICB and CMB.
Then we apply some spectral decompositions for the velocity, temperature and magnetic field. By the Galerkin method, using these decompositions, the equations of the six-cells model will be obtained.

The convective part of the model
As the velocity mode, one of the eigenmodes of available oscillations in the layer of a viscous liquid is used. The vertical component of this mode has no zeros along the radial boundary between the boundaries of the layer, i.e. it ensures a complete transfer of fluid between the boundaries. The temperature mode is also one of the eigenmodes of free temperature oscillations (eigenmodes of the Laplace operator), it is consistent with the spatial structure of the velocity mode.
It is important that, for any single-mode approximation, of veloocity to using the Galerkin procedure, the Coriolis term in the Navier-Stokes equation will vanish. Therefore, it is necessary to use a velocity, spatial structure, which contains information about the rotation of the layer. As such a method one can use one of the following: a mode, one of the solutions of the spectral problem For convenience, we retain the same dimensionless parameters as in the equation (1). Here the magnetic Prandtl number P m = ν/η.
The non-viscosity analogue of this problem has a skew-hermitian symmetry and is known as the Poincare problem. Therefore, we shall henceforth refer to the solutions of the problem (2) as Poincare modes. In the viscous case, the problem operator does not possess Hermitian or skew-Hermitian symmetry. This creates great problems with obtaining its explicit solutions [5].
We shall use the simplest approximations of these solutions with the using of the eigenmodes of free oscillations of a nonrotating shell. The structure of such modes is well known. Formally, they solving the spectral problem The operator of this problem (the Laplace operator) of Hermitian solutions decays into two orthogonal subspaces of toroidal and poloidal modes and form a complete orthogonal system.
The toroidal and poloidal eigenmodes have the form, respectively v P k,n,m = n(n + 1) where R T k,n (r) and R P k,n (r) some combinations of power functions and spherical Bessel functions , and Y m n (θ, ϕ) spherical harmonics. The entire set of solenoidal fields that are zero on the boundary is the direct sum of the linear orthogonal subspaces generated by the following sets of modes Each of these subspaces is invariant under the operator of the problem (2). This is easy to see if we write the matrix of the operator in the basis (4). Therefore, each of the Poincare modes lies entirely in one of the subspaces (5).
To simulate 6-jet convection, the space H P 2 , including the poloidal mode v P 0,4,±2 , is important. The formula (4) shows that the radial component of this mode in the (θ, ϕ) direction is determined by the spherical harmonic Y ±2 4 . In our model, we will use only one mode of velocity v 0 , which is the simplest approximation of the Poincare modes from the subspace H P 2 . We specify v 0 in the form where the coefficients β i are determined from the problem (2) by the Galerkin method.
It can be seen that the radial component v 0 has in (θ, ϕ)-structure according to the structure Y ±2 4 .
Now we consider the eigen modes of temperature free dissipation in an external core. They have the form T k,n,m = X k,n (r)Y m n , where X k,n (r) are linear combinations of spherical Bessel functions of the first and second kind.
In our model, the basic temperature mode T must be consistent with the radial velocity component. Therefore, we define it by the formula where β 3 and β 4 are the same as in (6). Finally, we choose the following velosity and temperature representation for our model

A kinematic dynamo model
For the magnetic field, we use some of the free decay modes, i.e. solutions of the spectral problem where B T and B P are the toroidal and poloidal magnetic field components, respectively.
The toroidal and poloidal solutions of this problem are

Solar-Terrestrial Relations and Physics of Earthquake Precursors
where j n (·) denotes the spherical Bessel functions of the first kind. The eigenvalues η T kn and η P kn are determined from the boundary conditions, and a T kn and a P kn the normalization coefficients. Further we assume these coefficients such that the r.m.s. of the modes is equal to unity. All these modes form an orthogonal system in the core.
To select magnetic modes, we used the following scheme borrowed from [9]. First, we order these modes as ascending eigenvalues (dissipation rate). The following sequence was obtained, where underlining means doublets with the same eigenvalues Next, we selected several modes with the lowest eigenvalues, denoting them B k (r). The magnetic field is written as B (r, t) = k g k (t)B k (r) and substituted into the induction equation (the third equation in (1)). The velocity, at the same, time can be written as (8), where u(t) = 1. Now applying the Galerkin method yields the system where η i is an eigenvalue of B i . In fact, we have obtained a small-mode approximation in kinematic dynamo with six convective cells. Let us denote by λ i the eigenvalues of the matrix of the system (12). The dynamo "works" if max λ i > 0. The mode with its eigenvalue is the leader, it grows faster than others. It is important that if the eigenvalue of the leading mode is complex, the magnetic field will oscillate.
Thus, we leave eight magnetic modes out of (11) that perform a working dynamo: P 0 01 , P ±2 03 , P 0 11 , T ±2 04 , P ±2 05 . The scheme of generation of magnetic modes of large-scale convection and the α-effect is shown in fig.4. To get acquainted with this particular choice of modes it is possible in [4].

A combined model of magneto-convection
In the above-described kinematic dynamo, only an unlimited growth of the field is possible. To obtain a stable generation of a limited field, a suppression mechanism must be introduced. Physically, this is realized due to the influence of the field through the Lorentz force on the structure of the flows. We also introduce into the system the mechanism of algebraic suppression of the α-effect.
Recall that in our construction Re m is a fixed value of the amplitude u(t) of the velocity mode v 0 . Therefore, we further consider a model of magnetohydrodynamic convection with variable velocity. To do this, the above-described decomposition of velocity, magnetic field and temperature must be substituted into all equations (1) and Galerkin method applied.
We introduce the noise components into the model in the α -generator and in the term before the Rayleigh number. We obtain a nonlinear dynamical system describing  self-consistent 6-jet magneto-convection in the core of the Earth

Solar-Terrestrial Relations and Physics of Earthquake Precursors
where μ and ς -eigenvalues of the velocity mode v and the temperature mode T , respectively, and S, L ij , H, W ij , A ij -are Galerkin coefficients. Noises ξ(t) and ζ(t) is a stochastic process with zero mean.

Simulation results
A numerical study was conducted in a model with eight modes. When performing a numerical simulation with the model (13), we used molecular values of the dissipation coefficients k = 10 −5 m 2 /s, η = 1m 2 /s, following the estimations made in [13]. The outer radius of the core r o = 3480 km and the angular velocity Ω = 7.29 × 10 −5 rad/s. The relevant parameter values are Re = 10 8 , P m = 10 5 and Ro = 3 · 10 −6 . Thus, there remain three control parameters: the Rayleigh number, the magnetic Reynolds number, and the amplitude of the α effect, depending on the values of these parameters, different dynamo modes are obtained. Note that by fixing the Reynolds number and varying the magnetic number of Reynolds, we actually vary the kinematic viscosity ν, which is the most uncertain parameter in the dynamo system. Estimates of the kinematic viscosity for the Earth according to estimates of various studies vary by many orders of magnitude.  0.549e-10 0.942e-10 0.501e-5 0.308e-10 0.371e-9 0.131e-8 0.304e-7

Solar-Terrestrial Relations and Physics of Earthquake Precursors
First, calculations of the solutions of the (13) system were carried out for identically zero. The model used four Rayleigh numbers Ra = 10 10 , Ra = 6.5 · 10 10 , Ra = 10 11 and Ra = 5 · 10 11 . Different R α and Re m numbers values were 0.1 ÷ 100.
First of all, an experiment was conducted for several combinations of selected parameters, the purpose of which was to determine the real dimension of the magnetic subspace of the model. The points of the calculated trajectories in the 8-dimensional space of magnetic modes were considered as an ellipsoid of scattering. It turned out that one, maximum two semi-axes are much more than the others. Those. one can say that the magnetic subspace is no more than two-dimensional. And the coordinates of the vectors of these largest axes determine the combined magnetic modes. In practice, these calculations were performed according to the standard scheme for analyzing the main components [12]. Examples of calculations are given in the table 1.
The analysis of the directions of the largest semiaxes showed that the initial magnetic modes are mainly divided into two groups corresponding to these directions, see the table 2. Therefore, it can be said that the model implements the generation of a magnetic field by mutual pumping of two magnetic structures corresponding to two combined modes. Physically, in axisymmetric works on a dynamo, this is expressed by the mutual pumping of toroidal and poloidal fields. In this paper, this separation does not work, but the general principle of the two magnetic structures that exchange energy continues to be preserved. Thus, it is proposed to make a two-mode combination of eight modes.  Finally, we choose the following representation of the magnetic field for our model B 1 = k 1 P 0 01 + k 2 P 2 03 + k 3 P −2 03 + k 4 P 0 11 + k 5 T 2 04 + k 6 T −2 04 + k 7 P 2 05 + k 8 P −2 05 B 2 = m 1 P 0 01 + m 2 P 2 03 + m 3 P −2 03 + m 4 P 0 11 + m 5 T 2 04 + m 6 T −2 04 + m 7 P 2 05 + m 8 P −2 05 (14) Thus, the system (13) acquires a view

Solar-Terrestrial Relations and Physics of Earthquake Precursors
Two noise processes ξ(t) and ζ(t) have been added to the system (15). These processes are the result of synchronization of the higher discarded velocity and magnetic field modes. Noise excitation ζ(t), which simulates the fluctuations of the heat flux from the inner core. Noise excitation ξ(t) models the spontaneous synchronization of small-scale modes that are not used in the model [14]. This spontaneous synchronization effect is a well-known phenomenon in the theory of turbulence. The structure of processes is as follows. We define on the time axis a random sequence of points 0 < τ 1 < θ 1 < τ 2 < θ 2 < . . . < τ k < θ k < . . . We assume that the of the coherent structure number k is formed at the time τ k and is destroyed at the time θ k . Then T est k = τ k − θ k−1 is the latensy time for the formation of the next structure, and T k = θ k − τ k is the time of its existence. During the waiting time, the process ξ(t) is zero, and for the time of existence ξ(t) = ξ k . Here ξ k are independent random variables with zero mean. The laws of distribution of these quantities, as well as T est k and T k , are chosen for numerical modeling.
The value of the Rayleigh number Ra = 10 10 does not give us stable solutions, it is obvious that it is not enough to launch the dynamo. The value of the Rayleigh number Ra = 5 · 10 11 immediately goes to infinity, it turns out to be too large. The two average values of the Rayleigh number Ra = 6.5 · 10 10 and Ra = 10 11 give an approximately identical picture of the solutions. We demonstrate the results of the study of this model. On the right are detailed fragments of calculations.
In fig.2. it can be seen that the velocity retains its sign, and in a magnetic field quasi-periodic reversals are observed with individual bursts. In fig.3. system goes out stationary. In fig.4. we see regular fluctuations in speed around a single level. And the magnetic field works in the mode of bursts. Situations with a change in the sign of speed are not physically plausible because assumes a reversal of convection in the opposite direction, but if we average the velocity over the period (i.e. go to a larger time scale), the reversals in the field are preserved, and the speed will become constant. In fig.5. we see the same velocity oscillations as in the previous case, but there are no reversals for the magnetic field, but the solution is in the nature of bursts.
In this way, as a result of the calculations, we found that for small values of the control parameter Re m < 40, we have mainly quasi-periodic and stationary solutions. With a large value of Re m > 50 and there are regular oscillations of the magnetic field and dynamo-burst with reversals. This is reflected in the time graphs of the amplitudes of the velocity and the dipole part of the field.   In the simulation, we used the exponential distribution law for waiting times and existence, and these quantities themselves were independent. Average values of T est k = 5 and T k = 30, i.e. the characteristic time of existence of coherent structures is much less than their waiting time. The jumps of ζ k and ξ k were uniformly distributed over the interval [-1;1].
The same values of the control parameters were used, with the imposition of noise. This is reflected in the time graphs of the amplitudes of the velocity and the dipole part of the field ( fig.6.)-(fig.9.).
By introducing noise components into the system, we expected to transfer the system from one parameter space to another, expecting to switch dynamos from one mode to another. But such a switch, as can be seen from the figures does not occur, the reason may be that we push the system too often. In the development of the  model, if the average waiting time and the existence of these structures were large, then it should already be infinitely large. As an option, the power of waiting time.

Conclusions
In this paper, a low-mode dynamo model is developed. The model contains, it seems to us, the minimum possible number of modes (four). For the fields of speed, temperature, and magnetic induction, the representation is used, using modes combined from the eigenmodes of the Laplace operator. For speed, this combination approximates one of the modes of the Poincare operator. The temperature mode is combined to match the speed. The two magnetic modes are combined on the basis of a computational experiment in a model containing eight modes of free damping of the magnetic field. For this, the ellipsoid of scattering was calculated. By the magnitude of the   eigenvalues, it was clear that sometimes there are signs of a single-mode approximation in the system, but basically it is two-mode. It is known that in the axisymmetric case of a dynamo this is expressed by the mutual pumping of toroidal and poloidal fields. In this paper, this separation does not work, but the general principle of the two magnetic structures that exchange energy continues to be preserved.

Solar-Terrestrial Relations and Physics of Earthquake Precursors
Numerical simulation shows that the model can be implemented in a stable mode of generating a magnetic field with and without reversals. The field reversals in the model have a regular character, or the character of dynamo-bursts. A wide variety of modes can also be provided by different distributions of the waiting times for pulses in noise. For simplicity, we chose an exponential distribution for modeling. It seems that the introduction of other types of distributions, for example, a power distribution, will give more complex statistics of reversals.