Mean-field analysis on large-scale magnetic fields at high Reynolds numbers

Solar magnetic fields comprise an 11-year activity cycle, represented by the number of sunspots. The maintenance of such a solar magnetic field can be attributed to fluid motion in the convection zone, i.e. a dynamo. This study conducts the mean-field analyses of the global solar dynamo simulation presented by Hotta et al. (2016). Although the study succeeds in producing coherent large-scale magnetic fields at high Reynolds numbers, the detailed physics of the maintenance of this field have not been fully understood. This study extracts the alpha-tensor and the turbulent magnetic diffusivity tensor through mean-field analyses. The turbulent magnetic diffusivity exhibits a significant decrease towards high Reynolds numbers. The decrease in the turbulent magnetic diffusivity suppresses the energy conversion of large-scale field to small-scale field. This implies that the decrease in the turbulent magnetic diffusivity contributes to the maintenance of a large-scale magnetic field at high Reynolds numbers. A significant downward turbulent pumping is observed; it is enhanced in the weak phase of the large-scale field. This study proposes a cyclic reversal process of a large-scale field which is dominantly driven by the alpha-effect and is possibly triggered by downward pumping.


INTRODUCTION
Solar magnetic fields comprise an 11-year activity cycle, represented by the number of sunspots. They are also characterized by their large-scale coherent structures and polarity reversals. The prominent examples are the 11-year polar field reversals and the sunspot parity rules (Hale et al. 1919). The maintenance of such a solar magnetic field can be attributed to fluid motion in the convection zone, i.e. a dynamo. Gilman (1983) and Glatzmaier (1984Glatzmaier ( , 1985 conducted studies using global three-dimensional (3D) magnetohydrodynamics (MHD) simulation. Furthermore, various 3D MHD simulations have been conducted in the past decade (Racine et al. 2011;Nelson et al. 2013;Masada Corresponding author: Ryota Shimada shimada ryota@eps.s.u-tokyo.ac.jp Fan & Fang 2014;Simitev et al. 2015;Duarte et al. 2016;Guerrero et al. 2016;Hotta et al. 2016;Käpylä et al. 2017;Strugarek et al. 2018;Hotta & Kusano 2021). Although realistic Reynolds numbers are large values of Re 10 13 and Rm 10 10 at the base of the convection zone (Ossendrijver 2003), most of the calculations used those of lower values of approximately 100-300 (Charbonneau 2020), i.e. they used large viscosity and magnetic diffusivity ( 10 12 cm 2 s −1 ) to suppress the small-scale fluctuating field and obtain a coherent large-scale magnetic field. The main drawback found by previous studies was that the energy and coherence of the large-scale magnetic field were destroyed at higher Reynolds numbers. In this context, Hotta et al. (2016) succeeded in producing a coherent largescale magnetic field at higher Reynolds numbers using high-resolution calculations. They reported that an efficient small-scale dynamo suppresses the small-scale flow, which consequently maintains the large-scale magnetic [10 4 erg cm −3 ] 16.5 9.7 9.8 (1) The number of grid points are counted on the spherical coordinate projected from Yin-Yang grid (Kageyama & Sato 2004).
(2) Imposed explicit diffusivities at the top boundary are presented in the tables. The explicit viscosity and explicit diffusivity have identical values. All cases also include the artificial viscosity reported by Rempel (2014). (3) The magnetic Reynolds number (Rm) and magnetic Prandtl number (Pm) are estimated from the energy spectra (Hotta et al. 2016). (4) The energy densities of the small-scale magnetic field (E turb ) and large-scale field (Emean) at the base of the convection zone (0.71R < r < 0.73R ) are presented. These values are considered as temporal average over 5. 5-27.4 yr. field. However, the specific physics involved in maintaining the large-scale magnetic field and reversing its polarity must be further analyzed.
Other approach to study the solar dynamo is the mean-field electrodynamics (Steenbeck et al. 1966), in which the large-scale magnetic field is considered as the mean field. This approach is useful in analyzing the physics of the 3D MHD calculations due to its simplicity. Racine et al. (2011) developed a method which extracts the mean-field parameters, that is, the α-tensor from the calculation using singular value decomposition. The mean-field approach was used to identify the major factors contributing to the induction of a large-scale magnetic field (see Brown et al. 2010;Nelson et al. 2013;Augustson et al. 2015).
This study aims to understand the physics involved in maintaining a coherent large-scale magnetic field at high Reynolds numbers. The mean-field analyses of the results of the 3D MHD simulations (Hotta et al. 2016) are conducted for this purpose. Particularly, the αtensor and β-tensor are extracted, and each term of the induction equation of the mean-magnetic field is estimated. Section 2 presents the basic settings of the simulation and the realization of a large-scale field. Section 3 presents the procedures and results of the mean-field analyses. Section 4 presents the maintenance and polarity reversal of a global-scale magnetic field at high Reynolds numbers. Lastly, Section 5 presents the conclusion.
In the course of our analysis and discussion, the radial pumping (γ r ) grabs the spotlight. The effect of the radial pumping on the solar cycle is first examined by Brandenburg et al. (1992) using the mean-field simulation. The background of this study is the existence of radial pumping in the convection zone suggested by Nordlund et al. (1992), which conducts local Cartesian simulation. Brandenburg et al. (1992) artificially select the value of γ r in their simulation at that time. Estimation of the radial pumping in 3D MHD simulation started by Ossendrijver et al. (2002) using local Cartesian simulation. In our work, the global 3D MHD simulation is used to obtain the global distribution of γ r . This allows us to estimate the distribution of induction by radial pumping and gain insight into polarity reversal. Furthermore, it is noteworthy that γ r under the several Reynolds numbers are extracted and discussed in this work.

OVERVIEW OF THE SIMULATION
This section presents a detailed description of the basic setting of the simulation and the realized large-scale field.

Basic settings
The data are calculated using a similar approach as in Hotta et al. (2016). The simulation adopts the reduced speed of sound technique (RSST) (Hotta et al. 2012) and solves the 3D MHD equations in spherical geometry (r, θ, ϕ) with gravity and rotation as follows: ρT where ρ, p, s, T , v, and B represent the density, gas pressure, specific entropy, temperature, fluid velocity, and magnetic field, respectively. The radial extent is restricted to r 1 < r < r 2 (r 1 = 0.71R , r 2 = 0.96R ). The subscript, 0, denotes the background spherically symmetric stratification; the solar standard model (Model S: Christensen-Dalsgaard et al. 1996) is employed for this calculation. The subscript, 1, denotes the fluctuation from the background stratification. ξ is the parameter in the RSST, which reduces the effective speed of sound by a factor of ξ. The gravitational acceleration, g, and radiative diffusivity, κ r , are adopted from Model S. The rotation rate is set as the solar rotation rate, |Ω 0 |/(2π) = 413 nHz (Thompson et al. 2003). Γ(r) is the cooling term, which is effective only near the surface. D and e ij represent the viscous stress tensor and strain rate tensor, respectively. A strong thermal conductivity κ = 2 × 10 13 cm 2 s −1 at the top boundary is adopted to obtain a solar-like differential rotation. The thermal conductivity, κ, explicit viscosity, ν, and explicit magnetic diffusivity, η, are set to have a radial dependence of 1/ √ ρ 0 (Fan & Fang 2014).
The analyzed cases are listed in Table 1. These cases are distinguished by the number of grid points and the imposed diffusivities. Case High contains twice as many grid points as the other cases in each direction. The explicit diffusivities, i.e., magnetic diffusivity and viscosity, are imposed in case Low, and only numerical diffusivities exist in the other cases.
The estimated magnetic Reynolds number (Rm) and magnetic Prandtl number (Pm) are listed in Table 1. The magnetic Reynolds and Prandtl numbers are evaluated using the energy spectra in cases Mid and High since no explicit diffusivities are used and the dissipations are all produced by numerical ones (see Hotta et al. 2016). The magnetic Reynolds number increases slightly from Low to Mid and significantly increases toward High.
The main differences between this simulation and Hotta et al. (2016) are the artificial diffusivity and the number of grid points. The artificial diffusivity suggested by Rempel (2014) is implemented on all physical quantities in ours, whereas one suggested by Rempel et al. (2009) is used exceptionally on density in Hotta et al. (2016) . The number of grid points for each direction in our case High is half as many as the counterpart. Figure 1 presents the snapshots of radial velocity and longitudinal magnetic field in case High. Because the obtained large-scale magnetic field obviously has axisymmetric feature, we define mean-field component of given physical quantity Q by: Q (t, r, θ, ϕ) = Q(t, r, θ, ϕ) − Q (t, r, θ),

Large-scale field
for discussions. Hereafter, and denote the operations described in Equations (6) and (7). The large-scale magnetic field is concentrated around the base of the convection zone (r < 0.8R ) (see Figure  15 in Appendix A). The temporal evolution of the mean toroidal magnetic field B ϕ at the base of the convection zone is shown in Figure 2. All the cases have spatially coherent structures at lower latitudes (|Θ| 30 • , Θ = 90 • − θ). These cases also exhibit irregular polarity reversals of the large-scale magnetic field. The reversal occurs every 5-10 yr in case Low and every 2-5 yr in the other cases.
The strengths of the large-scale and small-scale magnetic fields depend on the cases. The energy densities of the large-scale and small-scale magnetic fields are calculated as follows: Emean and the small-scale magnetic field E turb at the base of the convection zone (0.71R < r < 0.73R ) are shown. E turb is considered as the temporal average over a period of 5.5-27.4 yr, and Emean is averaged over the period in which Emean exceeds 5 × 10 4 erg cm −3 from 5.5-27.4 yr. Error bars indicate standard deviations.
where V denotes the integrated volume. A volume localized at the base of the convection zone (r < 0.73R ) is used for V in the following discussion. Figure 3 presents E mean and E turb . The energy density of the large-scale magnetic field decreases from case Low to Mid and remains at the same level from case Mid to High. Conversely, the energy density of the small-scale magnetic field monotonically increases toward case High with the increase in the magnetic Reynolds number. E mean is averaged over the period in which E mean exceeds 5 × 10 4 erg cm −3 in 5.5 yr < t <27.4 yr. The purpose of this operation is to avoid the misreading of the energy density of the large-scale magnetic field by including the period of polarity reversal. The large-scale magnetic field observed in case Mid is incoherent in time and does not vanish when the polarity is reversed, whereas that of cases Low and High is coherent in time and vanishes for  a significant amount of time ( 1 yr) when the polarity is reversed (Figure 2). Figure 4 depicts the mean flow field. It exhibits solarlike differential rotation such that the equator region rotates faster than the polar region. The shear of the differential rotation decreases in case Mid and recovers in case High. The mean meridional flow indicates the existence of two-cell circulation in each hemisphere and elongated structures parallel to the rotation axis at low latitudes. The basic structures are almost identical in all the cases.
The relative importance of the magnetic field and the flow field depends on the scale, and it can be measured using the energy spectra. The kinetic energy spectra (E KE ( )) is compared with the magnetic energy spectra (E ME ( )). Our normalization satisfies the following relations: (11) where denotes the spherical harmonic degree, v and B represent both the mean-field and fluctuating components, and ρ 0 represents the density of the background stratification. Figure 5 depicts these quantities at r = 0.72R from cases Low, Mid, and High. For cases Low and Mid, the kinetic energy exceeds the magnetic energy in almost all the spatial scales, whereas the magnetic energy exceeds the kinetic energy at a smaller scale for case High. This tendency in case High is also reported by Hotta et al. (2016), in which they concluded that the small-scale dynamo efficiently suppresses the turbulent flow, and a large-scale magnetic field is realized as a result.

MEAN-FIELD ANALYSIS
The induction equation of the mean magnetic field is given by: where E represents the so-called turbulent electromotive force (EMF), described as: The EMF can be described by an expansion of the mean magnetic field component and its first derivatives as follows: where a denotes a rank-two tensor and b denotes a rankthree tensor. To obtain a more physically meaningful form, Equation (14) can be rewritten as: where Schrinner et al. 2007). (∇ B ) (sym) represents the symmetric part of the gradient tensor, which is defined as: (17) where T is a transpose operation. The signs of β ij and δ i differ from those in Schrinner et al. (2007). α and β are symmetric rank-two tensors, γ and δ are vectors, and κ is a rank-three tensor. The α-term is used to describe the α-effect, which represents the contribution from helical flow. γ is the virtual velocity known as turbulent pumping, which advects, shears, and compresses the mean-field as if it is the physical velocity. The β term works as the turbulent diffusion.

α-tensor and β-tensor
The α-tensor, γ-vector, and β-tensor are extracted from the 3D simulation results to understand mean-field induction. The method reported in Simard et al. (2016) is employed in this procedure. This method is an extension of Racine et al. (2011) to include the first derivatives of the mean field in the fitting procedure. It is based on a linear least-squares fit of the temporal variation of the EMF to that of the mean magnetic field component, along with the component of its first derivatives. The relation between the EMF and the mean magnetic field, which is used for the linear least-squares fit, given by: a andb are pseudo-tensors and are assumed to be timeindependent. The following procedures are performed at each grid point, (r b , θ c ), and for each component, (m = r, θ or ϕ), of the EMF. We define and where k = r, θ, or ϕ. Subsequently, the fitting formula for Equation (18) can be written as: A merit function χ 2 is defined as: where N t denotes number of time steps t i of the data. The best value for Φ k is obtained by minimizing the merit function, χ 2 , using singular value decomposition (SVD). The design matrix, A, in the SVD is formed as: This matrix can be decomposed as: where U is an N t × 9 orthogonal matrix, w is a 9 × 9 diagonal matrix containing the singular values, and V is a 9 × 9 orthogonal matrix. The solution, Φ = (Φ 1 , · · · , Φ 9 ), is given by: where where σ 2 denotes the variance of the merit function, and V kl and w ll represent the elements of the V and w matrices. The relations between these pseudo-tensors (ã andb) and the true-tensors (a and b) considering the curvilinear nature of the spherical coordinate system (see Schrinner et al. 2007, for further details), are given as: The α-tensor, γ-vector, and β-tensor are automatically determined through Equations (16) and (28) afterã and b are estimated by the fitting procedure. δ-vector and κ-tensor are also obtained in this procedure. Firstly, the results of the α-tensor and γ-vector obtained through this procedure are explained for each Figure 6. Components of the a-tensor extracted from case High plotted in meridional plane. The signs of γ θ , α rθ , and α θϕ are different from those in Racine et al. (2011) because of the definition of the coordinates. case using the entire simulation interval. Figure 6 presents the result from case High. The diagonal components, α rr , α θθ , and α ϕϕ are all antisymmetric about the equator and have different signs between the surface and the base. Radial pumping, γ r , is downward (negative) in most of the convection zone and upward (positive) in the subsurface. The latitudinal pumping, γ θ , is equatorward (positive/negative in the northern/southern hemisphere) in the lower convection zone and polarward (negative/positive in the northern/southern hemisphere) in the subsurface. These properties are consistent with those of previous studies (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016). It can be observed that the global structures of the α-tensor are similar to those of the other cases. The resultant α-tensor and γ-vector do not significantly vary from those estimated using the method reported in Racine et al. (2011), which does not include the first derivatives of the mean magnetic field in the fitting procedure (see Appendix B).
The importance measure, ij , introduced by Augustson et al. (2015), is used to measure the relative importance of each component. It is defined as: and {v · v } represents the sum of the diagonal elements of the Reynolds stress tensor averaged over the duration of the simulation and over all the longitudes. ij ∼ aij vrms , E ∼ a vrms , and v rms ∼ v i v i where a is a norm of atensor, v rms is the root mean squared velocity. r 1 and r 2 denote the radii of the lower and upper boundaries of the simulation domain, respectively. Importance measure, This indicates that the upper two-by-two matrix formed by α rr , α rθ , α θθ , and γ ϕ are dominant over the others, which is consistent with Augustson et al. (2015). However, our results exhibit a larger contribution than those of γ r , which exceeds α rr , α θθ , and α rθ . In their study, the importance measures take the value, 0.054, in γr , 0.355 in αrr , 0.103 in α θθ , and 0.124 in α rθ . The relative contribution of each component of the a-tensor does not vary significantly in all the cases. Warnecke et al. (2018) extract the α-tensor and γvector based on test-field methods. These quantities in our result ( Figure 6) have smaller structure than those in Warnecke et al. (2018). The amplitudes of the structures, whose scales are smaller than 0.01R are within the range of 1σ and only the structures larger than 0.01R are reliable. γ r extracted by Warnecke et al. (2018) changes the sign between inside and outside the tangential cylinder, whereas our γ r has the same sign between these two regions. The cause of this difference might come from the difference in rotational constraint on the simulation, because the rotation rate of the simulation analyzed in Warnecke et al. (2018) is five times larger than that of ours.
Secondly, the results of the β-tensor obtained through this procedure are explained for each case using the entire simulation interval. Figure 7 shows the β-tensor extracted from case High. Negative values observed in the diagonal elements are within the range of 2σ. Note that estimation of the α-tensor is not significantly affected by the inclusion of the β-tensor (see Appendix B). To obtain the physical meaning from this β-tensor, the sum of the diagonal elements is calculated as: Figure 9. Radial plot of the turbulent magnetic diffusivity, β, explicit magnetic diffusivity, η, and explicit viscosity, ν, averaged over the latitude from cases Low, Mid, and High.
where β corresponds to the turbulent magnetic diffusivity introduced in the widely known mean-field dynamo model. Figure 8 presents the spatial structure of β from case Low to High. The result indicates that a positive value is observed through most of the convection zone and increases toward the surface. Figure  9 presents the radial plots of β averaged over the latitudes. β reaches its maximum value of 7×10 11 cm 2 s −1 for the Low, 4×10 11 cm 2 s −1 of case Mid, and 1×10 11 cm 2 s −1 of case High at depth r/R ∼ 0.91. It is clear that the turbulent magnetic diffusivity decreases with an increase in the magnetic Reynolds number. Note that we omit the results of δ-vector and κ-tensor because the values of some components have less than 1σ noise level. Unfortunately, we can not compare these results with those in the previous work by Simard et al. (2016), because there were no description of the results nor related discussion while these vector and tensor were shown in the fitting formula (see their equation 5).

Induction of the Mean Magnetic Field
A coherent large-scale magnetic field is observed in all the cases (see Figure 2 and Subsection 2.2). The elemental processes required to maintain the mean field are analyzed to determine the construction mechanism of the large-scale magnetic field. The induction equation can be decomposed as follows: (Brown et al. 2010; Equation (34) is used to search for important terms to induce the mean field. The terms on the right-hand side of Equation (34) represent the production of the mean magnetic field by the mean shear (MS), fluctuating shear (FS), mean advection (MA), fluctuating advection (FA), mean compression (MC), fluctuating compression (FC), and resistive diffusion (RD). In our results, the spatial distributions of these terms do not vary significantly from the cases, but some fluctuating structures that have a much smaller scale than the large-scale magnetic field are found in cases Mid and High. The result from case Low is presented in the following section to discuss the large-scale distribution. A time period stretching over approximately 5 yr, is selected to analyze the coherent large-scale magnetic field, when the mean magnetic field has constant polarity. The average was taken over this period. Figure 10 presents the mean magnetic field of case Low. The mean toroidal magnetic field, B ϕ , is concentrated at the base of the convection zone (r/R 0.8) and low latitudes (|Θ| 30 • ). The toroidal magnetic vector potential, A ϕ , depicted in Figure 10  B P = ∇ × ( A ϕ e ϕ ). B θ is concentrated at the base of the convection zone (r/R 0.8) and low latitudes (|Θ| 30 • ). B θ pointing in the opposite direction to the bottom one is distributed at the middle depth of the convection zone (0.8 r/R 0.9). These spatial structures of the mean magnetic field are common in the other cases, despite the difference in the amplitudes.    Figure 11 (b) demonstrates that the fluctuating shear produces negative B θ at low latitudes (|Θ| 20 • ) against the fluctuating advection. At high latitudes (|Θ| 20 • ), positive B θ is observed around the base of the convection zone (r/R 0.8), and the result presented in Figure 11 (b) is the radial average over r/R ≤ 0.8, due to which the signs of these terms are reversed against the signs at low latitudes.
The contribution of the fluctuating field, that is, EMF, is important in inducing B θ at the base of the convection zone and low latitudes, as explained earlier. The EMF is decomposed by each component of the α-tensor and γ-vector and their contribution is evaluated as: Figure 12 presents the dominant terms, that is, I αϕϕ , I α θϕ , and I γr , on the right-hand side of this formula. This indicates that B θ is produced by the effect of α ϕϕ against the counter-effects of α θϕ and γ r at the base of the convection zone and at low latitudes.

DISCUSSION
In Section 3, the mean-field parameters such as the a-tensor ( Figure 6) and the turbulent magnetic diffusivity β (Figure 9) are obtained. The structure of a large-scale magnetic field is also obtained at a steady polarity period (Figure 10), and its induction ( Figures  11 and 12) is analyzed. These results are used to infer the mechanism which maintains a large-scale field at high Reynolds numbers, which reverses their polarities. The cause of the enhancement of γ r in comparison to the previous studies (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016) is also discussed.

The Large-Scale Magnetic Field in High Magnetic
Reynolds Number The magnetic energy of the large-scale field decreases from case Low to case Mid and remains at the same level in case High (see Figure 3). This is explained by using the explicit diffusivities and turbulent magnetic diffusivities (Figure 9).
The general effects of the explicit diffusivity and turbulent magnetic diffusivity on large and small-scale magnetic fields are discussed here before discussing the results of this study. Explicit viscosity and magnetic diffusivity contribute to the smoothening of the spatial and temporal fluctuations of the velocity and magnetic field. Consequently, turbulence and small-scale dynamo are suppressed by explicit diffusivities, and a large-scale magnetic field is maintained. Conversely, the suppression of the turbulence and small-scale dynamo, which distract the large-scale magnetic fields, is less efficient when the explicit diffusivities become smaller. Consequently, the large-scale magnetic field is distracted. This implies that the energy of the large-scale magnetic field decreases while that of the small-scale magnetic field increases. The turbulent magnetic diffusion is expressed in the induction equation of a large-scale field as: whereas it is expressed in the induction equation of the small-scale field with a plus sign, as: This indicates that the energy dissipated from the largescale field by the turbulent diffusion, is converted into a small-scale field. Note that, not all the energy lost from the large-scale magnetic field due to turbulent diffusion is converted to that of the small-scale magnetic field. A lesser amount of energy present in the large-scale magnetic field is converted into a small-scale magnetic field when the turbulent magnetic diffusivity decreases. The simulation results demonstrate a decrease in the energy of the large-scale magnetic field and an increase in that of the small-scale one from case Low to case Mid (Figure 3), where the Reynolds number increases between these cases. The explicit viscosity, ν, and the magnetic diffusivity, η, between these cases demonstrate larger differences than the turbulent magnetic diffusivity (Figure 9). This result can be explained by the effect of the explicit viscosity and magnetic diffusivity, as described above. Case Mid has smaller explicit diffusivities than case Low (Figure 9), and the magnetic field of case Mid, thus forms a less large-scale structure. Consequently, the energy of the large-scale magnetic field decreases, and that of the small-scale magnetic field increase in case Mid.
However, the energy of the large-scale magnetic field remains at the same level from case Mid to case High, and the energy of the small-scale field increases. Only the turbulent magnetic diffusivity, β changes between these cases, although there is no change in the explicit diffusivities ( Figure 9). This can be attributed to the effect of the turbulent magnetic diffusivity, as explained earlier. The decrease in the turbulent diffusivity indicates the suppression of energy conversion toward the small-scale field. The maintenance of the magnetic energy of the large-scale field from case Mid to case High ( Figure 3) can thus be explained by the decrease in the turbulent diffusivity, which suppresses the energy converted from a large-scale magnetic field.
The decrease in the turbulent magnetic diffusivity in case High, results in the large-scale magnetic field being maintained. Here we refer to the previous study on the large-scale magnetic field with a high magnetic Reynolds number. Cattaneo & Tobias (2014) calculate kinematic dynamo with high magnetic Reynolds number, which gives 2.5D (two-dimensional and threecomponents) fluid velocity. They presented the suppression principle, which states that the construction of a large-scale magnetic field for a higher magnetic Reynolds number requires the suppression of the small-scale dynamo. Otherwise, the dynamo scale shifts towards a smaller one for higher magnetic Reynolds numbers, resulting in the destruction of the large-scale fields. Cattaneo & Tobias (2014) conducted kinematic calculations, in which the feedback from the magnetic field to the flow field was not considered. Our results for the magnetic energies of the large-scale and small-scale fields do not correspond to the results of the suppression principle, in that the large-scale magnetic field in case High is maintained without suppressing the small-scale dynamo. This is because the magnetic field exceeds the equipartition strength with velocity in the small-scale field in case High (Figure 5), and kinematic approximation is no longer applicable. In this case, the Lorentz force feedback from the small-scale magnetic field to the velocity becomes important. This effect is thought to be incorporated into the turbulent magnetic diffusivity through the small-scale velocity in the turbulent EMF (Equations (13) and (15)). Therefore, our study states that a largescale magnetic field can be maintained even when the small-scale dynamo works efficiently with high magnetic Reynolds numbers. Hotta et al. (2016) reported that an efficient small-scale dynamo suppresses the smallscale flow which destroys the large-scale magnetic field. Our analysis proposes that the destruction can be incorporated into the turbulent magnetic diffusivity in the mean-field electrodynamics. The estimation of the turbulent magnetic diffusivities quantitatively supports the effect proposed by Hotta et al. (2016).
In our simulation, the turbulent magnetic diffusivity significantly decreases with the number of grid points. We cannot confirm that this tendency continues in the further large numbers of grid points. Whether the suppression of turbulent magnetic diffusivity and consequent maintenance of the large-scale magnetic field occur at further high Reynolds numbers are needed to be investigated.

Polarity Reversal
All the cases in our simulations show the polarity reversal of the large-scale magnetic field (see Figure  2). The large-scale magnetic field comprises the spatial structures presented in Figure 10 during the interreversal period. Both B ϕ and B θ are concentrated at the base of the convection zone, where B θ pointing in the opposite direction, is distributed at the middle depth of the convection zone. B ϕ is produced by the Ω-effect and fluctuating advection (FA) against the destruction effect from the fluctuating shear (FS), whereas B θ is produced by the effect of α ϕϕ against the effects of α θϕ Figure 14. Radial plot of the terms which contribute to γr in Equation (38). These are averaged over the entire latitudes. The average period is from 5.5 yr to 27.4 yr. and γ r at the base of the convection zone and low latitudes (see Subsection 3.2).
The fluctuating shear reverses the sign of B ϕ (Figure 11) based on the obtained relation among the field induction effect. B θ which has reversed sign is dominantly induced by α θϕ at the base of the convection zone and at low latitudes. This brings the conclusion that α θϕ reverses the sign of B θ .
The Ω-effect can amplify B ϕ , and α ϕϕ can amplify B θ at the base of the convection zone after the signs of both B ϕ and B θ are reversed by these effects. Our analysis of the cycle phase dependence of the αtensor and γ-vector provides additional insights into the triggers of the reversal of B θ . The α-tensor and γvector are estimated using whole periods of the maxima and minima through the procedure described in Subsection 3.1 to estimate the cycle phase dependence. The definitions of the minima and maxima are presented in Appendix D. The radial pumping in the magnetic minima is observed to be stronger than the maxima ( Figure  13). α ϕϕ and α θϕ which contribute to the induction of B θ do not demonstrate a clear difference between the maxima and minima (see Appendix D). Consequently, the strengthened downflow produced by γ r during the minima may trigger the polarity reversal of B θ .

The Enhancement of Radial Turbulent Pumping
Our results for the α-tensor and γ-vector differ from those of the previous studies (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016) in relative contributions from γ r (Subsection 3.1). The contribution of γ r exceeds the components of the α-tensor. The model of the γ-vector in terms of the small-scale velocity and magnetic field under the stratification proposed by Rädler et al. (2003) as: is employed to interpret this difference, where τ 0 represents the correlation time. Figure 14 presents each term in Equation (38). This indicates that the amplitude of γ r is primarily enhanced by the contribution from the small-scale magnetic field. The model of the α-tensor and γ θ is also verified, and the results demonstrate that the amplitudes of these terms are mainly determined by the contribution from the small-scale velocity and do not exhibit significant differences between the cases (see Appendix E). Our enhancement of γ r when compared to the previous studies (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016) can be attributed to the fact that only γ r is affected by the small-scale magnetic field since the small-scale magnetic field is significantly enhanced in our simulation.

CONCLUSION
Mean-field analyses are conducted on the results of the global 3D MHD calculations to understand the physics which maintain a coherent large-scale magnetic field at high Reynolds numbers (Hotta et al. 2016). Three cases with increasing Reynolds numbers are considered (Figure 2). The energy of the large-scale magnetic field decreases from case Low to case Mid and remains at the same level in case High (Figure 3 (a)).
The turbulent magnetic diffusivity in each case is obtained through mean-field analysis (Racine et al. 2011;Simard et al. 2016). The results (see Figure 8 and 9) demonstrate that the turbulent diffusivity monotonically decreases with an increase in the Reynolds numbers. The maximum value of the turbulent magnetic diffusivity is 7×10 11 cm 2 s −1 in case Low and 1×10 11 cm 2 s −1 in case High, at the depth, r ∼ 0.91R .
The reduction of the turbulent diffusivity in the high-Reynolds number case leads to inefficient transformation since the turbulent diffusivity transforms the magnetic energy of the large-scale field into that of the smallscale field (see Subsection 4.1). This is considered to be the primary cause for the maintenance of the large-scale magnetic fields at high Reynolds numbers (case High).
The α-tensor and turbulent pumping are also obtained through mean-field analysis. The spatial structures of these components are consistent with those of the previous studies (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016) for instance: (1) The diagonal component of the α-tensor is antisymmetric about the equator and has different signs between the surface and the base.
(2) Radial pumping, γ r , is downward (negative) in most of the convection zone and upward (positive) in the subsurface. (3) The latitudinal pumping, γ θ , is equatorward (positive/negative in the northern/southern hemisphere) in the lower convection zone and polarward in the subsurface.
Our results for the α-tensor and γ-vector differ from those of previous studies (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016) in the relative contributions from γ r . The model of the α-tensor and γ-vector under the stratification proposed by Rädler et al. (2003) in our simulation indicates that only γ r is significantly affected by the small-scale magnetic field. Our enhancement of γ r in comparison to that of the previous studies can be attributed to this effect since the small-scale magnetic field is significantly enhanced in our simulation.
B ϕ is concentrated at the base of convection zone (r/R 0.8) and low latitudes (|Θ| 30 • ) during the period in which the polarity of the large-scale field is constant. B θ is concentrated at the base of the convection zone (r/R 0.8) and low latitudes (|Θ| 30 • ), where B θ pointing in the opposite direction, is distributed in the middle of the convection zone (0.8 r/R 0.9) during this period. Our estimation of the induction equation demonstrates that the Ω-effect (MS) and advection by fluctuating flow (FA) produce B ϕ , and shearing due to the fluctuating flow (FS) works against the production at low latitudes in the base of the convection zone. B θ in this area is produced by shearing due to the fluctuating flow (FS), and advection due to the fluctuating flow (FA) works against it. The contribu-tions from the α-tensor and turbulent pumping in producing B θ are also verified. The results demonstrate that the effect of α ϕϕ produces B θ at low latitudes in the base of the convection zone against the effects of α θϕ and γ r .
Based on the obtained relation among the field induction effect, the fluctuating shear (FS) reverses the sign of B ϕ (Figure 11), and α θϕ plays a significant role in reversing the sign of B θ in the base of the convection zone (Figure 12). Ω-effect and α ϕϕ can amplify both B ϕ and B θ after the sign of B ϕ and B θ are reversed. The cycle phase dependence of γ r ( Figure  18) indicates that the strengthened downflow due to γ r during the minima may trigger the polarity reversal of B θ . This effect may be incorporated into the fluctuating advection (FA) in Equation (34) at low latitudes (|Θ| 20 • ). We plan to conduct the mean-field dynamo simulation to quantitatively estimate the effect of α ϕϕ and γ r on the polarity reversal.  Figure 15 depicts the temporal evolution of the radial distribution of the mean toroidal magnetic field, B ϕ , at latitude, Θ = 20 • , where the field is concentrated. The large-scale magnetic field is observed to be concentrated at the base of the convection zone (r < 0.8R ) in all the cases.

B. ALPHA TENSOR OBTAINED WITHOUT USING FIRST DERIVATIVES OF MEAN MAGNETIC FIELD IN THE FITTING PROCEDURE
Subsection 3.1 presents the result of the α-tensor and the γ-vector ( Figure 6) obtained from the fitting method, which includes the first derivatives of the mean magnetic field in the fitting procedure (Equation 18). The result obtained by the method which does not include the first derivatives of the mean magnetic field in the fitting procedure is presented in this instance, as originally proposed by Racine et al. (2011). In this method, the fitting formula for Equation 18 is described as follows: E i (t, r, θ) =ã ij (r, θ) B j (t, r, θ), The fitting procedure is almost identical to that in Subsection 3.1 except for the designed matrix in SVD (see Racine et al. 2011, for further detail). The α-tensor and γ-vector are automatically determined through Equations (16) and (28) afterã is estimated by the fitting procedure. Figure 16. Components of the a-tensor extracted from case High plotted in meridional plane. The signs of γ θ , α rθ , and α θϕ differ from those in Racine et al. (2011) owing to the definition of the coordinates. The color bar and scaling are uniform across all the panels. Figure 16 presents the obtained α-tensor and γ-vector, and it can be observed that there is no significant difference between Figure 6 and Figure 16.  The maxima and minima are defined as the period in which the magnetic energy is larger and smaller than its temporal average over the entire simulation interval. The minima includes the reversal phase in this definition. Only the northern hemisphere is considered for this definition along with the related arguments in Sub-section 4.2 due to the weak parity synchronization in our simulation results.
The radial component of the turbulent pumping γ r during the maxima and minima are shown for the High case in Figure 13 (b) and for the other cases in Figure 18. All the cases exhibit the trend that downward pumping is stronger in the minima than in the maxima by approximately 1σ level. Figure 19 shows α ϕϕ and α θϕ during the maxima and the minima for the case High. Neither α ϕϕ and α θϕ do not indicate a clear difference between the maxima and minima.

E. MODEL OF THE MEAN-FIELD COEFFICIENTS
The models of the α-tensor are used in terms of the small-scale velocity and magnetic field proposed by Rädler et al. (2003), in order to interpret the difference between the mean-field coefficients in our study and those of previous works (Racine et al. 2011;Augustson et al. 2015;Simard et al. 2016),as follows: where τ 0 is the correlation time. We also use their model for the γ-vector, that is, Equation (38). The results for γ r are presented in Figure 14 and show a significant contribution from the small-scale magnetic field (Subsection 4.3). Conversely, the results for α and γ θ (Figure 20) exhibit a weak contribution from the small-scale magnetic field.