Algorithmic Aspects of Simulation of Magnetic Field Generation by Thermal Convection in a Plane Layer of Fluid

We present an algorithm for numerical solution of the equations of magnetohydrodynamics describing the convective dynamo in a plane horizontal layer rotating about an arbitrary axis under geophysically sound boundary conditions. While in many respects we pursue the general approach that was followed by other authors, our main focus is on the accuracy of simulations, especially in the small scales. We employ the Galerkin method. We use products of linear combinations (each involving two to five terms) of Chebyshev polynomials in the vertical Cartesian space variable and Fourier harmonics in the horizontal variables for space discretisation of the toroidal and poloidal potentials of the flow (satisfying the no-slip conditions on the horizontal boundaries) and magnetic field (for which the boundary conditions mimick the presence of a dielectric over the fluid layer and an electrically conducting bottom boundary), and of the deviation of temperature from the steady-state linear profile. For the chosen coefficients in the linear combinations, the products satisfy the respective boundary conditions and constitute non-orthogonal bases in the weighted Lebesgue space. Determining coefficients in the expansion of a given function in such a basis (for instance, for computing the time derivatives of these coefficients) requires solving linear systems of equations for band matrices. Several algorithms for determining the coefficients, which are exploiting algebraically precise relations, have been developed, and their efficiency and accuracy have been numerically investigated for exponentially decaying solutions (encountered when simulating convective regimes which are spatially resolved sufficiently well). For the boundary conditions satisfied by the toroidal component of the flow, our algorithm outperforms the shuttle method, but the latter proves superior when solving the problem for the conditions characterising the poloidal component. While the conditions for the magnetic field on the horizontal boundaries are quite specific, our algorithms for the no-slip boundary conditions are general-purpose and can be applied for treating other boundary-value problems in which the zero value must be admitted on the boundary.