Scaling relations in quasi-static magnetoconvection with a strong vertical magnetic field

Abstract The scaling law for the horizontal length scale $\ell$ relative to the domain height $L$, originating from the linear theory of quasi-static magnetoconvection, $\ell /L \sim Q^{-1/6}$, has been verified through two-dimensional (2-D) direct numerical simulation (DNS), particularly at high values of the Chandrasekhar number ($Q$). This relationship remains valid within a specific flow regime characterized by columnar structures aligned with the magnetic field. Expanding upon the $Q$-dependence of the horizontal length scale, we have derived scaling laws for the Nusselt number ($Nu$) and the Reynolds number ($Re$) as functions of the driving forces (Rayleigh number ($Ra$) and $Q$) in quasi-static magnetoconvection influenced by a strong magnetic field. These scaling relations, $Nu \sim Ra/Q$ and $Re \sim Ra Q^{-5/6}$, have been successfully validated using 2-D DNS data spanning a wide range of five decades in $Q$, ranging from $10^5$ to $10^9$. The successful validation of the relations at large $Q$ values, combined with our theoretical analysis of dissipation rates and the incorporation of the horizontal length scale's influence on scaling behaviour, presents a valid approach for deriving scaling laws under various conditions.


Introduction
The convective motion of electrically conducting fluids in the presence of an externally imposed magnetic field plays a crucial role in both natural phenomena and practical flows.In astrophysical environments, fluid motions are intricately coupled with magnetic fields, such as in the outer layers of the Sun and other late-type stars.The interaction between thermal convection and magnetic fields in these scenarios has been extensively studied [Proctor and Weiss, 1982, Cattaneo et al., 2003, Schüssler and Vögler, 2006].In industrial systems, magnetic damping is commonly employed in metallurgical applications, while magnetic levitation, pumping, and heating of liquid metal are vital aspects of nuclear engineering [Davidson, 2001].Despite its prevalence, understanding the dynamics of magnetoconvection in extreme parameter regimes, particularly in geophysical and astrophysical systems, remains a formidable challenge.Therefore, it is crucial to identify relevant flow regimes and establish corresponding scaling relations to characterize these complex flows.By discerning the points at which flow transitions occur, we can carefully extrapolate scaling laws to a reasonable range of input parameters, thereby providing practical utility in understanding and predicting these phenomena.
The planar Rayleigh-Bénard convection (RBC) is a useful model for investigating the fundamental dynamics of magnetic fields on convection.The canonical RBC consists of a fluid layer confined between a pair of plane, parallel boundaries separated by a distance L.
The boundaries are maintained at a constant temperature difference ∆.In the absence of a magnetic field, the Rayleigh number Ra = gβ∆L 3 /νκ and the Prandtl number P r = ν/κ characterize the thermal forcing and the ratio of the diffusivities present in the system, where β is the coefficient of thermal expansion, g the acceleration due to gravity, and ν and κ being the momentum and thermal diffusivity respectively.In the presence of an external magnetic field of strength B 0 , the damping effect by virtue of the Lorentz force must be considered.
The strength of the Lorentz force to the viscosity is quantified by the Chandrasekhar number , where ρ is the density, µ the magnetic permeability, and η magnetic diffusivity.
Early studies on magnetoconvection, both experimental as well as theoretical [Thompson, 1951, Chandrasekhar, 1961, Busse and Clever, 1982, Nakagawa, 1957, Jirlow, 1956], are primarily devoted to the linear stability analysis.In recent decades, much focus has been on the scaling relations between the dimensionless heat transport as a function of Ra and Q.The dimensionless convective heat transport is characterized by the Nusselt number N u, the ratio of the total heat flux to the conduction heat transfer.In recent decades, experiments by [Cioni et al., 2000, Aurnou and Olson, 2001, Burr and Müller, 2001, King and Aurnou, 2013] have suggested different scalings.Aurnou and Olson [2001] found N u ∼ (Ra/Q) 1/2 in contrast with the findings of Burr and Müller [2001] who report a scaling law N u ∼ (Ra/Q) 2/3 .Both of these studies consider 0 < Q ≲ 10 4 .Yu et al. [2018] conducted numerical simulations over a similar range of values of Q, and found the scaling relations between N u and Ra/Q are in good agreement with the empirical formulas obtained in the experimental studies of Aurnou and Olson [2001].Experimental investigations by Cioni et al. [2000] pushed the range of Chandrasekhar number by two decades up to Q ≈ 10 6 , and comprehensively characterized two major regimes influenced by the magnetic effects, with an additional transitional regime between them.The first regime, interpreted by them as a condition of marginal stability for the thermal boundary layer, the heat transport scaling N u ∼ Ra/Q holds.For the second major regime, which is characterized by higher thermal forcing, they report N u ∼ Ra 0.43 Q −0.25 .Recently, direct numerical simulations of quasi-static magnetoconvection by Yan et al. [2019] with magnetic strengths up to Q = 10 8 have shown the scaling laws of Aurnou and Olson [2001], Burr and Müller [2001] and Yu et al. [2018] to be limited to low values of Q.The proposed convective regimes are distinguished by flow characteristics, with the first regime being reminiscent of the linear convection, characterized by laminar, cellular structures.The heat transport scaling law for the second regime, characterized by the existence of quasi-laminar columnar structures, is suggested as N u ∼ (Ra/Q) γ , where Based on an exponential fit for γ, a value of 0.95 is suggested at Q = 10 16 .However, as the authors mention, it must be noted that the exponential fit has no physical basis and is shown only to provide a guide for the behavior at large values of Q.
Based on the scaling approach introduced by Grossmann and Lohse [2000], Zürner et al.
[2016] decomposed the dissipation rates in magnetoconvection into their bulk and boundary layer counterparts and identified four distinct regimes distinguished by the strength of the external magnetic field and the level of turbulence in the flow.Studies by Zürner et al. [2020], Akhmedagaev et al. [2020] report a detailed analysis of the spatial structure of magnetoconvective flows in cylindrical geometry with no-slip walls.Wall modes similar to those in RBC with rotation are observed near the linear stability limit.Furthermore, in these studies, the scaling relations for normalized Nusselt and Reynolds numbers reveal that the global transport properties approach a universal power law at larger degrees of supercriticality.
Following Grossmann and Lohse [2000], Bhattacharyya [2006] obtained N u ∼ N u(Ra, Q) scaling relations for different regimes characterized by weak and strong magnetic fields.It must be noted here Bhattacharyya [2006], Zürner et al. [2016] assumed the length scale associated with the magnetic induction to be the domain height L. The presence of a magnetic field opposes fluid motions perpendicular to the field lines while leaving the parallel component unopposed, thus the classical RBC under the effect of an externally applied magnetic field becomes highly anisotropic.Due to anisotropy, the length scale of the flow is significantly altered; the scaling of which remains vital in the scaling laws associated with the heat and momentum transport and Ohmic dissipation.The primary objective of our study is to explore the impact of the horizontal length scale on the scaling relations previously discussed, focusing specifically on the framework of 2D quasi-static assumptions.Although the real-world applications of magnetoconvection are three-dimensional (3D), 2D simulations are substantially cheaper in terms of computational time and thus can be utilized to understand certain aspects of the flow, including scalings.Goluskin et al. [2014] and Wang et al. [2020] are some noteworthy 2D studies that have explored the ability of the convection to drive vertically sheared, large-scale horizontal flow.Building on these works, recently, quasi-2D simulations of magnetoconvection have been used to study the magnetic damping of jet flows which finds applications in probing the mechanisms relevant to damping of large scale azimuthally directed jets on Jupiter [Aggarwal et al., 2022].Here, our 2D investigation extends beyond previous studies by considering even larger values of Q, reaching up to 10 9 .This range of values is particularly relevant for geophysical and astrophysical systems, such as the Earth's outer core, where Q is estimated to be on the order of 10 16 [Gillet et al., 2010].By identifying scaling relations that converge at values of Q approaching realistic estimates, we aim to enhance our understanding of the fundamental aspects of convection in these systems.

Physical model and numerical settings 2.1 Governing equations
The dimensionless equations governing the flow of a conducting fluid in Boussinesq approximation, under the action of a quasi-static, vertically imposed external magnetic field where u, p and θ represent the velocity, pressure and temperature fields.The vertical unit vector anti-parallel to the direction of gravity is e x , and Φ and J represent the electric potential and the current density respectively.At all places in the paper, x and y represent the vertical and horizontal directions respectively.The input parameters Ra, P r, and Q are defined in section 1.The rest of the symbols have their usual meanings.In equation 1, domain height L, free-fall velocity u f f = (βg∆L) 1/2 , free-fall time scale t f f = L/u f f , imposed temperature difference ∆, free-fall velocity-based dynamic pressure ρ 0 u 2 f f , and the magnitude of the imposed magnetic field B 0 are adopted for non-dimensionalization.All lateral boundaries are set to periodic.At the top and bottom walls, stress-free boundary conditions for u, isothermal for θ, and an insulating boundary condition on Φ are imposed, u x = ∂u y /∂x = 0 at x = 0, 1; θ = 1(0) at x = 0(1); ∂Φ/∂x = 0 at x = 0, 1.To compute J, the constraint of charge conservation, given by equation 1e, is invoked to obtain a Poisson equation for Φ, The Reynolds number, Re, referred to in the subsequent sections is defined as, Re = UL/ν, where U is the characteristic velocity.Depending on if the characteristic velocity , we further define Re, Re x and Re y to characterize the total, vertical and horizontal Reynolds numbers respectively.⟨•⟩ A and (•) denote area-average over the whole 2D domain, and time-averaging respectively.

Numerics
A conservative, second-order centered spatial discretization is employed on a staggered grid [van der Poel et al., 2015].Time marching is performed with a fractional-step third-order Runge-Kutta (RK3) scheme, in combination with the Crank-Nicholson scheme for the implicit terms.In the following section, the implementation of the subroutines for Φ and the validation of the code are described in detail.

Code validation
In this work, we have used AFiD [van der Poel et al., 2015, Zhu et al., 2018]  The solver was validated by comparing the Nusselt number against the values reported in Yan et al. [2019] at diffferent Rayleigh numbers away from the onset.Table 1 shows the Nusselt number from 3D test runs at Q = 10 6 to reproduce results from the literature.The resolutions employed for these cases in order, from smallest to largest Ra, are: (96×180×180), respectively.
In addition to that, the dominant length scale at Ra ⪆ Ra c (very close to the onset) was also calculated from the energy spectra of the vertical velocity, and validated against the predictions from the linear theory.In figure 1, the normalized kinetic energy spectra of the ] computed at mid-depth are shown.u x and u * x are the Fourier transform of u x and its complex conjugate respectively.The dominant integer wavenumber represented by the highest peak of spectra, indicates the number of (near-onset) critical horizontal wavelengths (λ c ) present in the domain.The aspect ratio for Q = 10 3 , 10 5 and Q = 10 7 cases was set to Γ = 11.0,Γ = 5.0 and Γ = 2.25 respectively which admit ≈ 10λ c .For Q = 10 6 case, Γ = 3.5, which allows for ≈ 11λ c to be resolved.The number of critical horizontal wavelengths observed in our simulations match reasonably well with the linear theory predictions for a given aspect ratio.3 Results and discussion

Parameter range and flow structure
We have considered five sets of simulations with Q ranging from 10 5 − 10 9 .For each set, at a fixed Q, the Rayleigh number is varied up to Ra ≳ 100Ra c , Ra c being the critical Rayleigh number at a given Q, which is calculated from the expression of the Rayleigh number corresponding to the marginal stability of the horizontal wavenumber.For details, the interested reader is referred to Chandrasekhar [1961].In all of our simulations, we have used P r = 1.More details on grid resolution and the parameter ranges explored are listed in Appendix A.
Figure 2 shows the visualization of the effect of increasing Ra at a fixed Q (Q = 10 6 ).
The flow structures at all Ra are stabilized by the Lorentz force, resulting in the convection directed along the externally imposed magnetic field.For the cases in the first row from the bottom, flow is dominated by cellular structures which eventually morph into columns as Ra is increased.The transport of heat in the cellular regime is primarily via conduction.In the columnar regime, the length scale of the structures stays more or less the same, until a third regime of magneto-plumes emerges.Upon further increasing the Rayleigh number, the plumes begin to merge, characterized by stronger buoyancy relative to the Lorentz force.We will not be dealing with the scaling laws in this merging-plume regime.With increasing Ra, the thermal boundary layer thickness also decreases.This has important consequences on the heat transfer scaling within the columnar regime, which we will discuss later.The flow in cellular, columnar and plume regimes is strictly steady.Unsteadiness kicks in when the stabilizing effect of the Lorentz force is weakened at around Ra/Ra c ∼ 80, which also marks the onset of the merging plume regime.Figure 2r shows such a case.Upon comparisons with the flow structures in 3D simulations, we observe that the inter-column distance in 2D is identical in the horizontal direction.However, no such uniformity is seen in 3D simulations of Yan et al. [2019].To illustrate the differences, in Figure 3, we have presented instantaneous temperature contours computed in a 3D simulation in the columnar regime.Given insulating boundary conditions (∂ x Φ| x=0,1 = 0) (non-conducting walls), electric currents induced by the horizontal motion (u × e x ) will close within the convective roll itself.In 3D, this gives rise to currents in the horizontal plane (y − z) near the boundaries.In the case of 2D, since one horizontal component of velocity is zero (u z in the present case), the current associated with that component also ceases to exist, and only one of the components, u y , remains.This gives rise to the Lorentz force which is confined to the only horizontal degree of freedom (y) existing in the 2D flow.Since everything is constrained to a single horizontal dimension, there is a tendency of the flow structures to adhere to a specific equidistant pattern, until the stabilizing effect of the Lorentz force is weakened at relatively higher Ra.At that point, the interplume distance starts to vary.
The anisotropy in the convective flow can be characterized by examining the Reynolds number based on the vertical and horizontal components of the rms velocity.Figure 4e clearly shows that the vertical component of the velocity dominates over the horizontal component.
As Ra is increased, the vertical component increases significantly relative to the horizontal component.
In figures 4a and 4b, vertical profiles of time-and horizontally-averaged mean temperature ⟨θ⟩ y , and horizontally-averaged rms fluctuations of temperature ⟨θ rms ⟩ y at Q = 10 6 are presented.For lower values of Ra, the rms fluctuations of temperature are small, and the mean temperature follows the conduction profile.At larger values of Ra, the mean temperature tends to be an isothermal profile in the bulk, and the near-wall peaks in the rms profiles indicate the presence of well-developed thermal boundary layers.The darkest shade in figure 4 represents Ra = 7 × 10 7 , which lies in the columnar regime (see figure 2k).
Figures 4c and 4d show the depth-dependence of the time-and horizontally-averaged viscous (ϵ ν = (P r/Ra) 1/2 ⟨|∇u| 2 ⟩ y ) and Ohmic dissipation (ϵ η = Q(P r/Ra) 1/2 ⟨|J| 2 ⟩ y ).Since we have adopted stress-free boundary conditions, the viscous dissipation is larger in the bulk of the domain, away from the boundaries.Ohmic dissipation, on the other hand, is dominant near the boundaries and decreases in the interior of the domain.At larger Ra, an increase in the Ohmic dissipation is expected due to larger currents generated by enhanced convection.1.94, 2.43, 2.92, 3.89, 4.86, 5.84, 6.81, 8.75, 9.73, 14.59, 19.45, 24.32, 38.91, 77.82, 97.28.The colormap values have been adjusted to [0.2, 0.8] for better visibility.N u (a) simulations.This is because, for 2D simulations, consistent with the mass conservation, we observed a higher Re x compared to that in 3D at the same Ra.To see a comparison of Re x in 2D and 3D, refer to figure 4e.Higher Re x is responsible for enhancing the buoyant fluxes in 2D, which ultimately translates into an increased N u, at least for the regime that we are most interested in (the 2D columnar regime), e.g.see figure 5a.The differing transition locations are evident.In addition to that, we can also clearly observe the columnar regime in both 2D and 3D, albeit with slightly different N u.The transition locations and the magnitude of the heat flux and flow velocity may differ; however structural similarities exist, and for 2D and 3D alike, the flow is characterized by cellular, columnar, and plume-like structures for the considered Ra range at different values of Q.This is important for scaling since a specific regime in 3D is expected and shown to have a counterpart in 2D.

Scaling analysis
The existence of the columnar and plume regimes can be explained very well by invoking the predictions from marginally stable thermal boundary layer analysis, which predicts and N u ∼ L/λ θ ∼ Ra 1/3 at smaller values of Q λ θ due to thinner thermal boundary layer thickness in the plume regime [Bhattacharjee et al., 1991].See figure 5b.
In the following sections, for the columnar regime, we will explain the scaling laws concerning not only the heat transport, but also the horizontal length scale, flow velocity, and Ohmic dissipation by using energetic arguments.To begin, we start with the exact relations for time-and volume-averaged (area-averaged in 2D) viscous dissipation rate ϵ ν = ν⟨|∇u| 2 ⟩, Ohmic dissipation rate ϵ η = (ρσ) −1 ⟨|J| 2 ⟩, and the thermal dissipation rate ϵ θ = κ⟨|∇θ| 2 ⟩ [Grossmann and Lohse, 2000, Ahlers et al., 2009, Bhattacharyya, 2006, Song et al., 2023], The convective part of ϵ θ can be expressed as, Introducing u, θ and ℓ as the representative scales for convective velocity, temperature, and the length scale, the convective part of thermal dissipation rate can be expressed as εθ ∼ uθ 2 /ℓ.Using this, after rearranging, we can rewrite equation 4 as [Song et al., 2023], N u can also be expressed in terms of convective heat transfer q ∼ uθ, as, From equations 5 and 6, we get ℓ/L ∼ θ/∆, leading to, Assuming the induction b, the total field can be written as B = B 0 + b.From Ohmic dissipation ϵ η = (ρσ) −1 ⟨|J| 2 ⟩, and J = (∇ × B)/µ, and assuming the length scale is ℓ we can write, where b represents the typical magnitude of the induced field.Under the quasi-static assumption, the induction equation can be simplified to, which upon using the scale estimates for b ∼ b, u Using this in equation 8, after rearranging and grouping relevant terms, we finally obtain a scaling law for the Ohmic dissipation with the horizontal length scale dependence, With Re and ℓ dependence of N u and ϵ η , we will proceed to derive the scaling relation of Re ∼ Re(Ra, Q, P r, ℓ/L) and ℓ/L ∼ Q γ in the columnar regime.It is to be noted the length scale of the structures in the columnar regime varies negligibly with Ra, which is why its Ra dependency is not considered.In figure 7a, in the columnar regime, we observe ϵ η >> ϵ ν .
Hence from the equation 2, we can consider ϵ η ≈ (ν 3 /L 4 )(N u − 1)RaP r −2 .Combining this together with equation 11 and 7, we derive the scaling law for Re in terms of the input parameters and the horizontal length scale, Using this in equation 7, length scale and Prandtl number dependence drops out, yielding the Nusselt number scaling, This result is consistent with the scaling law obtained from the marginally stable thermal boundary layer analysis discussed above.For non-magnetic convection, a regime in which the entire fluid layer is turbulent (at asymptotically large values of Ra) is referred to as the ultimate regime.This regime, predicted by Kraichnan [1962], is characterized by the heat transport becoming independent of the diffusion coefficients.With additional constraining forces, like the Lorentz force in the present case, a heat transport regime independent of the diffusion coefficients (ν and κ) can be observed at relatively low Rayleigh numbers.It is noteworthy to point out that in the columnar regime (N u−1 ∼ Ra/Q), an inherently different flow state in comparison to the ultimate state of conventional RBC, the heat transport N u × κ∆/L in the limit of an asymptotically large Q, is independent of ν and κ.It does, however, depend on the magnetic diffusion coefficient η.It is also evident that the scaling of the length scale has no bearing upon the scaling law for the dimensionless convective heat transport in the columnar regime.Equation 13 is validated against our data in figure 8a.
One way to obtain the dependence of the horizontal length scale on the Chandrasekhar number is to obtain the scaling for the critical horizontal wavenumber k c .From linear theory, it is well known that at large Chandrasekhar, 1961].
In figure 7b, the Q-dependence of the length scale found in the data is compared to this result from the linear theory.We find a good collapse of data for all the simulation sets according to ℓ/L ∼ Q −1/6 in the columnar regime.Using this in equations 11 and 12, we obtain the following relations for the Reynolds number and Ohmic dissipation scaling, Re ∼ RaQ −5/6 P r −1 (14) Data from our simulations show a good collapse for the above relations.See figure 8.
In figure 8a, at large Q, the N u data collapse in the columnar regime.Similarly, in figure Ra/Q Ra/(Q 5/6 P r)  8b, the Re data at large Q collapse, validating equation 14.It must however be noted for Reynolds number scaling, the exponent differs from one suggested by Yan et al. [2019] for 3D quasi-static magnetoconvection.They derive the scaling relation for Re by balancing the Lorentz force with the buoyancy in the columnar regime, which is evident from their 3D data.
Our derivations are based on energetic arguments involving the dissipation rates.Different assumptions invoked in these arguments is expected to lead to the differences observed in the exponent.The Ohmic dissipation data, shown in figure 8c, at all Q show a solid trend towards the scaling given by equation 15, collapsing for Q in the asymptotic limit, Q → ∞, in the columnar regime.

Concluding remarks
The effect of the horizontal length scale on the scaling behavior of the response parameters in quasi-static magnetoconvection has not been explored in the literature.Previous studies concerning the energetic scaling arguments have considered the domain height L as the sole length scale present in the system.In the present study, under the quasi-static assumption, we consider the horizontal width of the columns as a relevant length scale and explore the effect of its dependence on the Chandrasekhar number on the scaling of the dimensionless heat transport, flow velocity, and Ohmic dissipation.Based on our DNS data, the derived scaling laws are validated successfully for Nusselt number, Reynolds number, and Ohmic dissipation.
The Nusselt number scaling obtained from our analysis is shown to be consistent with the marginal stability predictions, which further validates our theoretical analysis.The scaling law for the Reynolds number and Ohmic dissipation also shows a reasonable collapse in our data at large Q, implying the validity of our scaling arguments.
In the present work, we have considered periodicity in the horizontal direction and a fixed Prandtl number of unity.Compared to periodic boundaries, the presence of sidewalls has been shown to lead to convection at a much lower Ra [Busse, 2008].Numerical simulations of planar RBC by Liu et al. [2018] confirmed the existence of wall-modes at Ra ≲ Ra c .Akhmedagaev et al. [2020] also confirmed the existence of these wall-modes at Ra ≳ Ra c .
The effect of sidewalls at even higher Q and Ra away from the onset will be something that we will explore in the future editions of this work.The presence of sidewalls may have implications in the design of liquid metal batteries [Kelley and Weier, 2018] and nuclear fusion technology [Ihli et al., 2008], wherein the geometries are often confined by walls.
Magnetoconvection is usually studied in liquid metals, which have very small Prandtl numbers and very high electrical conductivities.In the numerical and experimental low P r studies of Akhmedagaev et al. [2020] and Zürner et al. [2020], unlike the slender columns and magneto-plumes observed in our simulations, the dominant structures were observed to be a system of ascending and descending planar jets originating at the sidewalls.It must be noted, however, that the aforementioned studies consider convection in a cylindrical container in which sidewalls play a dominant role in organizing the flow.Hence, it remains to be seen if the slender columnar or plume structures can emerge in 2D simulations by virtue of a low Prandtl number or exclusively in the presence of sidewalls.Whether or not the flow structures are different, we believe that it does not affect the scaling argument for the conditions that we have considered.That can be seen from the low P r (P r = 0.025) simulations of Yan et al.
[2019], who performed the comparative calculations at P r = 1 and P r = 0.025.In their work, it was observed that the scaling behavior of the dimensionless heat transport seems to follow the same trend as in the case of P r being unity, i.
as the base solver, and implemented the subroutines for solving a Poisson equation for the electric potential Φ, and an explicit Lorentz body force term.Temperature and velocities are computed on cell faces, with θ and u x located at the same face (to avoid interpolation errors in calculating the term θ e x ).The scalar variables p and Φ are computed at the cell centers.Using the cell-centered values of Φ, the components of the current density are then constructed on the cell faces, which are used to finally compute the components of the Lorentz force.

Figure
Figure 5a gives a brief overview of the parameter space characterized by dimensionless heat transport.Three distinct scaling regimes with different slopes are observed at large values of Q: (1) a cellular regime close to the onset, (2) N u ∼ Ra regime, characterized by the existence of elongated, vertically aligned columns, and (3) N u ∼ Ra 1/3 regime, which is characterized by the existence of plume-like structures with thin stems and broadened heads.Some select cases representing these regimes are shown in figure 6.The corresponding scaling regimes in (N u, Ra) space are represented by black dashed (columnar) and dotted (plume) lines in figure 5a.It is important to note here that the states where transitions occur from

Figure 8 :
Figure 8: Scaling of the (a) Nusselt number, (b) Reynolds number, and (c) Ohmic dissipation in the columnar regime.

Table 1 :
Nusselt number from 3D test runs at Q = 10 6 to reproduce results from the literature.

Table 2 :
e. N u ∼ Ra/Q at asymptotically large values of Q; since in the columnar regime, there is no P r dependence of N u.Although, Re and ϵ η show P r dependence in this regime, in this work, our main goal was to test the scaling argument that we proposed for aymptotically large Q.Re and ϵ η scaling relations show a good agreement with what we observed in our DNS data (for Q and Ra dependence).Exploration of the P r dependence of the response parameters and the flow structures is a full-fledged campaign in itself that we propose to study in future.Simulation details.Q is the Chandrasekhar number, Ra is the Rayleigh number, N u m is the mean of the Nusselt numbers computed using the thermal dissipation, kinetic energy dissipation via viscous and Ohmic counterparts, plate Nusselt number, and the Nusselt number based on volume-averaged heat flux.Error represents the maximum deviation (in percentage) between the smallest and the largest individual Nusselt numbers.Re is the Reynolds number, Γ is the domain aspect ratio, N x and N y represent the resolution of the grid in vertical and horizontal directions.