© ESO, 2016

1. Introduction

Our understanding of the solar dynamo is based on theory, helioseismology, and the observational clues that we can derive from the emergence of the magnetic flux through the solar surface and its subsequent evolution. A well-studied theoretical model is the flux-transport dynamo (FTD) which explains the equatorward migration of the activity belts in terms of advective transport of toroidal flux; this model assumes a deep return flow of the observed near-surface poleward meridional flow. In order for this mechanism to work, the turbulent magnetic diffusivity affecting the toroidal field needs to be low enough to be overcome by the flow. Since simple estimates of the turbulent magnetic diffusivity in the convection zone based on the mixing length model give high values of the order of 1000 km² s^-1, most FTD models therefore assume that the toroidal flux is generated mainly by radial differential rotation associated with the tachocline and stored in a region of much lower magnetic diffusivity at or below the bottom of the convection zone (see reviews, e.g., Dikpati & Gilman 2009; Charbonneau 2010, 2014; Karak et al. 2014). Low diffusivity, possibly together with downward convective pumping (e.g. Jiang et al. 2013; Lawson et al. 2015), also prevents rapid loss of toroidal flux by radial diffusion (Zita 2010). Alternatively, it has been suggested that the required low values of the magnetic diffusivity affecting the toroidal field could be provided within the convection zone through diffusivity quenching by a sufficiently strong magnetic field (e.g., Muñoz-Jaramillo et al. 2011). In any case, the toroidal flux needs to be associated with a magnetic diffusivity that is orders of magnitude lower than the estimates for the convection zone based on mixing-length models.

Fig. 1 Time-latitude diagram (butterfly diagram) of the rate of sunspot group emergence based on the RGO/SOON data sets.

Fig. 2 Latitudinal cuts through the butterfly diagram shown in Fig. 1 for the years 1903 (during the weak cycle 14, solid red curve) and 1957 (during the strong cycle 19, solid black curve). In the left panel, the dashed (dash-dotted) curve represents a Gaussian fit to the number of spot groups in the southern (northern) hemisphere. In the right panel, the two hemispheres were fitted together with a single Gaussian.

In this paper, we use observations to infer the turbulent magnetic diffusivity, η_t, affecting the toroidal flux bands that provide the source of the magnetic flux emerging to form sunspot groups. We base our analysis on the properties of the time-latitude diagram of sunspot emergence (butterfly diagram) derived from the Royal Greenwich Observatory and Solar Optical Observing Network (RGO/SOON) sunspot records for solar cycles 12–23. We find that all cycles decline in the same manner once the centers of their activity belts have approached the solar equator to within their full width at half maximum (Sect. 2). Interpreting this result in terms of toroidal flux diffusing across the equator and cancelling with opposite-polarity flux, we determine a range of values for η_t that is within an order of magnitude with the mixing-length result for the convection zone (Sect. 3). This is in contrast to the requirement of FTD models that the toroidal field should be affected by a much lower magnetic diffusivity.

The high value of η_t would lead to severe spreading of the activity belts before they have propagated into the vicinity of the equator, contrary to what is observed. Such spreading could be prevented by a latitudinal inflow to the toroidal field bands. We estimate the structure and amplitude of such a flow in Sect. 4. It turns out that the properties of the inferred flow are similar to those of the observed shallow inflows toward the activity belts (Gizon et al. 2001; Gizon 2004; Zhao & Kosovichev 2004).

2. Properties of the butterfly diagram

Figure 1 shows the emergence rate of sunspot groups as a function of time and latitude, based on the RGO/SOON sunspot group records. We constructed this butterfly diagram by binning the sunspot groups in yearly and 1-degree bins, considering each group at the time of its maximum area. Similar to the analysis of Ivanov et al. (2011), we determined yearly values of the central latitude and the full width at half maximum of the butterfly wings, i.e., the migrating activity belts. We fitted Gaussians to the latitudinal profiles in both hemispheres using the curve_fit routine of NumPy (version 1.8.1). Examples of the original profiles and the Gaussian fits for the years 1903 (during the early phase of the weak cycle 17) and the year 1957 (during the early phase of the strong cycle 19) are shown in the left panel of Fig. 2. We also fitted the emergence rates in both hemispheres together (as a function of unsigned latitude); the corresponding fit is shown in the right-hand panel. In most of the following, we show the results for the fits based on the unsigned latitude (i.e., considering both hemispheres together) because the noise is lower. The main deviation from the Gaussian shape occurs near the equator, where fewer sunspots are observed than indicated by the Gaussian fits. Such an asymmetry is expected if the two activity belts (or rather the underlying oppositely directed toroidal flux systems) interact destructively and cancel each other across the equator.

We determine three parameters describing the evolution of the activity belt: amplitude A(t), central latitude λ_c(t), and its full width at half maximum, w(t). In addition, we consider the yearly international sunspot number, S(t) (version 2.0¹), all for cycles 12 through 23. Because cycles overlap in time by a few years, which is not accounted for in our Gaussian fitting, we follow Solanki et al. (2008) and consider the parameters only from years 2 to 9 during each cycle. Because we are interested in the relationship between activity and the latitudinal distribution (described by λ_c(t) and w(t)), we can choose to consider either S(t) and A(t) × w(t) as measure of the activity level. We tested both choices and found the results to be similar. Here we present the results based on S(t) because as a measure of activity it has been well studied and calibrated (Clette et al. 2014). In the following we therefore consider the relation between the three quantities λ_c(t), w(t), and S(t), which describe the evolution of the activity belts during the cycles covered by the data.

Figure 3 shows the sunspot number, S(t), as a function of the central latitude of the activity belt, λ_c(t), separately for cycles 12 to 23. This plot is similar to Fig. 6d of Ivanov & Miletsky (2014), except that we use the central latitude obtained from the fits for each year while they use an empirically derived model of the central latitude. This difference has no substantial impact. We see that each cycle begins with a different activity level at high latitude. As the activity belt propagates toward the equator, the sunspot number at first increases and then decreases sharply. The latitude at which the sunspot number begins to decrease depends on the strength of the cycle, with stronger cycles already beginning to decline when their activity belts are at higher latitudes. Because λ_c moves toward the equator in the same way for all cycles (Waldmeier 1955; Hathaway 2011), the decrease at high latitudes for stronger cycles is equivalent to a decrease earlier in the cycle than for weaker cycles. The fact that strong cycles begin to decline earlier is known as the Waldmeier effect. After they begin to decline, all the cycles show a similar sunspot number as a function of their mean unsigned latitude. The straight black line in Fig. 3 given by S = 25λ_c−200 gives a rough approximation of this common declining phase.

Fig. 3 Yearly International Sunspot Number S1 for cycles 12 to 23 as a function of the mean unsigned latitude of the sunspots each year. The line colors indicate the different cycles. The thin black curve is an approximate fit to the declining phase: S = 25λc−200. Error bars on the left side indicate standard deviations of λc from the Gaussian fitting procedure carried out for each year; the maximum value is shown in red while the mean over all years is given in black.

Next we consider the relationship between the full width at half maximum, w(t), and the central latitude, λ_c(t), of the activity belts for cycles 12 to 23, which is shown in Fig. 4. The belt width as a function of mean unsigned latitude behaves in analogy to the sunspot number shown in the previous figure. Stronger cycles begin with a high w that starts to decrease at higher latitudes than during weaker cycles. In their declining phases (at low latitudes), all cycles show a similar decrease of w. The common behavior in the declining phase can be approximated by the relation w = 1.1λ_c.

Fig. 4 Full width at half maximum (w) of the activity belts as a function of the mean unsigned latitude of the sunspots each year. Different colors represent different cycles as in Fig. 3. The thin black curve is an approximate fit to the declining phase: w = 1.1λc. Error bars on the top left indicate standard deviations of λc and w, respectively, from the Gaussian fit procedure carried out for each year; maximum values are shown in red while the mean values over all years are given in black.

An important point for the interpretation of the results is that the decline of activity once the sunspot belts have approached the equator to within their full width affects all latitudes within the belts in the same way; it is a universal decline, not just a decrease on its equatorward side. This is illustrated in Fig. 5, which shows the evolution of the activity belts for the strongest cycles (numbers 19 and 21) and the weakest cycles (numbers 14 and 16) in the data set. Relative to the year when λ_c was between 13° and 15° latitude (year zero, black lines), the plots show the latitudinal distribution of sunspot activity in the two preceding and up to five subsequent years. For the strong cycles, the decline begins at year zero while for the weak cycles the belts continue to migrate equatorward for another two years with no significant change of amplitude. The activity profiles for the following two years then show an overall decrease: activity declines at all latitudes of the activity belt while the shape stays roughly Gaussian.

Fig. 5 Latitudinal distribution of sunspots for several years during cycles 19 (upper left), 21 (upper right), 16 (middle and lower left), and 14 (middle and lower right). The black curves are for the years when the mean latitude of the sunspots was between 13° and 15°. The colors of the other curves indicate the year relative to this reference year (see legend in the top right panel). The curves were smoothed with a 3 degree boxcar filter. The late phases of cycles 14 and 16 are shown separately (bottom panels) to avoid overcrowding of lines.

At this point, we have seen that all of the cycles decline in a very similar way; the activity belts propagate toward the equator at the same rate and both their amplitude and width can be approximated by simple functions of their latitudinal distance from the equator. Cycles grow at different rates, with stronger cycles beginning to decline earlier (the Waldmeier effect). Once they begin to decline, however, cycles all decline in a similar way and the activity decreases over the whole width of the activity belts. We discuss the implications of these results in the next section.

3. Interpretation and estimate of η_t

The basic result obtained in the previous section is that all cycles decline in the same way. Figures 3–5 indicate that the activity belts propagate toward the equator showing increasing or approximately constant levels of activity until their distance from the equator becomes about equal to their full width at half maximum. Thereafter the activity strongly decreases across the whole latitude range covered by the belts. Stronger cycles begin with a larger width of their activity belts and thus start to decline earlier than weaker cycles. Once they begin to decline, all cycles follow the same pattern in terms of amplitude and latitudinal distribution. These properties suggest an interpretation of the declining phase in terms of cross-equator diffusion and cancellation of toroidal magnetic flux. We consider this possibility through a quantitative analysis, which provides an estimate of the magnitude of the relevant magnetic diffusivity.

For our analysis, we use radians as the unit for λ and w(t). We consider the distribution of sunspots as a function of latitude and time, N(λ,t), assuming that the activity belts in both hemispheres are symmetric with respect to the equator. We model this quantity as a Gaussian with full width at half maximum given by w(t), centered at λ_c and with a total area given by the sunspot number S/ 2, viz., $\begin{array}{ccc}\mathit{N}\mathrm{\left(}\mathit{\lambda ,t}\mathrm{\right)}\mathrm{=}\frac{\mathit{S}\mathrm{\left(}\mathit{t}\mathrm{\right)}}{\mathrm{2}\mathit{\sigma }\sqrt{\mathrm{2}\mathit{\pi }}}\exp \left(\mathrm{-}\frac{\mathrm{(}\mathit{\lambda }\mathrm{-}{\mathit{\lambda }}_{\mathrm{c}}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{2}{\mathit{\sigma }}^{\mathrm{2}}}\right)\mathit{,}& & \end{array}$(1)where $\mathit{\sigma }\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}\mathit{w}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathit{/}\sqrt{\mathrm{8}\ln \mathrm{2}}$.

We assume that the latitudinal profiles of the activity belts represent those of the underlying bands of toroidal magnetic flux from which the bipolar magnetic regions and sunspots emerge. This appears to be a reasonable assumption since the activity belts do not overlap and even at low latitudes it is uncommon to see sunspot groups that do not obey Hale’s polarity rules. If we assume further that the emerging flux (represented by the sunspot number) is proportional to the radially integrated toroidal flux density at each latitude, we obtain $\begin{array}{ccc}{\mathrm{\int }}_{\mathrm{0.7}{\mathit{R}}_{\mathrm{\odot }}}^{{\mathit{R}}_{\mathrm{\odot }}}\mathit{r}{\mathit{B}}_{\mathit{\phi }}\mathrm{\left(}\mathit{r,\lambda ,t}\mathrm{\right)}\mathrm{d}\mathit{r}\mathrm{=}\mathit{\alpha N}\mathrm{\left(}\mathit{\lambda ,t}\mathrm{\right)}\mathit{,}& & \end{array}$(2)where B_φ(r,λ,t) is the longitudinally averaged azimuthal component of the magnetic field, α is a constant of proportionality (which drops out of the final result), and R_⊙ is the solar radius.

We interpret the observational results illustrated in Figs. 3−5 in terms of diffusive cancellation of the oppositely oriented toroidal flux bands across the equator. This interpretation is motivated by the above result that the decline in activity each cycle begins only when the activity belts come within their full width of the equator, but from then on affects the whole activity belt. This behavior suggests that the loss of toroidal flux across the equator dominates over the losses due to cancellation at the axis of rotation as well as downward pumping of flux toward stable layers near the bottom of the convection zone. We consider the rate, D_eq, at which toroidal flux in the northern hemisphere is transported across the equatorial plane by diffusion. Applying Stokes’ theorem to the integrated toroidal flux over a hemisphere yields $\begin{array}{ccc}{\mathit{D}}_{\mathrm{eq}}\mathrm{=}{\mathrm{\int }}_{\mathrm{0.7}{\mathit{R}}_{\mathrm{\odot }}}^{{\mathit{R}}_{\mathrm{\odot }}}{\frac{{\mathit{\eta }}_{\mathrm{t}}\mathrm{\left(}\mathit{r}\mathrm{\right)}}{\mathit{r}}\frac{\mathit{\partial }{\mathit{B}}_{\mathit{\phi }}\mathrm{\left(}\mathit{r,\lambda ,t}\mathrm{\right)}}{\mathit{\partial \lambda }}|}_{\mathit{\lambda }\mathrm{=}\mathrm{0}}\mathrm{d}\mathit{r,}& & \end{array}$(3)where η_t(r) is the turbulent magnetic diffusivity. Since the depth dependence of the toroidal field (i.e., the storage location of the toroidal flux) is unknown, we rewrite Eq. (3) in the form $\begin{array}{ccc}{\mathit{D}}_{\mathrm{eq}}\mathrm{=}{\left(\frac{\mathit{\eta }}{{\mathit{R}}^{\mathrm{2}}}\right)}_{\mathrm{e}}\frac{\mathit{\partial }}{\mathit{\partial \lambda }}{\left({\mathrm{\int }}_{\mathrm{0.7}{\mathit{R}}_{\mathrm{\odot }}}^{{\mathit{R}}_{\mathrm{\odot }}}\mathit{r}{\mathit{B}}_{\mathit{\phi }}\mathrm{\left(}\mathit{r,\lambda ,t}\mathrm{\right)}\mathrm{d}\mathit{r})\right|}_{\mathit{\lambda }\mathrm{=}\mathrm{0}}\mathit{,}& & \end{array}$(4)where $\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}^{-1}$ is an effective diffusion time of the toroidal flux. If the flux is stored in a narrow radial layer around r = R₀, then $\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}\mathrm{=}{\mathit{\eta }}_{\mathrm{t}}\mathrm{(}{\mathit{R}}_{\mathrm{0}}\mathrm{)}\mathit{/}{\mathit{R}}_{\mathrm{0}}^{\mathrm{2}}$ and the relevant turbulent magnetic diffusivity affecting the toroidal field is given by η_t(R₀). For more general distributions, (η/R²)_e represents the effective quantity describing the diffusion of the toroidal flux. Using Eq. (2), we obtain $\begin{array}{ccc}{\mathit{D}}_{\mathrm{eq}}\mathrm{=}{\left(\frac{\mathit{\eta }}{{\mathit{R}}^{\mathrm{2}}}\right)}_{\mathrm{e}}\mathit{\alpha }{\frac{\mathit{\partial N}}{\mathit{\partial \lambda }}|}_{\mathit{\lambda }\mathrm{=}\mathrm{0}}& & \end{array}$(5)and substituting Eq. (1) we find $\begin{array}{ccc}{\mathit{D}}_{\mathrm{eq}}\mathrm{=}{\left(\frac{\mathit{\eta }}{{\mathit{R}}^{\mathrm{2}}}\right)}_{\mathrm{e}}\mathit{\alpha }\frac{{\mathit{\lambda }}_{\mathrm{c}}\mathit{S}}{\mathrm{2}{\mathit{\sigma }}^{\mathrm{3}}\sqrt{\mathrm{2}\mathit{\pi }}}\exp \frac{\mathrm{-}{\mathit{\lambda }}_{\mathrm{c}}^{\mathrm{2}}}{\mathrm{2}{\mathit{\sigma }}^{\mathrm{2}}}& & \end{array}$(6)for the rate at which the northern hemisphere loses toroidal flux by turbulent diffusion over the equator. Since opposite-polarity flux from the southern hemisphere diffuses into the northern hemisphere at the same rate, D_eq is actually half the rate at which the toroidal flux in the northern hemisphere decreases. Therefore we have $\begin{array}{ccc}{\mathit{D}}_{\mathrm{eq}}& \mathrm{=}& \frac{\mathrm{1}}{\mathrm{2}}\frac{\mathrm{d}}{\mathrm{d}\mathit{t}}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{\pi }\mathit{/}\mathrm{2}}{\mathrm{\int }}_{\mathrm{0.7}{\mathit{R}}_{\mathrm{\odot }}}^{{\mathit{R}}_{\mathrm{\odot }}}\mathit{r}{\mathit{B}}_{\mathit{\phi }}\mathrm{\left(}\mathit{r,\lambda ,t}\mathrm{\right)}\mathrm{d}\mathit{r}\mathrm{d}\mathit{\lambda }\\ & \mathrm{=}& \frac{\mathit{\alpha }}{\mathrm{2}}\frac{\mathrm{d}}{\mathrm{d}\mathit{t}}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{\pi }\mathit{/}\mathrm{2}}\mathit{N}\mathrm{\left(}\mathit{\lambda ,t}\mathrm{\right)}\mathrm{d}\mathit{\lambda }\\ & \mathrm{=}& \end{array}$(7)Equating Eqs. (6) and (7) yields $\begin{array}{ccc}\frac{\mathrm{d}\mathit{S}}{\mathrm{d}\mathit{t}}\mathrm{=}\mathrm{4}{\left(\frac{\mathit{\eta }}{{\mathit{R}}^{\mathrm{2}}}\right)}_{\mathrm{e}}\frac{{\mathit{\lambda }}_{\mathrm{c}}\mathit{S}}{\mathrm{2}{\mathit{\sigma }}^{\mathrm{3}}\sqrt{\mathrm{2}\mathit{\pi }}}\exp \frac{\mathrm{-}{\mathit{\lambda }}_{\mathrm{c}}^{\mathrm{2}}}{\mathrm{2}{\mathit{\sigma }}^{\mathrm{2}}}\mathrm{·}& & \end{array}$(8)To model the declining phase, we use the empirical relations obtained in Sect. 2, namely, $\begin{array}{ccc}{\mathit{\lambda }}_{\mathrm{c}}& \mathrm{=}& \mathrm{(}\mathit{S}\mathrm{+}\mathrm{200}\mathrm{)}\mathit{/}\mathrm{25}\mathrm{\times }\mathit{\pi }\mathit{/}\mathrm{180}\\ & \mathrm{=}& \mathrm{7}\mathrm{\times }{\mathrm{10}}^{-4}\mathrm{(}\mathit{S}\mathrm{+}\mathrm{200}\mathrm{)}\end{array}$(9)from Fig. 3 and $\begin{array}{ccc}\mathit{\sigma }& \mathrm{=}& \mathit{w}\mathit{/}\sqrt{\mathrm{8}\ln \mathrm{2}}\mathrm{=}\mathrm{0.47}{\mathit{\lambda }}_{\mathrm{c}}\\ & \mathrm{=}& \mathrm{3.3}\mathrm{\times }{\mathrm{10}}^{-4}\mathrm{(}\mathit{S}\mathrm{+}\mathrm{200}\mathrm{)}\end{array}$(10)from Fig. 4. Equation (8) then is a nonlinear ODE, which describes the evolution of the sunspot number during the declining phase of a cycle. We solve this ODE numerically. Figure 6 shows the observed decrease in S as a function of time (counted from activity minimum) for each cycle. The thick black curves show estimates for the evolution of S(t) based on Eq. (8) corresponding to two values of $\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}$: 450 km² s^-1 (solid black curve) and 300 km² s^-1 (dashed black curve). Since $\mathrm{1}\mathrm{\le }{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}\mathit{/}{\mathit{R}}^{\mathrm{2}}\lesssim \mathrm{2}$, we obtain values in the range of 150 km² s^-1 to 450 km² s^-1 for the effective turbulent diffusivity affecting the toroidal flux bands, where most of the uncertainty is due to the uncertainty in the depth at which the toroidal flux is located.

One may imagine alternative explanations of the observed evolution of the activity belts discussed in Sect. 2, which would allow the toroidal field to be stored in a low-diffusivity environment somewhere below the convection zone proper. One possibility would be magnetoshear instabilities in the tachocline (e.g., Cally et al. 2003; Gilman et al. 2007). However, these instabilities are not restricted to low latitudes and also would require a very peculiar nonlinear development to reproduce the uniform decline of activity throughout the whole activity belt and its onset depending on the width of the belt. The same problem is faced when assuming a field-dependent low-latitude upflow part of the meridional circulation; it needs to be of a very contrived field-dependent form to be able to reproduce the observed evolution of the activity belts, i.e., mimick a diffusion process.

Fig. 6 Sunspot number as a function of time of minimum activity. Different colors represent different cycles as in Fig. 3.

4. Inflows maintaining the activity belts

The high value for the turbulent magnetic diffusivity obtained here raises a problem. Our estimated range of diffusivities, 150–450 km² s^-1, corresponds to 5.7²–9.8² deg² yr^-1. If diffusion were the only process operating, then even with the lower diffusivity an initial Gaussian distribution with a FWHM of 14° degrees would increase its width by about 9.4° after only one year. This is inconsistent with the observed development of the activity belts during the rising phase of a cycle (cf. Fig. 4), which show an increase of their width only by roughly 2–3° per year. To restrict the spread of the wings to 3° or less would require the diffusivity to be less than 40 km²/s^-1. Higher values of the diffusivity therefore require the existence of an antidiffusive effect to counter spreading of the activity belts that is too strong. One possibility is to suppose that the toroidal field is affected by a latitudinal flow (inflow) toward the activity belt such that diffusive spreading is (partially) balanced by advection. Here we consider the simple kinematic problem with the aim of estimating the properties of an inflow that would balance the diffusive spreading.

To estimate the geometry and amplitude of such a flow field, we consider the induction equation for the azimuthally averaged magnetic field, B, $\begin{array}{ccc}\frac{\mathit{\partial }{B}}{\mathit{\partial t}}\mathrm{=}\mathrm{\nabla }\mathrm{\times }\mathrm{(}{U}\mathrm{\times }{B}\mathrm{)}\mathrm{-}\mathrm{\nabla }\mathrm{\times }{\mathit{\eta }}_{\mathrm{t}}\mathrm{\nabla }\mathrm{\times }{B}\mathit{.}& & \end{array}$(11)The rate at which the activity belts propagate toward the equator is about 1–2 degrees per year, while the expected diffusive spreading is substantially higher. Therefore we can make the approximation that the diffusion is balanced by the advective transport, i.e., $\begin{array}{ccc}\mathrm{\nabla }\mathrm{\times }\mathrm{(}{U}\mathrm{\times }{B}\mathrm{)}\mathrm{=}\mathrm{\nabla }\mathrm{\times }{\mathit{\eta }}_{\mathrm{t}}\mathrm{\nabla }\mathrm{\times }{B}\mathit{.}& & \end{array}$(12)Integration yields $\begin{array}{ccc}{U}\mathrm{\times }{B}\mathrm{=}{\mathit{\eta }}_{\mathrm{t}}\mathrm{\nabla }\mathrm{\times }{B}\mathit{,}& & \end{array}$(13)where we have set the arbitrary scalar potential to zero. The radial component of Eq. (13) reads $\begin{array}{ccc}{\mathit{U}}_{\mathit{\lambda }}{\mathit{B}}_{\mathit{\phi }}\mathrm{=}\frac{{\mathit{\eta }}_{\mathrm{t}}}{\mathit{r}\sin \mathrm{\left(}\mathit{\lambda }\mathrm{\right)}}\frac{\mathit{\partial }\sin \mathrm{\left(}\mathit{\lambda }\mathrm{\right)}{\mathit{B}}_{\mathit{\phi }}}{\mathit{\partial \lambda }}\mathrm{·}& & \end{array}$(14)On the left-hand side of this equation, we neglected the term −U_φB_λ, which describes the generation of toroidal field by differential rotation. This is justified since we show in Appendix A that the rate of flux change due to differential rotation in the low latitudes covers at most 10% of the required value derived from the observed butterfly diagram of sunspot emergence.

Similar to the determination of the effective diffusion time, $\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}^{-1}$, in Sect. 3 we derive an effective advection time, $\mathrm{(}{\mathit{U}}_{\mathit{\lambda }}\mathit{/}\mathit{R}{\mathrm{)}}_{\mathrm{e}}^{-1}$, of the toroidal flux owing to latitudinal inflows. Integrating Eq. (14) over radius and introducing the effective timescales we obtain $\begin{array}{ccc}\mathrm{(}{\mathit{U}}_{\mathit{\lambda }}\mathit{/}\mathit{R}{\mathrm{)}}_{\mathrm{e}}{\mathrm{\int }}_{\mathrm{0.7}{\mathit{R}}_{\mathrm{\odot }}}^{{\mathit{R}}_{\mathrm{\odot }}}\mathit{r}{\mathit{B}}_{\mathit{\phi }}\mathrm{(}\mathit{r,\lambda ,t}\mathrm{)}\mathrm{d}\mathit{r}& \mathrm{=}& \frac{\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}}{\sin \mathit{\lambda }}\\ & & \end{array}$(15)Using Eq. (2) then yields $\mathrm{(}{\mathit{U}}_{\mathit{\lambda }}\mathit{/}\mathit{R}{\mathrm{)}}_{\mathrm{e}}\mathrm{=}\frac{\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}}{\mathit{N}\mathrm{\left(}\mathit{\lambda ,t}\mathrm{\right)}\sin \mathit{\lambda }}\frac{\mathit{\partial }}{\mathit{\partial \lambda }}\left[\mathit{N}\mathrm{\left(}\mathit{\lambda ,t}\mathrm{\right)}\sin \mathit{\lambda }\mathrm{)}\right]\mathit{,}$(16)connecting the effective inflow with the time-latitude development of the sunspot emergence rate N(λ,t) shown in Fig. 1. If the toroidal flux is stored in a narrow radial layer around r = R₀, then (U_λ/R)_e = U_λ(R₀) /R₀ and the relevant inflow velocity is given by U_λ(R₀).

Fig. 7 Average activity belts based on aligning the butterfly wings shown in Fig. 1. The alignment was based on the year in which the average latitude of the wings is closest to 15 degrees. The red dots show the latitudes at which the activity belts are centered each year, and the red lines show the central latitude ± 12°.

We evaluated Eq. (16) on the basis of an avererage butterfly diagram to reduce the noise. All cycles in our data set were aligned in time with respect to the year during which the mean latitude of the sunspot groups in each wing is closest to 15°. The choice of 15° is somewhat arbitrary, but all cycles showed a reasonable number of sunspots near this latitude and the results are not strongly sensitive to this choice. The average butterfly diagram constructed in this way is shown in Fig. 7. We used the data for this average cycle and Eq. (16) together with Eq. (2) to obtain (U_λ/R)_eR_⊙ for an intermediate value of the magnetic diffusivity given by $\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}\mathrm{=}\mathrm{300}\hspace{0.17em}$km² s^-1. The result is shown in Fig. 8. Considering the ranges for $\mathrm{(}\mathit{\eta }\mathit{/}{\mathit{R}}^{\mathrm{2}}{\mathrm{)}}_{\mathrm{e}}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}\mathrm{=}\mathrm{150...450}\hspace{0.17em}$km² s^-1 and R/R_⊙ = 0.7...1, the corresponding amplitudes of U_λ become between 65% lower (≃ 3 m s^-1) and 50% higher (≃ 12 m s^-1).

We estimate U_λ under the assumption that advection exactly balances the outward diffusion. This is no longer valid near the equator in the declining phase of a cycle, when the activity belts approach each other. Since the inflows probably extend beyond the activity belts, the overlap between the flow systems of the two hemispheres reduces the inflows and thus allows for cross-equatorial diffusion of toroidal flux. The high velocities near the equator seen in Fig. 8 should therefore not be considered realistic.

The inflows derived here are similar in both extent and amplitude to the actually observed shallow inflows into the active region latitudes (Gizon et al. 2001; Gizon 2004; Zhao & Kosovichev 2004). Furthermore, Komm et al. (2015) showed that the inflows are part of the extended solar cycle, beginning at mid-latitudes (about 40 degrees). The inflows then propagate equatorward together with the activity belts. The extension in depth of these flows is unclear at the moment. Liang & Chou (2015) report that the helioseismic signature of the inflows is present down to 0.88–0.93 R_⊙, and partly present at depths from 0.75–0.88 R_⊙. On the other hand, Gizon & Rempel (2008) report a return flow at a depth of about 60 Mm (≃ 0.91 R_⊙). Such rather shallow inflows could result from cooling associated with the excess brightness of facular areas (Spruit 2003). Deep inflows could be driven by thermal perturbations caused by toroidal field amplification near the base of the convection zone (Rempel 2003). In fact, there are helioseismic indications of cycle-related changes of the sound speed in this region (Baldner & Basu 2008).

Since the inflows probably extend only over a limited depth range in the convection zone, we expect that part of the toroidal flux is spread over the convection zone by radial turbulent diffusion and an outward return flow. The effect of a deep return flow (in case of an inflow confined to the upper convection zone) is negligible since, owing to mass conservation and the strong increase of density with depth, its speed is much lower than that of the inflow. Radial turbulent diffusion cannot be neglected in this case, but the solution of a simple 2D advection-diffusion problem shows that a substantial part of the toroidal flux is still kept confined in the convergence zone of the inflows.

Fig. 8 Velocity required to prevent the activity belts from diffusing in latitude. The velocity is only shown in the latitudes within the red lines shown in Fig. 7.

5. Conclusion

Solar cycles begin to decline when the distance of the center of the activity belts from the equator becomes about equal to their (full) width (at half maximum). From this time onward, activity decreases across the whole latitude range covered by the belts and all of the cycles decline in the same way. Stronger cycles show wider activity belts and thus start to decline earlier than weaker cycles. This is the essence of the Waldmeier effect (Waldmeier 1955).

We interpret this result in terms of oppositely directed toroidal flux bands in each hemisphere that diffuse and cancel across the equator. The observational results are consistent with a diffusion process with effective value of the turbulent diffusivity acting on the toroidal magnetic field in the range 150 km² s^-1 to 450 km² s^-1. This is consistent (to order of magnitude) with the mixing length theory estimate, but inconsistent with the much lower values required by flux transport dynamos. Other interpretations of the observations in terms of magnetoshear instability or magnetic-field dependent upflows would require rather contrived properties of these processes, for which there is no evidence.

The flow field required to maintain the butterfly wings against the effects of diffusion is an inflow of a few meters per second, centered on the activity belts. Its spatial extent is at least 12° on both sides of the activity belt. Inflows with similar properties were observed in the upper part of the convection zone, but may extend also to deeper layers.

Acknowledgments

We acknowledge fruitful discussions with Matthias Rempel, Mark Cheung, Jie Jiang, and Emre Işik. NSO/Kitt Peak data used here are produced cooperatively by NSF/NOAO, NASA/GSFC, and NOAA/SEL. This work utilizes SOLIS data obtained by the NSO Integrated Synoptic Program (NISP), managed by the National Solar Observatory, which is operated by the Association of Universities for Research in Astronomy (AURA), Inc. under a cooperative agreement with the National Science Foundation. The synoptic magnetograms were downloaded from the NSO digital library, http://diglib.nso.edu/ftp.html. The sunspot numbers were obtained from WDC-SILSO, Royal Observatory of Belgium, Brussels, http://sidc.oma.be/silso. The RGO and SOON data sets were downloaded from http://solarscience.msfc.nasa.gov/greenwch.shtml.

References

Appendix A: Generation of toroidal flux by differential rotation

In order to estimate the effect of low-latitude generation of toroidal flux by differential rotation in the declining phase of an activity cycle, we consider the time derivative of the axisymmetric toroidal magnetic flux, ${{}^{\mathrm{\int }}\mathit{R}_{\mathrm{T}}}_{{\mathit{R}}_{\mathrm{\odot }}}^{}{\mathit{B}}_{\mathit{\phi }}\mathit{r}\mathrm{d}\mathit{r}$, radially integrated from the top of the tachocline, R_T, to the surface, R_⊙, as a function of latitude. This puts aside the question of where the toroidal flux is located in radius, but still requires us to make some assumptions on the internal structure of the magnetic field. Firstly, we assume that turbulent/convective pumping expels the horizontal components of the field from the near-surface shear layer (NSSL), so that the field there becomes purely radial. Secondly, we assume that the radial shear is zero between the NSSL and the top of the tachocline, and, thirdly, we assume that the magnetic field does not penetrate into the tachocline. Owing to the prolate shape of the tachocline (e.g., Antia & Basu 2011), the last assumption is somewhat questionable in higher latitudes, but here we are mainly interested in the flux generation in the low-latitude activity belts, where the main part of the tachocline is located below the convection zone. In addition, toroidal flux generation in the tachocline faces the problem of maintaining the required radial shear in the presence of a substantial Lorentz force and weak (or even absent) convection (Gilman & Rempel 2005; Spruit 2011).

In what follows we use a proper right-handed system of spherical coordinates, (r,θ,φ). The contribution of differential rotation, U_φ = Ω(r,θ) rsinθ, to the evolution of the azimuthally averaged toroidal field, B_φ, is given by $\begin{array}{ccc}\frac{\mathit{\partial }{\mathit{B}}_{\mathit{\phi }}}{\mathit{\partial t}}& \mathrm{=}& \frac{\mathrm{1}}{\mathit{r}}\left[\frac{\mathit{\partial }\mathrm{\left(}\mathit{r}{\mathit{U}}_{\mathit{\phi }}{\mathit{B}}_{\mathit{r}}\mathrm{\right)}}{\mathit{\partial r}}\mathrm{+}\frac{\mathit{\partial }\mathrm{\left(}{\mathit{U}}_{\mathit{\phi }}{\mathit{B}}_{\mathit{\theta }}\mathrm{\right)}}{\mathit{\partial \theta }}\right]\\ & \mathrm{=}& \frac{\mathrm{1}}{\mathit{r}}\left[\frac{\mathit{\partial }}{\mathit{\partial r}}\mathrm{\left(}\mathrm{\Omega }\hspace{0.17em}{\mathit{r}}^{\mathrm{2}}\sin \mathit{\theta }\hspace{0.17em}{\mathit{B}}_{\mathit{r}}\mathrm{\right)}\mathrm{+}\frac{\mathit{\partial }}{\mathit{\partial \theta }}\mathrm{\left(}\mathrm{\Omega }\hspace{0.17em}\mathit{r}\sin \mathit{\theta }\hspace{0.17em}{\mathit{B}}_{\mathit{\theta }}\mathrm{\right)}\right]\hspace{0.17em}\mathit{.}\end{array}$(A.1)Integrating over radius yields $\begin{array}{ccc}\frac{\mathit{\partial }}{\mathit{\partial t}}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}{\mathit{B}}_{\mathit{\phi }}\mathit{r}\hspace{0.17em}\mathrm{d}\mathit{r}& \mathrm{=}& \int \begin{array}{c}{\mathit{R}}_{\mathrm{\odot }}\\ {\mathit{R}}_{\mathrm{T}}\end{array}\sin \mathit{\theta }\hspace{0.17em}\frac{\mathit{\partial }\mathrm{\left(}\mathrm{\Omega }{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{\right)}}{\mathit{\partial r}}\hspace{0.17em}\mathrm{d}\mathit{r}\\ & & \mathrm{+}\frac{\mathit{\partial }}{\mathit{\partial \theta }}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{\Omega }\sin \mathit{\theta }\hspace{0.17em}{\mathit{B}}_{\mathit{\theta }}\mathit{r}\hspace{0.17em}\mathrm{d}\mathit{r}\mathit{.}\end{array}$(A.2)Since we assume B_r = 0 at r = R_T, we have ${\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}\sin \mathit{\theta }\hspace{0.17em}\frac{\mathit{\partial }\mathrm{\left(}\mathrm{\Omega }{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{\right)}}{\mathit{\partial r}}\mathrm{d}\mathit{r}\mathrm{=}{\mathrm{[}\sin \mathit{\theta }\hspace{0.17em}\mathrm{\Omega }{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{]}|}_{{\mathit{R}}_{\mathrm{\odot }}}\hspace{0.17em}\mathit{.}$(A.3)With our assumptions that B_θ = 0 in the NSSL and that Ω is independent of r between R_T and R_NSSL, we obtain $\begin{array}{ccc}\frac{\mathit{\partial }}{\mathit{\partial t}}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}{\mathit{B}}_{\mathit{\phi }}\mathit{r}\hspace{0.17em}\mathrm{d}\mathit{r}& \mathrm{=}& {\mathrm{[}\sin \mathit{\theta }\hspace{0.17em}\mathrm{\Omega }{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{]}|}_{{\mathit{R}}_{\mathrm{\odot }}}\\ & & \end{array}$(A.4)where R_NSSL indicates the bottom of the NSSL and we write Ω | _NSSL = Ω(R_NSSL,θ). The condition ∇·B = 0 for the azimuthally averaged field implies $\begin{array}{ccc}\mathit{r}\frac{\mathit{\partial }}{\mathit{\partial \theta }}\mathrm{\left(}{\mathit{B}}_{\mathit{\theta }}\sin \mathit{\theta }\mathrm{\right)}\mathrm{=}\mathrm{-}\sin \mathit{\theta }\frac{\mathit{\partial }}{\mathit{\partial r}}\mathrm{\left(}{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{\right)}\hspace{0.17em}\mathit{.}& & \end{array}$(A.5)Integration of this equation over radius yields $\begin{array}{ccc}\frac{\mathit{\partial }}{\mathit{\partial \theta }}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{NSSL}}}{\mathit{B}}_{\mathit{\theta }}\sin \mathit{\theta }\hspace{0.17em}\mathit{r}\mathrm{d}\mathit{r}& \mathrm{=}& \mathrm{-}\sin \mathit{\theta }{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}\frac{\mathit{\partial }\mathrm{\left(}{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{\right)}}{\mathit{\partial r}}\mathrm{d}\mathit{r}\\ & \mathrm{=}& {\mathrm{-}\sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}|}_{{\mathit{R}}_{\mathrm{\odot }}}\mathit{,}\end{array}$(A.6)where we used B_θ = 0 in the shear layer. By integrating over θ we obtain $\begin{array}{ccc}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{NSSL}}}{\mathit{B}}_{\mathit{\theta }}\sin \mathit{\theta }\hspace{0.17em}\mathit{r}\mathrm{d}\mathit{r}\mathrm{=}\mathrm{-}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{\theta }}{\sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}|}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{d}\mathit{\theta }\mathit{.}& & \end{array}$(A.7)We now substitute this relation into the last term of Eq. (A.4), using the assumptions that Ω is independent of r below the NSSL and B_θ = 0 within the NSSL, to obtain $\begin{array}{ccc}\frac{\mathit{\partial }}{\mathit{\partial t}}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}{\mathit{B}}_{\mathit{\phi }}\mathit{r}\mathrm{d}\mathit{r}& \mathrm{=}& {\mathrm{(}\sin \mathit{\theta }\hspace{0.17em}\mathrm{\Omega }{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{)}|}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{-}\frac{\mathit{\partial }}{\mathit{\partial \theta }}\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{NSSL}}}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{\theta }}\sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}{\mathrm{|}}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{d}\mathit{\theta }\\ & \mathrm{=}& {\mathrm{(}\sin \mathit{\theta }\hspace{0.17em}\mathrm{\Omega }{\mathit{r}}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{)}|}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{-}\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{NSSL}}}{\mathrm{(}\sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}\mathrm{)}|}_{{\mathit{R}}_{\mathrm{\odot }}}\\ & & \mathrm{-}\left(\frac{\mathit{\partial }\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{NSSL}}}}{\mathit{\partial \theta }}\right){\mathrm{\int }}_{\mathrm{0}}^{\mathit{\theta }}\sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}{\mathrm{|}}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{d}\mathit{\theta }\mathit{.}\end{array}$(A.8)Finally, we obtain $\begin{array}{ccc}\frac{\mathit{\partial }{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{T}}}^{{\mathit{R}}_{\mathrm{\odot }}}{\mathit{B}}_{\mathit{\phi }}\mathit{r}\mathrm{d}\mathit{r}}{\mathit{\partial t}}& \mathrm{=}& \sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}{\mathrm{|}}_{{\mathit{R}}_{\mathrm{\odot }}}\left(\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{-}\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{NSSL}}}\right)\\ & & \mathrm{-}\frac{\mathit{\partial }\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{NSSL}}}}{\mathit{\partial \theta }}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{\theta }}\sin \mathit{\theta }\hspace{0.17em}{\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}{\mathit{B}}_{\mathit{r}}{\mathrm{|}}_{{\mathit{R}}_{\mathrm{\odot }}}\mathrm{d}\mathit{\theta }\mathit{.}\end{array}$(A.9)The first term on the right-hand side of this equation corresponds to toroidal flux generation by radial differential rotation in the NSSL, while the second term describes the generation by latitudinal differential rotation. We know the value of both Ω | _RNSSL and Ω | _R⊙ from observations. For the surface rotation rate we use the synodic rate for magnetic fields of Hathaway & Rightmire (2011), $\begin{array}{ccc}\mathrm{\Omega }{\mathrm{|}}_{{\mathit{R}}_{\mathrm{\odot }}}& \mathrm{=}& \mathrm{14.30}\mathrm{-}\mathrm{1.98}{\cos }^{\mathrm{2}}\mathit{\theta }\mathrm{-}\mathrm{2.14}{\cos }^{\mathrm{4}}\mathit{\theta }{\mathrm{[}}^{\mathrm{◦}}\mathit{/}\mathrm{day}\mathrm{]}\mathit{.}\end{array}$(A.10)For the rotation rate at the base of the NSSL we add a latitude-independent value of 0.53°/day to the near-surface helioseismic result of Schou et al. (1998). This is motivated by the results of Barekat et al. (2014), indicating that the strength of the NSSL is independent of latitude.

The contributions to the generation of toroidal flux from both the radial and latitudinal shear can then be evaluated using azimuthally averaged KPNO/VTT and SOLIS synoptic maps; these are shown in Fig. A.1 (upper panels). The lower panels of this figure compare the combined flux generation by radial and latitudinal shear with the change of integrated toroidal flux as estimated from the sunspot emergence rate according to Eq. (2) in the main part of the paper. Here we chose the constant of proportionality, α, such that the result matches the total toroidal flux at activity maximum (≃ 5 × 10²³Mx/hemisphere) as estimated according to Cameron & Schüssler (2015) after integration in latitude. Considering the different scale ranges of the lower panels in Fig. A.1, it is clear that low-latitude generation of toroidal flux by differential rotation falls short of the required flux changes in the activity belts by at least a factor of 10. Incidentally, this result also indicates that the latitudinal propagation of the sunspot zone cannot be ascribed to a dynamo wave along the NSSL.

Fig. A.1 Generation rate of radially integrated toroidal magnetic flux as a function of time and latitude based on KPNO/VTT and SOLIS synoptic magnetograms. Upper left: flux generation by latitudinal shear. Upper right: flux generation by radial shear. Lower left: flux generation by latitudinal and radial shear combined. Lower right: change of integrated toroidal flux estimated from the sunspot emergence rate. We note the different scaling for the color bar of this plot. Black contours outline the region where sunspots were observed to emerge.

All Figures

	Fig. 1 Time-latitude diagram (butterfly diagram) of the rate of sunspot group emergence based on the RGO/SOON data sets.
In the text

Fig. 2 Latitudinal cuts through the butterfly diagram shown in Fig. 1 for the years 1903 (during the weak cycle 14, solid red curve) and 1957 (during the strong cycle 19, solid black curve). In the left panel, the dashed (dash-dotted) curve represents a Gaussian fit to the number of spot groups in the southern (northern) hemisphere. In the right panel, the two hemispheres were fitted together with a single Gaussian.

In the text

Fig. 3 Yearly International Sunspot Number S1 for cycles 12 to 23 as a function of the mean unsigned latitude of the sunspots each year. The line colors indicate the different cycles. The thin black curve is an approximate fit to the declining phase: S = 25λc−200. Error bars on the left side indicate standard deviations of λc from the Gaussian fitting procedure carried out for each year; the maximum value is shown in red while the mean over all years is given in black.

In the text

Fig. 4 Full width at half maximum (w) of the activity belts as a function of the mean unsigned latitude of the sunspots each year. Different colors represent different cycles as in Fig. 3. The thin black curve is an approximate fit to the declining phase: w = 1.1λc. Error bars on the top left indicate standard deviations of λc and w, respectively, from the Gaussian fit procedure carried out for each year; maximum values are shown in red while the mean values over all years are given in black.

In the text

Fig. 5 Latitudinal distribution of sunspots for several years during cycles 19 (upper left), 21 (upper right), 16 (middle and lower left), and 14 (middle and lower right). The black curves are for the years when the mean latitude of the sunspots was between 13° and 15°. The colors of the other curves indicate the year relative to this reference year (see legend in the top right panel). The curves were smoothed with a 3 degree boxcar filter. The late phases of cycles 14 and 16 are shown separately (bottom panels) to avoid overcrowding of lines.

In the text

	Fig. 6 Sunspot number as a function of time of minimum activity. Different colors represent different cycles as in Fig. 3.
In the text

	Fig. 7 Average activity belts based on aligning the butterfly wings shown in Fig. 1. The alignment was based on the year in which the average latitude of the wings is closest to 15 degrees. The red dots show the latitudes at which the activity belts are centered each year, and the red lines show the central latitude ± 12°.
In the text

	Fig. 8 Velocity required to prevent the activity belts from diffusing in latitude. The velocity is only shown in the latitudes within the red lines shown in Fig. 7.
In the text

Fig. A.1 Generation rate of radially integrated toroidal magnetic flux as a function of time and latitude based on KPNO/VTT and SOLIS synoptic magnetograms. Upper left: flux generation by latitudinal shear. Upper right: flux generation by radial shear. Lower left: flux generation by latitudinal and radial shear combined. Lower right: change of integrated toroidal flux estimated from the sunspot emergence rate. We note the different scaling for the color bar of this plot. Black contours outline the region where sunspots were observed to emerge.

In the text