© The Authors 2026

1. Introduction

The theory of how magnetic fields appear in turbulent flows of conductive fluid has wide applications, particularly in explaining the magnetic fields observed in many astrophysical objects (see, e.g., Brandenburg & Subramanian 2005; Jouve et al. 2008; Moss et al. 2013; Brandenburg et al. 2024). There are two possible mechanisms to generate a small-scale magnetic field, that is, a field with a characteristic scale much smaller than the integral scale of turbulence (Karak & Brandenburg 2016). First, there is a cascade mechanism of generation that is initiated at large scales under the influence of global shear and convective flows. Second, there is a faster generation process caused by small-scale turbulence. The latter has been studied in numerous papers (see, e.g., Zeldovich et al. 1990; Falkovich et al. 2001; Brandenburg et al. 2012) and is the subject of this paper.

Conductive fluids can be classified by their magnetic Prandtl number,

$$\begin{array}{c}\hfill Pm=\nu /\eta ,\end{array}$$(1)

where ν is viscosity, and η is magnetic diffusivity. For Pm ≫ 1, the resistive scale lies deep inside the viscous range of turbulence (Batchelor 1950; Zel’dovich et al. 1984; Chertkov et al. 1999). This limit is realized, for instance, in the interstellar medium and in processes of star formation (Han 2017). We concentrate on the other case of fluids with low and intermediate Pm. This means that the resistive scale exceeds the viscous scale, so that the excitation of magnetic fluctuations is driven by the inertial-range or bottleneck velocity fluctuations. Low-Pm processes take place in the interior of solar-type stars and planetary convective envelopes (Petrovay & Szakaly 1993).

The theoretical description of the process in the incompressible fluid is based on the Kazantsev equation, which relates the evolution of the magnetic field pair correlator and the velocity structure function. Experimental possibilities are rather limited in terrestrial conditions. Direct numerical simulations (DNS) of low Pm, unlike the case of high Prandtl numbers, face difficulties in modeling magnetic field advection in the inertial range of turbulence and demonstrate some discrepancies between different simulations (see, e.g., Warnecke et al. 2023, Fig. 2).

However, during the past decades, significant progress has been made in the theoretical understanding of the properties of turbulence (L’vov et al. 1997; Donzis & Sreenivasan 2010; Biferale et al. 2011; Iyer et al. 2020) and of the environment and technics of DNS (Iskakov et al. 2007; Warnecke et al. 2023; Rempel et al. 2023; Brandenburg et al. 2023; Yeung et al. 2025), which has resulted in high resolution and advanced to lower values of Pm. This progress gives us hope that we can proceed from a qualitative correspondence to a more or less thorough quantitative comparison at least between theory and DNS.

In the study of turbulence, after the theoretical model formulation (for decaying turbulence) in the classical work (Kolmogorov 1941) and its development for stationary turbulence (Novikov 1965), the precision testing of the model results was performed in experiments (Moisy et al. 1999). The coincidence of theoretical predictions and experimental measurements has, citing Moisy et al. (1999), demonstrated ‘the relevance of isotropic homogeneous turbulence state approximation, used in almost all theoretical approaches to turbulence’. For a low-Pm dynamo, such a program seems hardly possible because of experimental difficulties. However, in the absence of experimental approbation, an accurate quantitative comparison with DNS becomes even more important. It could validate both of the approaches: it would verify the resolution of DNS, and verify the applicability of theoretical assumptions of δ–time correlation (Tobias et al. 2012) and Gaussianity of the turbulent velocity field. It could also test nontrivial effects resulting from non-Gaussianity, for instance, from the weakening of the magnetic field generation, and corrections for the magnetic energy spectrum (Kopyev et al. 2022a,b, 2024). Finally, if the theory were confirmed for moderate to low magnetic Prandtl numbers, it might be scaled to extremely low Pm, which are unattainable by current facilities and typical for astrophysical objects.

The correct interpretation of numerical data for velocity correlators for their comparison with the theory is one of the problems in this path. In the Kazantsev theory, velocity correlators are assumed to be δ-correlated in time, and the velocity distribution acts in the Kazantsev equation by means of the multiplier b(ρ) in the structure function,

$$\begin{array}{c}\hfill {⟨\delta {\mathit{v}}_{\Vert }(\mathit{\rho },t)\delta {\mathit{v}}_{\Vert }(\mathit{\rho },t\prime )⟩}_{\mathrm{Kaz}}=2b(\rho )\delta (t-t\prime ),\end{array}$$(2)

where the longitudinal velocity increments in two near points are

$$\begin{array}{c}\hfill \delta {\mathit{v}}_{\Vert }(\mathit{\rho },t)=(\mathit{v}(\mathit{r}+\mathit{\rho },t)-\mathit{v}(\mathit{r},t))\frac{\mathit{\rho }}{\rho }.\end{array}$$(3)

The correlator is independent of t, r, and of the direction of ρ because turbulence is homogeneous and isotropic.

To restore the amplitude b(ρ) from a given correlator, we can use the expression (Vainshtein 1980; Kichatinov 1985)

$$\begin{array}{c}\hfill b(\rho )=\frac{1}{2}{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}}d\tau ⟨\delta {\mathit{v}}_{\Vert }(\mathit{\rho },\tau )\delta {\mathit{v}}_{\Vert }(\mathit{\rho },0)⟩.\end{array}$$(4)

For the δ -correlated process, there is no difference whether δv_∥ is calculated at a fixed point r of space or at an arbitrarily moving point r(t); the correlator is the same. The same result is also obtained for the quasi-Lagrangian reference frame, that is, r(t) tracing one arbitrary fluid particle (L’vov et al. 1997). This frame is distinguished, for example, by the fact that the Kazantsev predictions must be consistent with the results obtained from the evolutionary approach (Zel’dovich et al. 1984; Chertkov et al. 1999) wherever both of the models are applicable; and the evolutionary approach explicitly considers the quasi-Lagrangian frame.

In real turbulent flows, including those in DNS, all structure functions have a nonzero correlation time, and the values of b(ρ) calculated in different frames (i.e., for different r(t) in Eq. (3)) differ essentially. The analysis of data from John Hopkins database shows that at high Reynolds numbers, for all scales ⟨δv_∥(ρ, τ)δv_∥(ρ, 0)⟩^Eul ∝ 1/τ (Kopyev et al., in prep.). If this is so, the logarithmic divergence of Eq. (4) corresponds to a growth of the magnetic field that is faster than exponential; this would contradict the Oseledets theorem (Oseledec 1968). We note that, unlike in the Eulerian case, the quasi-Lagrangian b(ρ) converges well: ⟨δv_∥(ρ, τ)δv_∥(ρ, 0)⟩^Lagr decays exponentially with respect to τ (Biferale et al. 2011).

Consider now some experiment or simulation. We determine below which b(ρ) should be substituted in the Kazantsev equation in order to compare the theory results with the experiment. In other words, we determine which b(ρ) simulates the correct amplitude for the effective δ-correlator. The choice of the Eulerian frame (r = const) in Eq. (3), although it might seem most natural, is in fact not preferable and is probably just incorrect. This might have been the reason that prevented Mason et al. (2011) from obtaining a good correspondence between theory and DNS. A consistent solution of the problem must be based on an accurate dynamical analysis and is far from being performed yet. We assume, however, that the correct answer is to choose the quasi-Lagrangian frame. The grounds for this assumption are as follows:

– First, from the renewing model (Zel’dovich et al. 1984), it follows (Zeldovich et al. 1990; Rogachevskii & Kleeorin 1997) that if the correlation time is shorter than all other characteristic timescales, the quasi-Lagrangian frame gives the correct result: the Kazantsev approach in this case can be verified by evolutionary models. This is even the case near the generation threshold: the characteristic time for the magnetic field evolution is long, so the δ approximation for velocity correlations works well.

– Second, for an arbitrary correlation time but only for high magnetic Prandtl numbers (which restricts the consideration to the viscous range of scales), the choice of the quasi-Lagrangian frame is also correct (Vainshtein 1980; Kichatinov 1985). Recently, Il’yn et al. (2022) have shown the quasi-Lagrangian correlators to be equivalent to effective δ-correlators in long-time approximation.

By analogy with the limit cases, it is reasonable to imply the quasi-Lagrangian velocity increments in Eqs. (3) and (4). The relevance of this assumption is validated in this paper by means of a direct comparison of the predictions of the theory based on the Kaznatsev equation with the velocity structure function b(ρ) determined in this way with the results of DNS. The good agreement is the proof of our choice of the quasi-Lagrangian structure function.

We calculate the critical magnetic Reynolds number and the increments of the magnetic field correlator near the critical regime of generation. For this purpose, we solve the Kazantsev equation with quasi-Lagrangian velocity structure function (4). The basic notations and equations are introduced in Sect. 2.

Unfortunately, there are not many data on quasi-Lagrangian velocity structure functions. We therefore considered two different cases: the Taylor Reynolds number Re_λ = 140 (Sect. 3), for which DNS data on the magnetic field (Schekochihin et al. 2007) and velocity correlators (Biferale et al. 2011; Donzis & Sreenivasan 2010) are available, and the limit of the infinite Reynolds number (Sect. 4), where the shape of the velocity correlator is chosen based on theoretical reasons. For this case, we use different approximations and consider the effects produced by the declination of the velocity structure function exponent from the Kolmogorov scaling in the inertial range of turbulence (L’vov et al. 1997; Iyer et al. 2020) and the properties of the transition range between the inertial and largest-eddy scales. We compare the results with the numerical data Iskakov et al. (2007), Schekochihin et al. (2007), Brandenburg et al. (2018), Warnecke et al. (2023) and find a quite good correspondence. We also show that the Reynolds-dependent intermittency of velocity structure functions (Iyer et al. 2020) can explain the observed decrease in the critical Reynolds number at small Prandtl numbers (Warnecke et al. 2023). We calculate the magnetic field correlator growth rate in the vicinity of the generation threshold and compare the results with the DNS (Sect. 5).

Finally, we summarize the results of the paper and comment on further prospects and possible improvements in comparison of the Kazantsev theory predictions with experiments and simulations. In particular, we indicate the properties of quasi-Lagrangian turbulence that can be found from DNS together with the magnetic generation properties, to make the comparison with the theory quantitative (not only qualitative) and precise.

2. Basic equations and parameters

The evolution of the magnetic field B(r, t) is governed by the induction equation,

$$\begin{array}{c}\hfill \frac{\partial \mathit{B}(\mathit{r},t)}{\partial t}=\mathrm{\nabla }\times [\mathit{v}(\mathit{r},t)\times \mathit{B}(\mathit{r},t)]+\eta {\mathrm{\nabla }}^2\mathit{B}(\mathit{r},t).\end{array}$$(5)

Since the back-reaction of the magnetic field on the flow is quadratic in field strength and the initial (seed) magnetic field is weak, we can neglect its effect on the velocity dynamics. In this kinematic regime, the magnetic field acts as a passive vector field advected by the flow. The solenoidal velocity field v(r, t) is treated as a prescribed stochastic field with stationary statistics.

The correlation function of the magnetic field

$$\begin{array}{c}\hfill G(\rho ,t)=⟨B_{\Vert }(\mathit{r}+\mathit{\rho },t)B_{\Vert }(\mathit{r},t)⟩,B_{\Vert }=\mathit{B}·\mathit{\rho }/\rho \end{array}$$(6)

is assumed to be independent of r in a homogenous isotropic flow. Its evolution under the assumption (2) is described by the equation (Kazantsev 1968)

$$\begin{array}{c}\hfill \frac{\partial }{\partial t}G(\rho ,t)=2S(\rho )(G{″}_{\rho \rho }+\frac{4G{\prime }_{\rho }}{\rho })+2S\prime G{\prime }_{\rho }+2(S″+4\frac{S\prime }{\rho })G,\end{array}$$(7)

where

$$\begin{array}{c}\hfill S(\rho )=\eta +\frac{1}{2}b(\rho ),\end{array}$$

and the function b(ρ) is defined in (4). We are interested in the exponential behavior of the magnetic field, and we therefore searched for a solution in the form

$$\begin{array}{c}\hfill G=e^{\gamma t}\psi (\rho )/({\rho }^2\sqrt{S}).\end{array}$$

Then, ψ satisfies the Schrodinger-type equation (Kazantsev 1968)

$$\begin{array}{cc}& \psi {″}_{\rho \rho }=\frac{\gamma }{2S(\rho )}\psi +U(\rho )\psi ,\hfill \end{array}$$(8)

$$\begin{array}{cc}& U=-\frac{1}{{\rho }^2}(\frac{3\sigma (\sigma +4)+1}{4}+\frac{\rho \sigma \prime }{2}),\hfill \end{array}$$(9)

$$\begin{array}{cc}& \sigma (\rho )=\frac{\mathrm{d}lnS}{\mathrm{d}ln\rho }-1\hfill \end{array}$$(10)

with the boundary condition ψ(0) = 0,  ψ(∞) < ∞. We solved this equation numerically by means of a modification of the stochastic quantization idea (Il’yn et al. 2021, Appendix A).

To compare the theory with the results of DNS, we normalized the dimensional parameters. Generally, an isotropic hydrodynamical flow is completely characterized by three dimensional parameters (Novikov 1965; Moisy et al. 1999):

$$\begin{array}{c}\hfill \nu ,\epsilon =\nu \langle \frac{\partial v_i}{\partial x_k}\frac{\partial v_i}{\partial x_k}\rangle ,u\prime =\sqrt{⟨v_i{v}_i⟩/3}=v_{\mathrm{rms}}/\sqrt{3},\end{array}$$

where ν is the viscosity, ε is the total energy flux from larger to smaller scales, and v_rms is the volume integrated root-mean-squared velocity. From these three parameters, we composed one universal dimensionless combination, for example,

$$\begin{array}{c}\hfill R{e}_{\lambda }\equiv \frac{u\prime \lambda }{\nu }=\frac{\sqrt{15}{u\prime }^2}{\sqrt{\epsilon \nu }},\end{array}$$(11)

where $\lambda =u\prime \sqrt{15\nu /\epsilon }$ is the Taylor scale; this parameter is called the Taylor Reynolds number. However, the classical Reynolds number is often used,

$$\begin{array}{c}\hfill Re=v_{\mathrm{rms}}L/\nu ,\end{array}$$

where L is the pumping scale, or the largest-eddies scale, and is not defined universally. This uncertainty produces difficulties in matching different experimental and/or DNS results and their comparison with the theory (e.g., Brandenburg et al. 2018 reported some obstacles when comparing their results with those of Schekochihin et al. 2007). We defined

$$\begin{array}{c}\hfill R{e}^{\mathrm{Sch}}=\frac{R{e}_{\lambda }^2}{30}=\frac{{u\prime }^4}{2\epsilon \nu }.\end{array}$$(12)

This corresponds to the relations between Re and Re_λ used by Schekochihin et al. (2007).

The magnetic properties of the fluid can be described by the dimensionless magnetic Prandtl number (1). The magnetic Reynolds number

$$\begin{array}{c}\hfill Rm=Pm·Re\end{array}$$(13)

is also generally used. Unlike Pm, it depends on the choice of L in the definition of the Reynolds number.

We are interested in the stability condition Pm_c(Re_λ), or Rm_c(Re) such that there is a solution of (8) with γ = 0 and there is no solution with γ > 0 for smaller Rm. We note that the second relation is conventional, but the first one is more universally defined.

In the vicinity of this stability curve, γ is known to depend log-linearly on Rm (Rogachevskii & Kleeorin 1997; Kleeorin & Rogachevskii 2012),

$$\begin{array}{c}\hfill \gamma \propto ln(Rm/R{m}_c).\end{array}$$

Based on Eq. (8), we calculated the proportionality coefficient numerically and compared it with the results of Warnecke et al. (2023). To solve Eq. (8), we assumed some model for b(ρ). This produced one more difficulty, since b(ρ) is the quasi-Lagrangian correlation function and is difficult to measure (Biferale et al. 2011). By definition of the Lagrangian correlation time τ_c(ρ) introduced by L’vov et al. (1997), we have

$$\begin{array}{c}\hfill b(\rho )=⟨{(\delta v_{\Vert })}^2⟩(\rho ){\tau }_c(\rho ).\end{array}$$(14)

In Biferale et al. (2011), Donzis & Sreenivasan (2010), the simultaneous structure function ⟨(δv_∥)²⟩ and the correlation time τ_c were found from DNS for Re_λ = 140; for much higher Re_λ, we used theoretical considerations (L’vov et al. 1997) that argue that both ⟨(δv_∥)²⟩ and τ_c are power-law functions of ρ inside the inertial range.

3. Re_λ = 140: Moderate Re and Pm

In this section, we refer to the results of Biferale et al. (2011), Donzis & Sreenivasan (2010) for data on velocity statistics and Iskakov et al. (2007), Schekochihin et al. (2007) for data on magnetic generation. Fortunately, all these papers contained the DNS performed for the same Reynolds number Re_λ = 140. The critical Prandtl number is not small for this Reynolds number, so the bottleneck region and even the viscous range of scales can affect the results and have to be taken into account. The main parameters of the DNS performed by Biferale et al. (2011) are listed in Table 1.

Biferale et al. (2011) calculated the dependence of the Lagrangian correlation time τ_c on ρ,

$$\begin{array}{cc}& {\tau }_c(\rho )=C_1{t}_{\nu }\frac{{(1+{(C_2\rho /r_{\nu })}^2)}^{1/3}}{{(1+{(C_3\rho /L)}^2)}^{1/3}},\hfill \end{array}$$(15)

where

$$\begin{array}{c}\hfill C_1=2.15,C_2=6.75·{10}^{-2},C_3=2.94,\end{array}$$

r_ν = ν^3/4ε^−1/4 is the Kolmogorov viscous scale, and $t_{\nu }=\sqrt{\nu /\epsilon }$ is the viscous characteristic time.

The second-order simultaneous velocity structure function is well investigated; in the viscous and inertial ranges, its behavior is universal and is very well described by the generalized Batchelor approximation (Donzis & Sreenivasan 2010),

$$\begin{array}{c}\hfill ⟨{(\delta v_{\Vert })}^2⟩=\frac{\epsilon }{15\nu }\frac{{\rho }^2}{(1+{(C_B\rho /r_{\nu })}^q{)}^{\frac{2-{\zeta }_2}{q}}},\rho \ll L,\end{array}$$(16)

$$\begin{array}{c}\hfill C_B=7.6\times {10}^{-2},{\zeta }_2=0.67,q=1.82.\end{array}$$

It is easy to see that in the ultraviolet regime, it coincides with Kolmogorov’s analytical formula, $⟨{(\delta v_{\Vert })}^2⟩=\frac{\epsilon }{15\nu }{\rho }^2$, ρ ≪ r_ν, while inside the inertial range, it demonstrates a power-law behavior, ⟨(δv_∥)²⟩∝ρ^ζ₂, L ≫ ρ ≫ r_ν/C_B. We note that the coefficients C_B in Eq. (16) and C₂ in Eq. (15) are both comparably small. This corresponds to the existence of the bottleneck transition region between the viscous and the inertial ranges: viscosity is essential at scales ρ ∼ r_ν/C₂ ≃ r_ν/C_B, which are significantly larger than the viscous scale.

For our purpose, we enlarged the limits of applicability of Eq. (16) for larger scales up to the integral scale, where ⟨(δv_∥)²⟩ = 2u′², ρ ≥ L. In accordance with the idea of extended self-similarity (ESS), this can also be done by means of a fractionally rational approximation¹,

$$\begin{array}{cc}& {⟨{(\delta v_{\Vert })}^2⟩}_{(\mathrm{main})}=\frac{\epsilon }{15\nu }\frac{{\rho }^2}{(1+{(C_B\rho /r_{\nu })}^q{)}^{\frac{2-{\zeta }_2}{q}}(1+{(C_m\rho /L)}^2{)}^{\frac{{\zeta }_2}{2}}},\hfill \end{array}$$(17)

where C_m = 1.6 follows from the matching with the integral scale.

The two transitional regions, the bottleneck and the range between the inertial and integral scales, are presented in Eq. (17) by the two brackets in the denominator. To estimate the accuracy of the approximation and to visualize the contributions of the two transition regions, we consecutively substituted a piecewise power-law function for each of the brackets,

$$\begin{array}{cc}\hfill {⟨{(\delta v_{\Vert })}^2⟩}_{(\mathrm{Sharp}-\mathrm{int})}& =\frac{\epsilon }{15\nu }\frac{{\rho }^2}{(1+{(C_B\rho /r_{\nu })}^q{)}^{\frac{2-{\zeta }_2}{q}}}\theta (L/a_{\mathrm{int}}-\rho )\hfill \\ \hfill & +2u\prime \theta (\rho -L/a_{\mathrm{int}})\hfill \end{array}$$(18)

and

$$\begin{array}{c}\hfill {⟨{(\delta v_{\Vert })}^2⟩}_{(\mathrm{Sharp}-\mathrm{bot})}=\frac{\epsilon }{15\nu }\frac{{\rho }^2}{(1+{(C_m\rho /L)}^2{)}^{\frac{{\zeta }_2}{2}}}\theta (\rho -a_{\mathrm{bot}}{r}_{\nu })\\ \hfill +\frac{\epsilon }{15\nu }{\rho }^2\theta (a_{\mathrm{bot}}{r}_{\nu }-\rho ).\end{array}$$(19)

The coefficients a_bot = 13.1 and a_int = 1.6 were also found from the matching condition, as were C_B and C_m. The value of a_int coincides to high accuracy with C_m, which indicates that the two approximations are rather close.

The results are presented in Table 2. The theoretical predictions are rather close to the numerical results; the use of the Eulerian structure function b(ρ) would give an estimate of Rm_c at least an order of magnitude lower. The fact that the theoretical prediction for Rm_c is lower than the experimental result may be caused by the assumption of Gaussianity of the flow: when the declination from Gaussianity is taken into account, that is, the third-order correlator of the velocity, the result can increase (Kopyev et al. 2024).

The critical values of the magnetic Reynolds number and Prandtl number are close for all the three models, but the results still differ by about 15%. This indicates the effect of the transitional regions on the near-threshold generation. The importance of the bottleneck region is determined by the fact that the magnetic Prandtl number is not small enough to neglect the effect of the viscous range. From the value of C_B in Eq. (17) and C₂ in Eq. (15), it follows that the effect of the viscosity reaches about 15r_ν, while for Pm = 0.25, the magnetic diffusivity scale is only about 3r_ν, so it lies deep inside the bottleneck range. It is natural to assume that for smaller Pm, the effect of viscosity on the magnetic field generation decreases. In contrast, the effect of the external transition between the inertial and the integral ranges can only become even stronger for higher Reynolds numbers.

4. Very high Reynolds numbers

The lack of experimental or DNS data on the hydrodynamic Lagrangian correlation time prevented us from determining b(ρ) for other Reynolds numbers, and, hence, from calculating Rm_c for them. However, for a well-developed turbulence with a wide inertial range, we assumed in accordance with the general theoretical approach to turbulence that b(ρ) has a universal shape regardless of the details of the flow. To normalize the function, we took the largest scales. For these scales,

$$\begin{array}{c}\hfill \underset{\rho \to \infty }{lim}⟨{(\delta v_{\Vert })}^2⟩=\frac{2}{3}{v}_{\mathrm{rms}}^2=2{u\prime }^2,\end{array}$$(20)

and the integral Lagrangian correlation time $T_{\mathrm{L}}=\underset{\rho \to \infty }{lim}{\tau }_c$ can be found from the approximation based on Sawford’s second-order stochastic model (Sawford 1991; Sawford et al. 2008),

$$\begin{array}{c}\hfill T_{\mathrm{L}}=\sqrt{\frac{\nu }{\epsilon }}\frac{2R{e}_{\lambda }}{\sqrt{15}{C}_0}=\frac{2}{C_0}\frac{u{\prime }^2}{\epsilon },\end{array}$$(21)

where C₀ is the Kolmogorov constant. It is rather difficult to find (Lien & D’Asaro 2002; Ouellette et al. 2006; Zimmermann et al. 2010; Uma-Vaideswaran & Yeung 2025), and we used the value obtained by Sawford & Yeung (2011),

$$\begin{array}{c}\hfill C_0=6.9\pm 0.2.\end{array}$$(22)

In what follows, we therefore consider different models for b(ρ) that satisfy the normalization condition consistent with Eq. (14)

$$\begin{array}{c}\hfill b_{\infty }\equiv b(\infty )=2u{\prime }^2{T}_{\mathrm{L}}.\end{array}$$(23)

According to Eqs. (21), (11),(1),

$$\begin{array}{c}\hfill \frac{b_{\infty }}{2\eta }=Pm\frac{R{e}_{\lambda }^2}{15}\frac{2}{C_0}\propto Rm.\end{array}$$(24)

In what follows, we use the dimensionless value

$$\begin{array}{c}\hfill X={\left(\frac{b_{\infty }}{2\eta }\right)}^{1/(1+s)}\end{array}$$(25)

to characterize the generation properties of a flow. Here, s is the scaling exponent of ⟨|δv_∥|⟩. The value of X describes the relative width of the range of scales that can affect the generation. We note that for any particular series of DNS or experiments, Rm ∝ X^1 + s, but because the definition of the integral scale L is uncertain, X is more accurately defined than Rm.

4.1. Sharp model

We started with the natural and simple theoretical model proposed by Vainshtein & Kichatinov (1986), which focuses on the main features of turbulence (hereinafter, model ‘Sharp’). It describes the power-law behavior of b(ρ) in the inertial range and its approximate constancy at large scales,

$$\begin{array}{c}\hfill b{(\rho )}_{\mathrm{Sharp}}=b_{\infty }\{\begin{array}{cc}{(\rho /\mathrm{\Lambda })}^{1+s},\hfill & \rho <\mathrm{\Lambda }\hfill \\ 1\hfill & \rho >\mathrm{\Lambda }\hfill \end{array},s=1/3.\end{array}$$(26)

This model was used in numerous papers (see, e.g., Rogachevskii & Kleeorin 1997; Vincenzi 2002; Boldyrev & Cattaneo 2004; Arponen & Horvai 2007; Schober et al. 2012; Kleeorin & Rogachevskii 2012; Kopyev et al. 2026). The viscous range is not of interest for a dynamo study at very high Reynolds numbers, since it is far below the magnetic diffusion scale and cannot contribute to the generation. The integral scale L used in the definition of Re, Rm was assumed to be proportional to Λ. The sharp break in Eq. (26) results in a kink in σ(ρ) and in a discontinuity in σ′. This discontinuity produces a δ-function in the potential (9). This is a peculiar property of the model and just enhances and emphasizes the maximum of U(ρ), but the term ρσ′ does not play an essential role. It changes the resulting Rm_c to only about 20% (Kopyev et al. 2026). A much larger contribution to Rm_c is produced by the sharp kink at the same ρ = Λ.

From Eq. (25), X is the ratio of the end of the inertial range scale Λ and the magnetic diffusion scale r_d,

$$\begin{array}{c}\hfill X={\left(\frac{b_{\infty }}{2\eta }\right)}^{1/(1+s)}=\frac{\mathrm{\Lambda }}{r_{\mathrm{d}}}\text{where}{r}_{\mathrm{d}}:b(r_{\mathrm{d}})=2\eta .\end{array}$$

Kopyev et al. (2026) reported that the generation threshold corresponds to a critical value of X equal to

$$\begin{array}{c}\hfill X_c=20.5.\end{array}$$

Then, according to Eq. (24) and Eq. (12),

$$\begin{array}{c}\hfill {R{m}_c}_{(\mathrm{Sharp})}^{\mathrm{Sch}}=\frac{C_0}{4}{X}_c^{4/3}\simeq 100.\end{array}$$

In what follows, we compare this result with the results of other models.

4.2. Intermittency: Dependence on s

The scaling exponent of b inside the inertial range is an important parameter. It is composed of the exponents of ⟨(δv_∥)²⟩ and τ_c. Let ⟨|δv_∥|ⁿ⟩∝ρ^ζ_n ; then, according to the bridge relations (L’vov et al. 1997; Biferale et al. 2011), τ_c ∝ ρ^{1 − ζ₂ + ζ₁} and, hence,

$$\begin{array}{c}\hfill b(\rho )\propto {\rho }^{1+{\zeta }_1}.\end{array}$$(27)

In the Kolmogorov phenomenological theory, we obtain s = ζ₁ = 1/3. However, the intermittency of the velocity structure functions implies that ζ_n > n/3 for n < 3 (Frisch 1995).

The values of ζ_n depend on Re; for very high Reynolds numbers (Re_λ ≳ 650), we obtain s = ζ₁ = 0.39 (Benzi et al. 2010; Iyer et al. 2020). Even this small change in s affects the generation properties significantly: for s = 0.39, we obtain X_c = 14.6, which results in smaller critical Reynolds number,

$$\begin{array}{c}\hfill P{m}_c\simeq \frac{2160}{R{e}_{\lambda }^2},R{m}_c^{\mathrm{Sch}}\simeq 70.\end{array}$$

This qualitative behavior is an important indication for the explaination of the decrease in Rm_c for high enough Reynolds numbers. Iyer et al. (2020) showed that ζ₁ increases as a function of Re. This is a sufficient reason for the decrease in Rm_c observed by Warnecke et al. (2023).

4.3. Smooth model: Effect of the transition region

Not only the inertial range, but also the transition scales from the inertial to the integral range affect the generation of high Reynolds numbers. To demonstrate this, we considered the ’Smooth’ model: it coincides with the Sharp model in the inertial and integral ranges, but provides a more accurate description in between,

$$\begin{array}{c}\hfill b{(\rho )}_{\mathrm{Smooth}}=b_{\infty }\frac{{(\rho /\mathrm{\Lambda })}^{1+s}}{{(1+{(\rho /\mathrm{\Lambda })}^2)}^{\frac{1+s}{2}}}.\end{array}$$(28)

Since at ρ ≪ Λ and at ρ → ∞ this coincides with b_Sharp, the definition of r_d and the relation X = Λ/r_d remain the same. The shape of b(ρ) for the Smooth and Sharp models is illustrated in Fig. 1.

Fig. 1. Shape of b(ρ) for the Sharp (Eq. (26)) and Smooth (Eq. (28)) models for the same values of the parameters s, b∞, and Λ.

We considered the Smooth model with a classical nonintermittent scaling s = 1/3 and with a more realistic s = 0.39. The resulting X_c and Rm_c^Sch are presented in Table 3. Comparing the results of the theoretical models with the DNS data, we found that the agreement is quite reasonable. All the models show a value of Rm_c lower than the upper limit Rm_c^Sch ≤ 350 found by Schekochihin et al. (2007). The data of Schekochihin et al. (2007), Iskakov et al. (2007), Warnecke et al. (2023) demonstrate the decrease in Rm_c at Re_λ ≳ 300, and the value of Rm_c obtained at the highest Reynolds number Re ≃ 5 ⋅ 10⁴ achieved by Warnecke et al. (2023)² is Rm_c^Sch ≳ 200. This agrees very well with the prediction of the Smooth model for s = 0.39, which is the most realistic of our models.

Comparison of the Smooth and Sharp models shows that the effect of the transition region is significant: the critical values of the parameters differ by almost a factor of two. We determine the reason of the difference below.

The velocity correlator profiles and the corresponding potentials U(ρ) drawn at the critical values Λ_c are rather close in both models everywhere except in the vicinity of Λ_c(Sharp) (see Fig. 2). This appears to be the region in which the transition zone has the strongest effect. On the other hand, this is the boundary of the region in which the generation occurs. The generation threshold at any given Reynolds number corresponds to the zero energy level; Fig 2 shows that the scales ρ₀  :  U(ρ₀) = 0 for the Smooth models differ by less than 20% from the corresponding scales ρ_0(Sharp) = Λ_c(Sharp) for the Sharp models. Hence, the generation region ρ ≲ ρ₀ is situated for the Sharp and Smooth models within or near the scale Λ_c(Sharp) and rather far from the integral scale Λ_c(Smooth). Although the Sharp model is less plausible and accurate than the Smooth model, it indicates the scales that contribute most to the generation better.

Fig. 2. Effective potential U(ρ) (Eq. (9)) that corresponds to the generation threshold for the two Sharp and two Smooth models. The length scale is normalized by the diffusion scale rd, which is taken to be the same for both models. The vertical arrows correspond to the δ functions in the Sharp model potential.

The difference in the values of Rm_c is therefore probably caused by different meanings of the scale Λ in the two types of models. While Λ_c(Smooth) corresponds to the integral scale of turbulence, Λ_c(Sharp) merely acts as a transition range marker, which is about twice smaller.

The function b(ρ) contributes to U(ρ) by means of the logarithmic derivative σ(ρ). It is therefore interesting to compare these functions in critical regimes for different models. In Fig. 3 the functions σ(ρ) that correspond to the generation threshold are plotted. The rising parts of the curves are almost identical, the whole difference being in the descending part of the curve. In accordance with Kazantsev (1968), the scales important for the generation satisfy the condition σ > 0. For every s, the graphs for both models intersect the line σ = 0 at ρ ≃ Λ_c(Sharp). This shows that the Kazantsev criterion agrees well with the criterion ρ ≲ ρ₀.

Fig. 3. Logarithmic derivative σ(ρ), Eq. (10), for the four models. The length scale is normalized by the same rd. The parameters of each model correspond to the generation threshold.

The transition region probably determines the generation boundary. On the other hand, in the absence of the self-similarity of velocity correlators at the largest scales, the properties of this transition region may not be uniquely related to the velocity properties at the integral scales. If this is so, the comparison of some characteristics related to the transition region, for example, ρ₀ or the scale at which σ = 0 for different DNS-s or experiments, would be more effective than the comparison of large-scale parameters such as the magnetic Reynolds numbers.

In other words, it may be more effective to describe the generation threshold by means of some combination of parameters that describe the surrounding of the point ρ₀ or ρ : σ = 0 rather than by means of integral characteristics such as Rm_c.

5. Increment near the generation threshold

We considered the growth rate $\gamma =\underset{t\to \infty }{lim}dlnG/dt$ of the magnetic field correlator near the generation threshold. Namely, we were interested in

$$\begin{array}{c}\hfill g\equiv {Rm\frac{d\gamma }{\mathit{dRm}}|}_{Rm=R{m}_c}={\frac{d\gamma }{dln(Rm/R{m}_c)}|}_{Rm=R{m}_c}.\end{array}$$

We note that since the expression contains the ratio Rm/Rm_c , it does not depend on the choice of normalization.

Rogachevskii & Kleeorin (1997), Kleeorin & Rogachevskii (2012) predicted this value to be constant along the critical curve (Rm = Rm_c) for small Prandtl numbers. The data of Warnecke et al. (2023) confirmed this prediction for a wide range of Prandtl numbers beginning with Pm ≲ 0.1. The value of the coefficient found by Brandenburg et al. (2018), Warnecke et al. (2023) is

$$\begin{array}{c}\hfill {\left(\frac{{\gamma }_{\mathrm{DNS}-\mathrm{W}}}{ln(Rm/R{m}_c)}\right)}_{\mathit{Pm}}=0.022v_{\mathrm{rms}}{k}_f,\end{array}$$

where k_f is the pumping wave number; taking the relation ε/(v_rms³k_f) = 0.041 given by Brandenburg et al. (2018) into account, we can express this in terms of universal variables,

$$\begin{array}{c}\hfill {\left(\frac{{\gamma }_{\text{DNS-W}}}{ln(Rm/R{m}_c)}\right)}_{\mathit{Pm}}=0.18\frac{\epsilon }{u{\prime }^2}.\end{array}$$(29)

To compare this coefficient with the Kazantsev theory, we were again forced to take the value Re_λ = 140, since for this Reynolds number, we have information on the Lagrangian statistics (Biferale et al. 2011). The Prandtl number corresponding to this Re_λ is Pm ≈ 0.3 and does not belong to the range in which the increment slope is constant. However, according to the data of Warnecke et al. (2023), the value is close to the boundary and the slope does not differ significantly.

Making use of the velocity statistics given by Biferale et al. (2011), Donzis & Sreenivasan (2010), we solved Eq. (8) numerically for values of γ close to zero and obtained (for Re_λ = 140)

$$\begin{array}{c}\hfill {\left(\frac{{\gamma }_{\mathrm{th}}}{ln(Rm/R{m}_c)}\right)}_{\mathit{Re}}\simeq 1.3\epsilon /u{\prime }^2.\end{array}$$(30)

To compare this with Eq. (29), we note that the increment derivatives in Eq. (29) are taken at constant Prandtl numbers, while in Eq. (30), they are taken at constant Reynolds numbers. These two derivatives are related by

$$\begin{array}{c}\hfill \frac{{(\partial \gamma /\partial lnRm)}_{\mathit{Pm}}}{{(\partial \gamma /\partial lnRm)}_{\mathit{Re}}}=-\frac{\mathit{Re}}{P{m}_c}\frac{dP{m}_c}{\mathit{dRe}}=-\frac{dlnP{m}_c}{dlnRe}.\end{array}$$

This follows from Eq. (13) and from the condition that γ = 0 along the critical curve Rm_c(Re). From the data presented by Schekochihin et al. (2007), Warnecke et al. (2023), we found

$$\begin{array}{c}\hfill -\frac{dlnP{m}_c}{dlnRe}={(1-\frac{dlnR{m}_c}{dlnPm})}^{-1}\simeq 0.5.\end{array}$$

The results for the theoretical calculation and for the DNS results are summarized in Table 5. We also added g obtained from the DNS data (Schekochihin et al. 2007) (see Table 4),

$$\begin{array}{c}\hfill {\left(\frac{{\gamma }_{\text{DNS-Sch}}}{ln(Rm/R{m}_c)}\right)}_{\mathit{Re}}\simeq \frac{0.49+0.22}{ln(220/110)}\simeq 0.7\frac{\epsilon }{u{\prime }^2}.\end{array}$$(31)

The results differ significantly. Even the results of the two numerical simulations are quite different; when this is taken into account, the fact that the theoretical and numerical results are on the same order is a good correspondence.

The significant difference between the theoretical and numerical results can be explained as follows, for example:

– The data for hydrodynamics are taken from different sources, and this affects the errors strongly. In particular, C₀ is difficult to determine and depends strongly on the isotropy of the flow and other factors. We used the value (22) found analytically (Sawford 1991) in the frame of a particular model. The difference given by various DNSs and experiments is up to 1.5 times in both directions (Lien & D’Asaro 2002; Uma-Vaideswaran & Yeung 2025).

– In Table 4 the estimate of the slope of a curve near the inflection point is taken by only two distant points, which may lead to an underestimation of the slope. In addition, these two points are themselves known to have some errors.

– The details of DNS, such as periodic boundary conditions, can affect the generation.

– Finally, the real slope may be smaller than the one predicted by the Kazantsev model because of the non-Gaussianity of the velocity flow. If this factor were taken into account, the slope would decrease (Kopyev et al. 2024).

In the case of high Reynolds numbers, we lack enough data for nonsimultaneous velocity statistics, and we restricted ourselves by dimensional considerations. For the Sharp and Smooth models, we obtained inside the inertial range of scales

$$\begin{array}{c}\hfill {b(\rho )|}_{r_{\nu }\ll \rho \ll L}={\rho }^{4/3}{\left(\frac{\rho }{\mathrm{\Lambda }}\right)}^{s-1/3}{\epsilon }^{1/3}\stackrel{\sim }{a}(R{e}_{\lambda }),\end{array}$$(32)

where $\stackrel{\sim }{a}$ is a dimensionless parameter. As Re → ∞, the function $\stackrel{\sim }{a}(R{e}_{\lambda })$ has a constant limit (Frisch 1995). Taking Eq. (23) into account, we obtain for the Sharp model from the condition b(Λ) = b_∞

$$\begin{array}{c}\hfill \mathrm{\Lambda }={\left(\frac{b_{\infty }}{{\epsilon }^{1/3}\stackrel{\sim }{a}}\right)}^{3/4}={\left(\frac{4}{C_0\stackrel{\sim }{a}}\right)}^{3/4}\frac{u{\prime }^3}{\epsilon }.\end{array}$$

For the Smooth model, the result differs no more than by a multiplier of ∼2. Then, the increment is

$$\begin{array}{c}\hfill g=\frac{c{X}_c^{s-1}}{1+s}\frac{2\eta }{r_{\mathrm{d}}^2}=\frac{c}{1+s}\frac{b_{\infty }}{{\mathrm{\Lambda }}^2}=\frac{c}{1+s}\sqrt{\frac{C_0}{4}}\frac{\epsilon }{u{\prime }^2}{\stackrel{\sim }{a}}^{3/2}.\end{array}$$

This confirms that the increment does not depend on the Reynolds number for very high Re. Kleeorin & Rogachevskii (2012) obtained this result for s = 1/3; we prove that the presence of intermittency practically does not affect it, so that the independence of g is preserved up to the values of Re at which the bottleneck effect is essential.

6. Discussion

We studied the magnetic field generation in a turbulent flow near the generation threshold. We compared theoretical predictions of the extended Kazantsev theory with results of DNS for two cases: the value of the local Reynolds number Re_λ = 140, which is distinguished by the presence of magnetohydrodynamical and Lagrangian hydrodynamical data, and the limiting case of very high Reynolds numbers. To calculate the time-integrated velocity structure function b(ρ) (see Eq. (4)), we used the quasi-Lagrangian statistics (in particular, the Lagrangian velocity correlation time). This provided a better concordance with the DNS results by some orders of magnitude than those found with the Eulerian statistics. The values of critical magnetic Reynolds and magnetic Prandtl numbers obtained from the theory and from DNS differ by no more than ∼10%.

An ideal way to compare the DNS with theory should be to find b(ρ) and Pm_c in the same DNS and to solve Eq. (8) with this particular b(ρ). Even when we neglect the term with σ′ (which is difficult to calculate numerically) in the potential, we would obtain a rather high accuracy of presumably better than 20%. In this sense, the Kazantsev equation is a type of bridge connecting quasi-Lagrangian turbulence to the parameters of the magnetic field generation.

In absence of such a combined DNS, we took the data of Biferale et al. (2011), Donzis & Sreenivasan (2010) for the velocity statistics for one particular Reynolds number and compared the result with the data of Schekochihin et al. (2007). The result is presented in Table 2, and the agreement is quite good.

One more difficulty for the comparison of the theory and the numerical/experimental data is the use of different normalizations by different authors. This creates ambiguity and the additional need to coordinate data. In order to avoid this, it would be better to provide the stability curve in terms of Pm_c as a function of Re_λ, since, unlike Rm and Re, they are independent of the pumping properties of a particular flow and are universal in the sense that they can be expressed in terms of the basic parameters ε, v_rms, ν, η. Another way is to use a more universal definition of the integral scale. We presented our results in terms of Re^Sch,  Rm^Sch that are defined by Eqs. (12) and (13).

The comparison of different models shows the importance of the transition region between the inertial and the integral ranges. While the bottleneck transition between the viscous and inertial ranges is only important for relatively small Reynolds numbers, this outer transition range was found to be essential for the magnetic field generation at all Reynolds numbers.

For extremely high Reynolds numbers, the piecewise power-law model (Sharp) represents the qualitative properties of the generation well, but lowers the generation threshold significantly; the Smooth model appears to be much more accurate and agreed well with DNS. For the velocity scaling exponent found by Iyer et al. (2020), the obtained critical magnetic Reynolds number fits the DNS results well.

The comparison of the Sharp and the Smooth models helped us to determine the region that is most important for the value of the generation threshold. It is about twice smaller than the parameter Λ of the Smooth model, which is comparable to the pumping scale. This scale belongs to the transition region from the inertial to the pumping range. It can be marked by the condition σ(ρ)≈dlnb/dlnρ = 0, or by the claim that the effective potential (9) is zero at the point.

Based on the comparison of our data for different scaling exponents, we propose an explanation of the decrease in the critical magnetic Reynolds number as a function of Re for 650 ≳ Re_λ ≳ 300. The earlier explanation referred to the effect of the bottleneck, but the downward tendency of Rm_c also continues for Pm ≲ 0.01, where the effect of the bottleneck is certainly negligible. We note that the critical magnetic Reynolds number depends essentially on the scaling exponent of the velocity structure function inside the inertial range. As shown by Iyer et al. (2020), the scaling exponent decreases as a function of Re from one-third at relatively small Re to 0.39 as Re → ∞. This decrease is enough to ensure the decrease in Rm_c.

We also studied the dependence of the growth rate exponent of the magnetic field correlator on the magnetic Reynolds number. The results are presented in Table 5. The theoretical prediction differs by two to four times from the numerical results, which in turn differ significantly from each other. On one hand, this is a good correspondence taking into account the technical uncertainties. On the other hand, the difference might be decreased when the non-Gaussianity of the velocity field is taken into account. As in the case of Rm_c, the corresponding correction would change the theoretical result in the right direction.

To conclude, we stress the importance of a precise comparison of the kinematic dynamo theory with DNS. It would quantitatively verify the key theoretical simplifications, such as a δ-correlated in time effective velocity field. We obtained encouraging results, which confirm this assumption, by an attempt of such a comparison using the existing numerical and theoretical data on turbulence. We hope that future DNSs will be able to combine kinematic dynamo and Lagrangian velocity properties, which would allow for an accurate comparison with the Kazantsev theory. The precise comparison would also allow for a quantitative assessment of the effect of non-Gaussian effects, such as dynamo suppression and changes in the magnetic energy spectrum. A successful confirmation of the theory at moderately low magnetic Prandtl numbers would provide a critical basis for the extrapolation of its predictions to extreme, astrophysically significant Pm values, which are beyond our current computational capabilities.

Acknowledgments

This work of AVK was supported by the Foundation for the Advancement of Theoretical Physics and Mathematics (BASIS).

References

All Tables

All Figures

	Fig. 1. Shape of b(ρ) for the Sharp (Eq. (26)) and Smooth (Eq. (28)) models for the same values of the parameters s, b∞, and Λ.
In the text

	Fig. 2. Effective potential U(ρ) (Eq. (9)) that corresponds to the generation threshold for the two Sharp and two Smooth models. The length scale is normalized by the diffusion scale rd, which is taken to be the same for both models. The vertical arrows correspond to the δ functions in the Sharp model potential.
In the text

	Fig. 3. Logarithmic derivative σ(ρ), Eq. (10), for the four models. The length scale is normalized by the same rd. The parameters of each model correspond to the generation threshold.
In the text