© The Authors 2026

1. Introduction

The amplitude of the approximately 11-year solar cycle varies over a wide range, from extended periods of very low activity, such as the Maunder minimum, to episodes of unusually high activity, such as the modern grand maximum in the 1960s (Petrovay 2020; Usoskin 2023). This variability significantly affects the space environment near Earth and numerous aspects of human activity (Clark 2006), motivating efforts to predict solar cycle amplitudes. Since the solar cycle arises from a dynamo process operating within the Sun, prediction schemes based on solar dynamo models with data assimilation are particularly appealing. The first such attempts were made at the end of cycle 23, when two groups predicted cycle 24 using flux-transport dynamo models (Dikpati et al. 2006; Dikpati & Gilman 2006; Choudhuri et al. 2007; Jiang et al. 2007). These efforts subsequently prompted debate regarding the feasibility and robustness of solar cycle prediction based on dynamo model outputs (Tobias et al. 2006; Bushby & Tobias 2007). Nevertheless, research along this line has continued. In addition to flux-transport dynamo models (Jiang et al. 2013; Sanchez et al. 2014; Macario-Rojas et al. 2018), α − Ω mean-field dynamo models have also been applied to solar cycle prediction (Kitiashvili & Kosovichev 2008; Jouve et al. 2011; Kitiashvili 2020). This body of work raises the key physical question of identifying the necessary physical conditions for data assimilation of the surface magnetic field to result in predictive dynamo simulations.

The polar magnetic field constitutes a major component of the solar surface field, particularly during cycle minimum, and exhibits a strong correlation with the amplitude of the subsequent solar cycle. Schatten et al. (1978) first attempted to interpret the polar-field precursor within the framework of a Babcock–Leighton (BL) dynamo (Babcock 1961; Leighton 1969). In this class of models, the toroidal magnetic field is generated from a preexisting poloidal field through shearing by differential rotation (the Ω effect). The resulting toroidal flux subsequently emerges to form sunspots and is therefore expected to be correlated with the poloidal field present at the preceding cycle minimum, which is predominantly concentrated in the polar regions. Based on this dynamo interpretation, a number of subsequent studies have attributed the polar-field precursor relationship to the BL mechanism, regarding it as observational support for the operation of the solar dynamo in a BL regime (Schatten et al. 1978; Muñoz-Jaramillo et al. 2013).

However, Charbonneau & Barlet (2011) (hereafter CB11) suggested that such an interpretation of the polar-field precursor is not necessarily conclusive. They investigated three distinct stochastic-forcing dynamo models and found that the polar field exhibits predictive skill in both a mean-field α − Ω model with meridional circulation and a BL dynamo model, with nearly identical correlation coefficients in each case. Notably, the predictive skill of the surface polar field disappeared in the α − Ω model when meridional circulation was omitted. Thus, they concluded that the surface poloidal field only shows predictive skill if it is coupled back into the dynamo loop. At least part of the polar field should be submerged and carried into the tachocline to contribute to the induction of the toroidal magnetic component that will give rise to sunspots of the subsequent cycle. The condition can be satisfied in both BL and mean-field dynamo models.

Moreover, interpretations outlined above are based on dynamo models in which the dynamo operates primarily in the tachocline. In such dynamo models, both mean-field and BL types, the surface magnetic field is transported by meridional circulation into the tachocline. There, it feeds back into the dynamo loop by acting as the source to be amplified through the radial shear to generate the toroidal magnetic field (Pipin & Kosovichev 2011; Charbonneau 2020). However, in recently developed distributed-shear BL dynamo models (Zhang & Jiang 2022; Zhang et al. 2024; Jiang & Zhang 2025), toroidal flux is generated in the majority of the convection zone by latitudinal shear on the internal latitudinal poloidal field. Owing to the imposed condition of vanishing radial diffusion of the surface magnetic field (Zhang & Jiang 2022; Luo et al. 2026), the polar field cannot be submerged by the meridional flow, and the surface radial magnetic field does not directly participate in the internal dynamo process. Consequently, based on this dynamo model, a natural question arises: whether the assimilation of the surface magnetic field into the BL dynamo can retain predictive skill, and if so, what physical mechanism underlies it.

Beyond the unresolved questions, the methodologies commonly employed to assess the predictability in dynamo models also warrant reevaluation. The conventional approach typically involves two steps (e.g., Yeates et al. 2008; Charbonneau & Barlet 2011; Kumar et al. 2021): first, introducing stochastic forcing into a dynamo model and performing long-term simulations; and second, statistically analyzing correlations between the polar field strength at cycle minimum and the amplitude of the subsequent cycle inferred from the model output. This approach has several limitations when applied to the study of dynamo predictability. First, it primarily quantifies the correlation between the polar field at cycle minimum and the following cycle amplitude, rather than directly addressing whether the observed surface magnetic field actively provides the source information for the subsequent dynamo process. Second, the method requires long-term simulations for each dynamo model, which are computationally demanding and model dependent. Third, the inferred predictive skill sensitively depends on the level of stochasticity and nonlinearity prescribed in the model (Charbonneau & Barlet 2011; Kumar et al. 2021). In addition, because there is no unified parameter that directly quantifies the solar cycle strength in dynamo models, different proxies are adopted across studies (Yeates et al. 2008; Charbonneau & Barlet 2011; Guo et al. 2021; Ghosh et al. 2024), making consistent comparisons between models operating in different regions of the convection zone difficult. These limitations highlight the need for a more general framework that can be applied uniformly across dynamo models while directly probing the capability of predictive skill.

The lack of a unified parameter to quantify solar cycle strength in dynamo models largely reflects the open question of how emerging flux forms from the toroidal field (Fan 2021). Motivated by the observed hemispheric coherence of active-region orientations over a solar cycle (Hale et al. 1919), with only a small fraction violating Hale’s polarity law, Cameron & Schüssler (2015) (CS15) proposed that the hemispheric net toroidal flux, defined as the azimuthal average of the toroidal field, is the relevant quantity for sunspot formation. They further apply Stokes’ theorem to the induction equation. Thus the net toroidal flux can be calculated from a closed line integral bounding the cross-sectional area of each hemisphere. Using solar observations, CS15 concluded that the net toroidal magnetic flux is determined by the emerged surface flux and thus can be calculated from the observed surface magnetic field distribution. In particular, the dominant contribution to the net toroidal flux arises from the roughly dipolar polar field, providing a natural explanation for the polar-field precursor relationship. Finley et al. (2024) applied this approach to a 3D MHD dynamo simulation of a cool star, in which the net toroidal flux in each hemisphere is directly accessible, and verified the results of CS15 that the surface observables capture the toroidal flux generation inside the convection zone in their dynamo model. In contrast, Pipin & Kosovichev (2024) applied the same approach to their mean-field dynamo model (Pipin et al. 2023) and found that the net toroidal flux can only partially be estimated from the surface magnetic field distribution. These differing results suggest that the extent to which surface magnetic fields encode information about the net toroidal flux may depend sensitively on the underlying dynamo regime.

Although the methodology proposed by CS15 does not identify the detailed physical origin or the location of the toroidal flux generation, it provides a quantitative measure of the net toroidal flux. Thus the application of Stokes’ theorem to the induction equation provides a general framework for linking the surface observables to the internal dynamo process, and hence for assessing the predictability of dynamo models. Within this framework, the net toroidal magnetic flux integrated over the entire convection zone can be adopted as a unified proxy for the amplitude of the solar surface activity. Moreover, the extent to which this quantity is predominantly determined by the surface line integral associated with the observed surface magnetic field distribution can serve as a criterion for evaluating the predictive skill of a given dynamo model.

The paper is organized as follows. In Sect. 2, we present the general framework for analyzing solar cycle prediction. In Sect. 3, we introduce five dynamo models with different types. In Sect. 4, we apply the proposed framework to these dynamo models to assess their predictive skill and to elucidate the underlying physical mechanisms. Finally, Section 5 summarizes the main results.

2. The general framework for understanding solar cycle prediction

Since all dynamo models considered in the following sections are axisymmetric, we restrict our attention to the azimuthally averaged large-scale magnetic field. Its temporal evolution is governed by the mean-field induction equation (Moffatt 1978; Krause & Raedler 1980), which can be written as follows:

$$\begin{array}{c}\hfill {\partial }_t\overline{\mathit{B}}=\mathbf{\nabla }\times (\overline{\mathit{U}}\times \overline{\mathit{B}}+\overline{\mathcal{E}}),\end{array}$$(1)

where $\overline{\mathit{U}}$ and $\overline{\mathit{B}}$ are the azimuthally averaged large-scale flow field and magnetic field, respectively, and $\overline{\mathcal{E}}=⟨\mathit{u}\times \mathit{b}⟩$ is the mean electromotive force (EMF) representing the azimuthally averaged correlation between the turbulent-scale fluctuating magnetic field b and velocity field u. The term $\overline{\mathcal{E}}$ gives rise to turbulent magnetic diffusivity η, the α-effect, and possible turbulent pumping. The α-effect describes the regeneration of the poloidal field from the toroidal field and can arise from two distinct physical mechanisms, which are commonly used to distinguish different classes of dynamo models. In the α − Ω mean-field model, the α-effect originates from kinetic helicity in rotating, stratified convection. In contrast, in BL dynamos it results from the emergence, systematic tilt, and subsequent surface decay of sunspots.

We follow the approach of CS15 to calculate the net toroidal flux rate integrated in the whole convection zone by applying Stokes’ theorem to the induction equation (Eq. (1)). We define a contour 𝒞 that encloses the area 𝒮 in a meridional plane, as illustrated in Fig. 1. The normal vector of the surface 𝒮 is oriented parallel to the direction of the positive azimuthal field, B_ϕ. The time derivative of the net toroidal magnetic flux, corresponding to the flux generation rate, in the northern hemisphere of the Sun, Φ_tor^N, can be calculated as follows:

$$\begin{array}{cc}\hfill \frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}& =\frac{\mathrm{d}}{\mathrm{d}t}{\displaystyle {\int }_{\mathcal{S}}}{B}_{\varphi }\mathrm{d}S,\hfill \end{array}$$(2a)

$$\begin{array}{cc}& ={\displaystyle {\oint }_{\mathcal{C}}}(\overline{\mathit{U}}\times \overline{\mathit{B}}+\overline{\mathcal{E}})·d\mathit{l},\hfill \end{array}$$(2b)

Fig. 1. Schematic depiction of the application of Stokes’ theorem. The white and gray zones represent the convection zone and the radiative zone, respectively. The blue contour 𝒞 shows the path of integration around the area 𝒮 in the northern hemisphere. And the contour 𝒞 is clockwise and consists of the solar surface (a, surface), a radial segment in the equatorial plane (b, equator), a circular arc slightly below the bottom of the convection zone (c, base), a part along the axis of rotation (d, axis).

where dS is the surface element of 𝒮, and dl is the line element along 𝒞. The net toroidal flux of the southern hemisphere of the Sun can be calculated using a similar approach.

Following CS15, we calculate the right-hand side of Eq. (2b) in a reference frame corotating with the solar equator at an angular rate Ω^E. Note that ${\overline{U}}_{\varphi }=0$ along the rotation axis, $\overline{{\mathit{B}}_p}=0$ along the base line segment since the poloidal magnetic field does not penetrate the low-diffusivity radiative zone. The line integrals in Eq. (2b) can be expanded into the following form:

$$\begin{array}{c}\hfill \begin{array}{c}\hfill \frac{d{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathit{dt}}={\displaystyle \underset{I_1}{\underbrace{{\int }_0^{\pi /2}({\overline{U}}_{\varphi }-{\overline{U}}_{0\varphi }){\overline{B}}_r{R}_{\odot }{d\theta |}_{r=R_{\odot }}}}}+{\displaystyle \underset{I_2}{\underbrace{{\int }_{r_i}^{R_{\odot }}({\overline{U}}_{\varphi }-{\overline{U}}_{0\varphi }){\overline{B}}_{\theta }{dr|}_{\theta =\frac{\pi }{2}}}}}\\ \hfill +{\displaystyle \underset{I_3}{\underbrace{{\int }_{r_i}^{R_{\odot }}({\overline{\mathcal{E}}}_r{|}_{\theta =0}-{\overline{\mathcal{E}}}_r{|}_{\theta =\frac{\pi }{2}})dr+{\displaystyle {\int }_0^{\pi /2}}{\overline{\mathcal{E}}}_{\theta }{R}_{\odot }d\theta {{|}_{r=R_{\odot }}-{\displaystyle {\int }_0^{\pi /2}}{\overline{\mathcal{E}}}_{\theta }{r}_{i}d\theta |}_{r=r_i}}}},\end{array}\end{array}$$(3)

where R_⊙ denotes the solar radius, ${\overline{U}}_{0\varphi }=R_{\odot }sin\theta {\mathrm{\Omega }}^E$ is the background rotation speed, and r_i is the radial location of the base line segment c in Fig. 1. The value of r_i is chosen as the first interior grid point above the bottom boundary and depends on the specific dynamo model considered in the following calculations. The integrals I₁ and I₂ represent the line integrals of $\overline{\mathit{U}}\times \overline{\mathit{B}}$ evaluated along the surface and equatorial segments of the integration contour, respectively. They quantify the magnetic induction produced by differential rotation acting on the large-scale magnetic field. Hereafter, I₁ and I₂ are referred to as the surface and equatorial induction integral, respectively. In parallel, the integrals I₃ represent the contributions from the turbulent EMF. Specifically, when an isotropic scalarization parameter is employed, $\overline{\mathcal{E}}$ in Eq. (3) can be calculated by the following equation:

$$\begin{array}{cc}\hfill {\overline{\mathcal{E}}}_r& =\alpha \overline{B_r}-\eta \frac{1}{rsin\theta }\frac{\partial }{\partial \theta }(sin\theta {\overline{B}}_{\varphi }),\hfill \end{array}$$(4a)

$$\begin{array}{cc}\hfill {\overline{\mathcal{E}}}_{\theta }& =\alpha \overline{B_{\theta }}-\gamma {\overline{B}}_{\varphi }+\eta \frac{1}{r}\frac{\partial }{\partial r}\left(r{\overline{B}}_{\varphi }\right).\hfill \end{array}$$(4b)

It includes turbulent diffusion, pumping effect, and field-generation processes associated with the α-effect, which originates from kinetic helicity in α − Ω mean-field dynamo models and sunspot emergence and evolution in BL dynamos. Note that pumping is assumed to be vertical, which means the density pumping is neglected (Petrovay 1994). Hereafter, I₃ is referred to as the total EMF integrals.

In CS15, the primary objective was to estimate the net toroidal flux within a hemisphere using only observational data: the surface differential rotation and the magnetic field distribution. To achieve this, they introduced key approximations. Specifically, the equatorial induction integral I₂ along the equator was assumed to vanish, justified by helioseismology observations indicating that the angular velocity is nearly constant along the equator (Schou et al. 1998). Furthermore, the contribution of the total EMF integrals I₃ was simplified to consider only diffusive decay, represented as an exponential decay term with e-folding time τ. Under these assumptions, $\frac{d{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathit{dt}}$ is determined almost exclusively by I₁. Consequently, CS15 suggested that the surface induction integral I₁ captures the majority of the net toroidal flux generation, and since the polar field dominates I₁, this provides the physical basis for the polar-field precursor.

However, CS15 inferred from observations that surface magnetic fields dominate the net toroidal flux generation, but this may not hold across all dynamo regimes. As noted by Pipin & Kosovichev (2024), in some dynamo configurations the contributions I₂ and I₃ may reach magnitudes comparable to the surface induction integral I₁. The variability motivates us to apply the framework to assess predictability in different dynamo models. Dynamo simulations provide full access to the internal toroidal flux, allowing the net toroidal flux generation rate to be decomposed into the three contributions I₁, I₂, and I₃. Within this framework, the assessment of predictive skill can be formulated as a clear physical criterion: does the surface induction integral I₁ dominate the total generation rate of the net toroidal flux? We will apply this criterion to assess the predictive skill of the different dynamo models presented in the following section.

3. Dynamo models

In this section, we introduce the five dynamo models. The first three models employ parameter ranges similar to those of the models considered by CB11. The remaining two models are distributed-shear Babcock–Leighton (BL) dynamos recently developed by Zhang & Jiang (2022) and Jiang & Zhang (2025); one includes meridional circulation, while the other does not.

3.1. Model equations and boundary conditions

The axisymmetric large-scale magnetic field we investigate in Eq. (1) from Sect. 2 can be expressed in spherical coordinates as

$$\begin{array}{c}\hfill \overline{\mathit{B}}(r,\theta ,t)=B_{\varphi }(r,\theta ,t){\widehat{\mathbf{e}}}_{\varphi }+\mathrm{\nabla }\times \left[A(r,\theta ,t){\widehat{\mathbf{e}}}_{\varphi }\right],\end{array}$$(5)

where $B_{\varphi }(r,\theta ,t){\widehat{\mathbf{e}}}_{\varphi }$ represents the toroidal field and $\mathrm{\nabla }\times \left[A(r,\theta ,t){\widehat{\mathbf{e}}}_{\varphi }\right]={\mathit{B}}_p$ represents the poloidal field. For kinematic dynamo, the large-scale flow fields $\overline{\mathit{U}}$ are assumed to be steady and can be prescribed as

$$\begin{array}{c}\hfill \overline{\mathit{U}}(r,\theta )=rsin\theta \mathrm{\Omega }(r,\theta ){\widehat{\mathbf{e}}}_{\varphi }+{\mathit{U}}_p(r,\theta ),\end{array}$$(6)

where U_p(r, θ) represents meridional circulation and possible radial pumping, and Ω(r, θ) is the angular speed. The α − Ω type dynamo equations, which describe the spatiotemporal evolution of the poloidal and toroidal fields, are derived by expanding the mean-field induction equation (Eq. (1)) in the azimuthal direction and the meridional plane, and can be written in the following form:

$$\begin{array}{c}\hfill \frac{\partial A}{\partial t}+\frac{1}{s}({\mathit{U}}_p·\mathrm{\nabla })(sA)=\eta ({\mathrm{\nabla }}^2-\frac{1}{s^2})A+S(r,\theta ,B_{\varphi }),\\ \hfill \end{array}$$(7)

$$\begin{array}{c}\hfill \begin{array}{c}\hfill \frac{\partial B_{\varphi }}{\partial t}+\frac{1}{r}[\frac{\partial \left(r{U}_r{B}_{\varphi }\right)}{\partial r}+\frac{\partial \left(U_{\theta }{B}_{\varphi }\right)}{\partial \theta }]=\\ \hfill \eta ({\mathrm{\nabla }}^2-\frac{1}{s^2})B_{\varphi }+s({\mathit{B}}_p·\mathrm{\nabla }\mathrm{\Omega })+\frac{1}{r}\frac{d\eta }{\mathit{dr}}\frac{\partial (r{B}_{\varphi })}{\partial r},\end{array}\end{array}$$(8)

where s ≡ r sin θ, U_θ represents the latitudinal component of the meridional circulation, while U_r represents the radial component plus possible radial pumping γ, and η is the turbulent diffusivity. The source term S(r, θ, B_ϕ) = αB_ϕ describes the α-effect in mean-field dynamo models, while it represents the regeneration of the poloidal field at the solar surface in BL dynamo models. All parameters including Ω(r, θ), U_p(r, θ), η, and the source term S(r, θ, B_ϕ) for the five dynamo models presented in Sect. 3.2 are summarized in the upper part of Table 1.

The computational domain for mean-field models is 0.6 R_⊙ ≤ r ≤ R_⊙ and 0 ≤ θ ≤ π, whereas BL models adopt a bottom boundary at 0.65 R_⊙. Regarding the boundary conditions, both classes of models share a perfect conductor condition at their respective lower boundaries, requiring A = 0 and ∂(rB_ϕ)/∂r = 0; at the poles, A = B_ϕ = 0 for both (Charbonneau 2020). Two kinds of outer boundary conditions are considered. For the first three models described in Sect. 3.2, a potential-field condition with no electrical currents at r = R_⊙ (Dikpati & Choudhuri 1995) is imposed. For the remaining two models, we adopt a vertical-field condition, defined by ∂(rA)/∂r = 0 and B_ϕ = 0 at r = R_⊙, together with near-surface magnetic pumping. This setup of the outer boundary ensures that the resulting surface magnetic field evolution is consistent with surface flux transport models (Cameron et al. 2012; Yeates et al. 2023).

The initial magnetic field was set with B_p = 0, retaining only a latitudinally dependent B_ϕ, and was then relaxed to cyclic solutions for each dynamo model presented in Sect. 3.2. The first three models were calculated using the code SURYA developed by A.R.Choudhuri and his colleagues (Dikpati & Choudhuri 1994; Chatterjee et al. 2004), while the last two were solved with the Crank-Nicolson scheme combined with an approximate factorization technique (van der Houwen & Sommeijer 2001) developed at Beihang University.

3.2. Different dynamo models

3.2.1. Model 1: Flux transport BL model

Flux transport dynamo (FTD) models were the earliest dynamo models used to explain polar-field precursor and are widely recognized for their predictive capability. The model developed by Dikpati & Charbonneau (1999) is a representative example of this class. The parameters in the dynamo equation (Eqs. (7) and (8)), are adopted from the first case listed in their Table 1. Note that the original model yields symmetric solutions within one hemisphere, we extend the solution to a global domain (0 ≤ θ ≤ π) to incorporate possible cross-hemispheric flux exchange by I₃, while maintaining this north-south symmetry. Model 2 and Model 3 are also extended to obtain global solutions for the same reason, compared with the original model. This setup is essentially equivalent to Model 3 in CB11 except for the computational domain.

In the model, the surface poloidal field is proportional to the toroidal field strength at the same latitude at the tachocline, with a classical algebraic α-quenching as the nonlinear mechanism. The meridional circulation transports the surface poloidal flux toward the poles and subsequently downward to the tachocline, where it is stretched by differential rotation to generate the toroidal field for the subsequent cycle. Previous studies have demonstrated a correlation between the polar field at cycle minimum and the subsequent cycle strength by introducing stochastic fluctuations in the poloidal source terms (e.g., Yeates et al. 2008; Charbonneau & Barlet 2011). Applying the general framework described in Sect. 2 to this dynamo model therefore provides a useful test of its ability to assess the predictive skill of a given dynamo configuration.

3.2.2. Model 2: α − Ω mean-field model with meridional circulation

Our second model is an α − Ω mean-field model that includes meridional circulation. The α-effect is confined to the radial range [0.7 R_⊙, 0.8 R_⊙] and to a mid-latitude band (π/4 ≤ θ ≤ 3π/4), with algebraic α-quenching nonlinearity. The α-effect is co-spatial with the region of toroidal field generation. All parameters forms are identical to those of Model 2 in CB11 to enable direct comparison, but the values are chosen to obtain an 11-year cycle resembling the real Sun. Specifically, we take α₀ = 0.25 m s⁻¹, η₀ = 0.5 × 10¹¹ cm² s⁻¹, u₀ = 15 m s⁻¹, which means $C_{\alpha }=\frac{{\alpha }_0{R}_{\odot }}{{\eta }_0}=3.48$, $R_m=\frac{u_0{R}_{\odot }}{{\eta }_0}=209$ instead of the value C_α = 0.5 and R_m = 2500 adopted by CB11. The presence of meridional circulation allows the surface poloidal field to feed back into the toroidal field generation region, and this model therefore exhibits predictive skill in CB11. Applying our framework enables us not only to quantify the contribution of the surface magnetic field to the net toroidal flux generation, but also to compare it directly with the contribution arising from the α-effect.

3.2.3. Model 3: Classical α − Ω mean-field model

Our third model is a classical α − Ω mean-field dynamo operating without meridional circulation. We adopt the same parameters as Model 1 in CB11, except for a slightly decreased α₀ amplitude. Specifically, we take C_α = 4.65, compared to the value 5.0 used in CB11. In the absence of meridional circulation, there is no mechanism to feed the surface poloidal field back into the dynamo, and the model therefore does not exhibit a polar-field precursor, as demonstrated in CB11. This model allows us to investigate its behavior and contrast it with Model 2 within our general framework.

3.2.4. Model 4: Distributed-shear BL model

Previous three dynamo models are consistent with the ones considered by CB11. Recently Zhang & Jiang (2022) and Jiang & Zhang (2025) developed a new type of BL dynamo, distributed-shear BL model. The model differs fundamentally from the FTD models in that the subsurface meridional circulation in such a model plays a minimal role. While the surface radial field originates from the evolution of sunspots at surface as the BL type, it does not exert a back-reaction process by meridional circulation due to the vertical outer boundary condition and near-surface pumping. A central objective of including this model is to determine whether the surface field in such a BL dynamo retains its predictive skill for the subsequent solar cycle. We adopt the solution of reference case in Zhang & Jiang (2022).

3.2.5. Model 5: Distributed-shear BL model without meridional circulation

Meridional circulation in dynamo models not only couples the surface poloidal field back into the dynamo loop but also concentrate the surface magnetic flux toward the poles. This polar dominance has been central to previous interpretations of the polar-field precursor, since a polar-concentrated field substantially enhances the surface induction integral I₁ in Eq. (3). If the surface magnetic field distribution is no longer polar dominated, whether predictive skill is retained becomes a key question to address.

In FTD models, the meridional circulation is a necessary ingredient. In contrast, distributed-shear BL dynamos can operate without meridional circulation, as it is not required for sustaining cyclic surface magnetic field evolution (Leighton 1964, 1969; Baumann et al. 2004). This distinction motivates the construction of Model 5, which is a distributed-shear BL dynamo without meridional circulation. The source term is identical to that used in Model 3, with parameters set to α₀ = 0.2 m s⁻¹, η_CZ = 3.2 × 10¹⁰ cm² s⁻¹, η_s = 1.15 × 10¹² cm² s⁻¹. All other parameters are kept the same as in the reference case of Zhang & Jiang (2022).

4. Results

In this section, we apply the general framework presented in Sect. 2 to the five kinematic axisymmetric dynamo models to assess their predictive skill. Using the output of each dynamo model, we directly calculate the net flux generation rate within the northern hemisphere $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ based on Eq. (2a). We also directly calculate the surface and equatorial induction integral I₁ and I₂, and the total EMF integrals I₃ based on Eq. (3). By combining the evolution of the surface magnetic field distributions with the poloidal field configurations at solar minimum, we examine the relative magnitudes and phase relationships between these integrals and $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$. This analysis enables a physical interpretation of the results and helps elucidate the mechanisms underlying the models’ predictive skill. The results for the five dynamo models are summarized in the lower part of Table 1.

Figure 2 illustrates the results of Model 1, the flux transport BL model. Panel a shows the temporal evolution of the surface induction integral I₁, equatorial induction integral I₂, the total EMF integrals I₃, and the net flux generation rate within the northern hemisphere $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ calculated based on the parameters and output of the dynamo model. As expected, the sum I₁ + I₂ + I₃ consistently reproduces $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ in all models investigated below, validating our numerical implementation. For visual clarity and ease of comparison, we therefore display only $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ in this and the following figures. Panel (a) shows that I₁ dominates $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ and remains strongly phase-coherent with it throughout the time. Both I₁ and $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ peak around solar minimum, reaching amplitudes of approximately 3.26 × 10²³ Mx/yr, while the contribution of I₂ and I₃ can be negligible. This demonstrates that, in Model 1, toroidal flux generation is overwhelmingly controlled by the surface induction term, explaining the strong predictive skill of this dynamo configuration.

Fig. 2. Results of Model 1, a flux transport BL dynamo model (Dikpati & Charbonneau 1999). Panel (a): temporal evolution of the surface induction integral I1 (orange), equatorial induction integral I2 (green), the total EMF integrals I3 (red), and the net flux generation rate within the northern hemisphere $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ (blue). Panel (b): snapshot of the poloidal field configuration at cycle minimum. Panel (c): time–latitude diagram of the surface radial magnetic field Br. The vertical dashed lines in panels (a) and (c) mark the timing of cycle minimum, as shown in panel (b).

This result can be understood from the configuration of the poloidal field at cycle minimum, shown in Panel b, together with the surface magnetic field distribution in Panel c. According to Eq. (3), the equatorial induction integral I₂ represents the flux generation captured by B_θ crossing the equator, while I₃ arises from the EMF. During solar minimum, there is essentially no B_θ crossing the equator, rendering I₂ negligible. Meanwhile, $\overline{\mathcal{E}}$ is dominated by the diffusion term, whose contribution is strongly suppressed by the efficient meridional circulation in this flux-transport dynamo model. Consequently, the contributions of both I₂ and I₃ to the toroidal flux generation are very small: I₃ reaches a maximum amplitude of only 1.63 × 10²² Mx/yr and acts predominantly as a loss term, while I₂ remains close to zero. In contrast, the strong polar field evident in Panel c dominates the surface poloidal field distribution, causing the surface induction integral I₁ to control the toroidal flux generation rate.

This outcome is consistent with a main result of CS15. I₁ effectively traces the peak in toroidal flux generation. This peak, which occurs around cycle minimum, largely determines the total net toroidal flux available for the subsequent cycle and thus governs the ensuing solar cycle. Through this physical linkage, the surface magnetic field, dominated by the polar field near solar minimum, acts as an effective precursor for predicting the strength of the next cycle, providing a physical explanation for the predictive skill of this dynamo model.

Figure 3 illustrates the result of Model 2, an α − Ω mean-field model with meridional circulation. As shown in Panel (a), the overall behavior is similar to Model 1: I₁ (4.59 × 10²³ Mx/yr) continues to dominate the sum, accounting for the majority of the net toroidal flux generation rate $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ (7.00 × 10²³ Mx/yr) and maintains phase coherence with it throughout the time. While I₂ remains near to 0, I₃ in this model is non-negligible and provides a moderate positive contribution to the total toroidal flux generation, marking a departure from Model 1.

Fig. 3. Same as Fig. 2, but for Model 2, an α − Ω mean-field dynamo model with meridional circulation (Charbonneau & Barlet 2011).

This behavior can still be understood from the configuration of the poloidal field at cycle minimum (Panel b) and the surface magnetic field distribution (Panel c). As in Model 1, there is essentially no B_θ crossing the equator during solar minimum, and the differential rotation along the equator is weak, resulting in I₂ near to 0. However, in this α − Ω mean-field model with meridional circulation, $\overline{\mathcal{E}}$ contains both inductive contributions from the α-effect and turbulent diffusion. Unlike Model 1, where diffusion acts primarily as a loss term, the inductive component here provides a net positive contribution to toroidal flux generation. The presence of meridional circulation ensures that only a very small portion of the poloidal field across the equator and, more importantly, maintains the polar field’s dominance in the surface distribution. Consequently, $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ is still primarily captured by I₁. This preserves the close correspondence between the evolution of the polar field and the amplitude of the subsequent cycle.

Crucially, these findings demonstrate that even when the poloidal field generation mechanism shifts from the BL mechanism to the internal α-effect, the physical linkage established by applying Stokes’ theorem to the induction equation remains valid, and the surface field dominated by polar field retains its predictive skill for the following solar cycle strength. This implies that the existence of polar-field precursor is not a definitive indicator of the BL mechanism’s operation. The polar-field precursor phenomenon, in itself, shows no unique relation to the dynamo type and therefore cannot be used to conclusively constrain the dynamo type. This conclusion is consistent with the statistical inferences from long-term simulations in CB11 and the analysis by Kumar et al. (2021), both of which similarly caution against using the precursor correlation as a model discriminator.

Given that the poloidal field generation mechanism is not the decisive factor, we now seek to identify the fundamental physical driver of the predictive skill by the surface field. Classical α − Ω mean-field models were widely recognized as lacking predictive skill (Bushby & Tobias 2007; Sanchez et al. 2014; Bhowmik et al. 2023), which leads us subsequently investigate Model 3, a classical α − Ω mean-field model, to verify its predictive skill within the general framework. In contrast to Models 1 and 2, the performance of Model 3 shown by Fig. 4 reveals a fundamentally different regime. As shown by Panel a, the net toroidal flux generation rate $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ is found to be dominated not by I₁, but primarily by I₂ and I₃. Furthermore, I₁ loses its phase coherence with the flux generation rate, peaking instead around solar maximum. This shift is physically consistent with the poloidal field configuration at solar minimum (Panel b). In the absence of meridional circulation, the configuration possesses a significant B_θ crossing the equator and a strong poloidal field concentrated near the lower boundary below 0.7 R_⊙. The former significantly enhances I₂, while the latter implies that $\overline{\mathcal{E}}$ is dominated by turbulent induction associated with the α-effect. As a result, both I₂ and I₃ become substantial and can no longer be negligible. As these two terms exceed I₁, the evolution of the surface magnetic field becomes decoupled from the production of interior toroidal flux, explaining why the surface poloidal field in this model lacks predictive skill for the subsequent cycle (Pipin et al. 2023).

Fig. 4. Same as Fig. 2, but for Model 3, a classical α − Ω mean-field model (Charbonneau & Barlet 2011).

Thus far, our analysis of the first three models obtain the consistent result with prior statistical approaches (Charbonneau & Barlet 2011). This agreement both confirms that the polar-field precursor cannot be used to constrain the dynamo regime and demonstrates the effectiveness of the present framework. By replacing purely statistical correlations with a direct evaluation of the physical boundary integrals responsible for toroidal flux generation, the framework provides a physically grounded and robust measure of dynamo predictability.

In addition, Model 3 generates multiple bands of toroidal flux of alternating polarity within each hemisphere. It has been questioned whether the net toroidal flux, $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$, vanishes due to cancellation between these bands (Finley et al. 2024), and thus cannot serve as a proxy for surface flux emergence. However, we show that it retains a substantial value even in a classical α–Ω mean-field dynamo. This result further supports the view that the net toroidal flux provides a physically meaningful link between the interior toroidal field and the emergence of magnetic flux at the solar surface.

The performance of the first three models highlights meridional circulation as a pivotal factor for the predictive skill of the relevant dynamo models. Traditionally, the primary role of meridional circulation is thought to be the transport of surface poloidal fields to the tachocline, where they are coupled back into the dynamo loop and act as physical source via radial shear (Charbonneau & Barlet 2011; Kumar et al. 2021). In our framework, this process manifests as the suppression of I₂ and I₃ by minimizing cross-equatorial B_θ. However, meridional circulation simultaneously reshapes the surface field into a polar-dominated configuration, thereby significantly enhancing I₁. The effect of the spatial distribution of surface field has not been isolated in previous studies. Furthermore, as discussed in Sect. 3.2.4, emerging distributed-shear BL models rely on latitudinal shear throughout the convection zone to generate toroidal flux. In these regimes, even with meridional circulation, the surface field may not undergo a traditional feedback process due to the specific setup of the upper boundary. The predictive skill of such distributed-shear BL dynamos remains uninvestigated. To address these gaps, we employ Models 4 and 5, based on the distributed-shear BL framework developed by Zhang & Jiang (2022) and Jiang & Zhang (2025), to isolate the role of surface field distribution and evaluate predictive skill in this dynamo regime.

The behavior of Model 4 within the framework is illustrated in Fig. 5. Similar to Models 1 and 2, the surface induction integral I₁ accounts for the majority of the net toroidal flux generation rate $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ (5.36 × 10²² Mx/yr out of 6.12 × 10²² Mx/yr) and remains strongly phase-coherent with it throughout the cycle. In contrast to the earlier models, however, both I₂ and I₃ provide non-negligible contributions. These two terms operate largely in opposite phase and therefore partially cancel each other, such that $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ remains dominated by I₁. This cancellation preserves the phase coherence between I₁ and $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$. Panel b shows that the poloidal field configuration at solar minimum presents an amount of cross-equatorial B_θ, allowing I₂ to contribute positively to the net toroidal flux generation. At the same time, cross-hemispheric flux exchange acts as an efficient loss mechanism for toroidal flux (Cameron et al. 2018; Finley et al. 2024). And the turbulent EMF is dominated by the diffusive component, yielding a negative contribution from I₃.

Fig. 5. Same as Fig. 2, but for Model 4, a distributed-shear BL dynamo model (Zhang & Jiang 2022).

The above results suggest that the predictive linkage between the surface field and the subsequent cycle amplitude remains effective in the distributed-shear BL dynamo model. In particular, the polar field also dominates the surface poloidal distribution in this regime as shown in Panel c. Therefore, even in a model where the surface radial field does not directly serve as a source for the toroidal field and is not coupled back into the dynamo loop (due to the imposed vertical outer boundary condition), the surface field retains its role for evaluating the subsequent cycle’s strength.

Having established that the direct feedback of surface flux into the dynamo loop is not the decisive factor for the predictive skill of a dynamo model, we proceed to Model 5. By removing meridional circulation from Model 4, we aim to isolate and examine the specific influence of the surface magnetic field distribution. Unlike conventional tachocline-based models, the distributed-shear BL model explains the equatorward migration of toroidal flux through latitudinal dependence of the toroidal field generation (Jiang & Zhang 2025). Therefore, this model can regenerate the butterfly diagram without invoking the subsurface return flow of meridional circulation (Babcock 1961; Jiang & Zhang 2025).

The behavior of Model 5 is shown in Fig. 6. Compared with Model 4, the equatorial induction integral I₂ is significantly enhanced, while I₃ remains of comparable magnitude. This enhancement introduces a small phase offset between I₁ and the net toroidal flux generation rate $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$. Nevertheless, I₁ continues to account for the dominant fraction of $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ (9.26 × 10²² Mx/yr out of 1.10 × 10²³ Mx/yr) and remains largely phase-coherent with it. As shown Panel (b), the absence of meridional circulation allows the cross-equatorial B_θ to accumulate, leading to the increase in I₂. In contrast, the cross-hemispheric flux exchange and the associated diffusive loss captured by I₃ remain comparable to those in Model 4. Since I₂ and I₃ contribute with opposite signs, their net effect does not significantly alter the dominant contribution of I₁. This result demonstrates that predictive skill does not require the polar field to be the dominant component of the surface magnetic field, as illustrated in Panel c.

Fig. 6. Same as Fig. 2, but for Model 5, a distributed-shear BL dynamo model without meridional circulation (Zhang & Jiang 2022).

A contrasting behavior emerges when comparing Models 3 and 5. Despite exhibiting nearly identical surface magnetic field distributions characterized by mid-latitude maxima and dominant equatorward branches, Model 3 shows no predictive skill, whereas Model 5 retains it. The difference arises from the spatial relationship between the surface field and the dominant induction region. In the classical α − Ω mean-field Model 3, the absence of meridional circulation leads to a strong accumulation of poloidal field near the lower boundary, which substantially enhances I₃ through positive α-effect induction and reduces the relative contribution of I₁. As a result, the surface poloidal field becomes spatially decoupled from the deep radial shear layer where toroidal flux is generated, preventing the surface field from serving as an effective precursor. In contrast, Model 5 operates in the majority of the convection zone, and its surface magnetic field represents the radial component of the interior poloidal field. The surface field is therefore intrinsically linked to the shear process occurring immediately beneath it, allowing I₁ to capture the net toroidal flux generation rate. We thus conclude that, in addition to the flux-transport feedback operating in Models 1 and 2 as proposed by CB11, a second physical mechanism can also ensure the predictive skill of a dynamo. Specifically, when the surface field directly reflects the interior poloidal field, it remains physically coupled to the dynamo process and serves as a reliable precursor.

5. Discussion and conclusions

In this study, we apply Stokes’ theorem to the magnetic induction equation, and establish a direct physical link between the surface magnetic field and the subsequent dynamo process. This direct connection builds a general framework that provides a quantitative criterion for assessing dynamo predictability by examining the magnitude and phase relationship between the surface induction integral I₁ and the net toroidal flux generation rate $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$. Grounded in robust physical principles, the framework is then applied to investigate five representative dynamo models, helping to elucidate the mechanisms underlying solar cycle predictability.

The net azimuthally averaged toroidal flux generation rate is adopted as the proxy for the solar surface magnetic field in the paper. This is motivated by observational arguments, as proposed by Cameron & Schüssler (2015). This assumption has been questioned in the context of classical α–Ω mean-field dynamo models, which often produce multiple bands of toroidal flux of alternating polarity within each hemisphere. Finley et al. (2024) suggested that there would be zero net flux in each hemisphere for such models. However, in our Model 3 (a classical α − Ω mean-field dynamo), we still obtain a substantial net toroidal flux, even larger than in the other models considered in the paper. In addition, the existence of multiple toroidal-flux bands is not unique to classical α − Ω mean-field dynamo. Flux-transport B–L dynamos, such as Model 1 and Fig. 11 of Dikpati & Gilman (2006) also present similar multiband structures. The presence of a nonzero net toroidal flux across different types of dynamo models suggests that it provides a physically meaningful link between the interior toroidal field and the emergence of magnetic flux at the solar surface.

Using the hemispheric net toroidal flux integrated over the entire convection zone as a proxy for the surface magnetic field does not account for the uncertain location of toroidal field generation (e.g., confined to the lower layers or distributed throughout), which remains an open question. The flux-emergence process is also not yet fully understood. Recent BL dynamo models (e.g., Cameron & Schüssler 2020; Cloutier et al. 2025; Zhang & Chen 2026) have begun to incorporate simplified treatments of flux emergence, enabling a more direct link between the net toroidal flux and the surface field.

Our results show that predictive skill in the surface field, such as the polar-field precursor, is not a definitive indicator of the B-L mechanism and does not provide a conclusive constraint on the dynamo type. In addition, these results highlight that, for tachocline-based models to exhibit predictive skill, meridional circulation is crucial in coupling the surface magnetic field back into the dynamo loop. Model 3 (classical α − Ω mean-field model) exhibits a regime lacking predictive skill, where the net toroidal flux generation rate is dominated by the contributions from I₂ and I₃. This is attributed to B_θ crossing the equator and strong poloidal field concentrated near the lower boundary below 0.7 R_⊙, which significantly enhances I₂ and I₃, respectively. However, upon the introduction of meridional circulation in Model 2 (α − Ω mean-field model with meridional circulation), the flow suppresses this bottom concentration and reduces the cross-equatorial B_θ. Consequently, I₂ is significantly diminished and I₃ is partially weakened, thereby endowing the model with predictive skill. It is anticipated that as the meridional flow speed is further increased toward a fully advection-dominated regime, I₂ and I₃ would be further suppressed, eventually converging to the behavior of Model 1 (FTD). This physical transition aligns with the findings of Sanchez et al. (2014), who find that the predictive skill is shown to decrease with the strength of α-effect, and to increase with the strength of the imposed meridional circulation. This may also illustrate the nearly identical correlation coefficients between FTD and mean-field models with strong meridional circulation investigated by Charbonneau & Barlet (2011).

Our results reveal a novel predictive mechanism distinct from the traditional process of physically coupling the surface magnetic field back into the dynamo loop. Our results show that when the surface magnetic field represents the radial component of the interior poloidal field, the predictive skill also arises. In the distributed-shear dynamo model (Zhang & Jiang 2022), meridional circulation does not directly transport the surface magnetic field back into the internal dynamo loop. Remarkably, both models including Model 5, which operates without any meridional circulation, exhibit predictive skills. This comparative setup effectively rules out the necessity of a polar-dominated surface field distribution, a geometric factor often overlooked in previous studies. The underlying reason lies in the nature of these dynamos: the poloidal field manifests as a large-scale dipole structure, and toroidal flux generation relies primarily on latitudinal shear. In this configuration, I₁ captures $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ by this natural physical correlation, endowing the models with inherent predictive skills.

Contrary to previous assumptions (Cameron & Schüssler 2015; Finley et al. 2024), our results reveal that moderate cross-hemispheric flux exchange actually facilitates the estimation of $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ by I₁. Specifically, when a cross-equatorial B_θ is captured by I₂, it simultaneously interacts with the latitudinal shear to generate toroidal fields of opposite polarity across the hemispheres. This also creates a strong diffusive gradient, which is captured by I₃ as a negative contribution. Consequently, I₂ and I₃ tend to operate largely in opposite phase and therefore partially cancel each other. This cancellation preserves the phase coherence between I₁ and $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$. Given that the solar magnetic field at cycle minimum is dominated by large-scale antisymmetric axial components (Luo et al. 2024), such cross-hemispheric flux exchange likely plays a critical role in the high predictive skill observed in reality.

Acknowledgments

We thank the anonymous referee for valuable comments and suggestions that have improved the manuscript. The research is supported by the National Natural Science Foundation of China (grant Nos. 12425305, 12350004, and 12173005) and China’s Space Origins Exploration Program.

References

All Tables

All Figures

Fig. 1. Schematic depiction of the application of Stokes’ theorem. The white and gray zones represent the convection zone and the radiative zone, respectively. The blue contour 𝒞 shows the path of integration around the area 𝒮 in the northern hemisphere. And the contour 𝒞 is clockwise and consists of the solar surface (a, surface), a radial segment in the equatorial plane (b, equator), a circular arc slightly below the bottom of the convection zone (c, base), a part along the axis of rotation (d, axis).

In the text

Fig. 2. Results of Model 1, a flux transport BL dynamo model (Dikpati & Charbonneau 1999). Panel (a): temporal evolution of the surface induction integral I1 (orange), equatorial induction integral I2 (green), the total EMF integrals I3 (red), and the net flux generation rate within the northern hemisphere $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{tor}}^N}{\mathrm{d}t}$ (blue). Panel (b): snapshot of the poloidal field configuration at cycle minimum. Panel (c): time–latitude diagram of the surface radial magnetic field Br. The vertical dashed lines in panels (a) and (c) mark the timing of cycle minimum, as shown in panel (b).

In the text

	Fig. 3. Same as Fig. 2, but for Model 2, an α − Ω mean-field dynamo model with meridional circulation (Charbonneau & Barlet 2011).
In the text

	Fig. 4. Same as Fig. 2, but for Model 3, a classical α − Ω mean-field model (Charbonneau & Barlet 2011).
In the text

	Fig. 5. Same as Fig. 2, but for Model 4, a distributed-shear BL dynamo model (Zhang & Jiang 2022).
In the text

	Fig. 6. Same as Fig. 2, but for Model 5, a distributed-shear BL dynamo model without meridional circulation (Zhang & Jiang 2022).
In the text