SOLAR MAGNETIC FIELD REVERSALS AND THE ROLE OF DYNAMO FAMILIES

The solar dipolar magnetic field can be analyzed in terms of its axial and equatorial harmonic components. As has long been known (Hoeksema 1984), the axial dipole component, having a magnitude of |B⁰₁|, is observed to be largest during solar minimum when there is a significant amount of magnetic flux located at high heliographic latitudes on the Sun. These two so-called polar caps possess opposite polarity, and match the polarity of the trailing flux within active regions located in the corresponding hemisphere that emerged during the previous sunspot cycle. Long-term observations of surface-flux evolution indicate that a net residual amount of such trailing-polarity flux breaks off from decaying active regions and is released into the surrounding, mixed-polarity quiet-sun network. This flux is observed to continually evolve as flux elements merge, fragment, and move around in response to convective motions (Schrijver et al. 1997), but the long-term effect is that the net residual flux is slowly advected poleward by surface meridional flows. Such poleward advection results in a net influx of trailing-polarity flux into the higher latitudes. At the same time, an equivalent amount of leading-polarity flux from each hemisphere cancels across the equator, as is necessary to balance the trailing-polarity flux advected poleward. Over the course of a sunspot cycle, this process is repeated throughout subsequent sunspot cycles, during which flux from the trailing polarities of active regions eventually cancels out the polar-cap flux left over from previous cycles. Once the leftover flux has fully disappeared, the buildup of a new polar cap having the opposite polarity occurs by the subsequent activity minimum.

In contrast to the axial dipole component, the equatorial dipole components, having magnitudes |B⁻¹₁| and |B¹₁|, are largest during maximum activity intervals and weakest during activity minima. Individual active regions on the photosphere each contribute a small dipole moment that, aside from the small axial component arising from the Joy's law tilt, is oriented in the equatorial plane. Together the equatorial dipole moments from the collection of active regions add vectorially to form the overall dipole moment. When many active regions are on the disk, it thus follows that the equatorial dipole mode is likely to have a higher amplitude. During periods of quiet activity with few active regions on disk, the equatorial dipole amplitude is minimal.

Because the WSO data span three sunspot activity cycles, a bit of historical perspective on the evolution of the dipole can be gained, as shown in Figure 1. Figure 1(a) shows the amplitude of the axial dipole moment since mid 1976 and its evolution as a function of the activity level, represented in the figure by the sunspot number⁷ (SSN). It is also evident that, during the most recent minimum prior to Cycle 24, the magnitude of the axial dipole component is much lower than during any of the three previous minima (i.e., those preceding Cycles 21–23). The connection between the axial dipole component and the flux in the northern hemisphere is illustrated in the time history of the mean flux density as integrated over the northern hemisphere, shown in Figure 2. We observe a slight lag for the northern hemisphere averaged magnetic flux with respect to the axial dipole coefficient due to the contribution of other axisymmetric modes possessing a different phase.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 1(b) illustrates the magnitude of the equatorial dipole since mid 1976. In step with the relatively lower number of active regions during Cycle 23 when compared with the maxima for Cycles 21 and 22, the equatorial dipole strength is found to be lower during the most recent maximum than during the maxima corresponding to Cycles 21 or 22.

Given the variation in sunspot cycle strengths throughout the past few centuries, we suspect that cycle-to-cycle variations in the magnitudes of the axial and equatorial modes are not unusual. Proxies of the historical large-scale magnetic field, such as cosmic-ray-induced variations of isotopic abundances measured from ice-core data (Steinhilber et al. 2012), also show such longer-term variation and thus seem to be consistent with this view. Interestingly, the range over which the variation in the ratio of the energies possessed by the equatorial versus the axial dipole components is about the same for the three sunspot cycles observed by WSO, as shown in Figure 1(c). Longer-term measurements of this ratio unfortunately are not available due to the lack of a sequence of long-term magnetogram maps to which the harmonic decomposition analysis outlined in Section 2 can be applied.

The process by which old polar caps are canceled out and replaced with new, opposite-polarity polar caps, as described in the previous section, manifests itself as a change in sign of the axial dipole amplitude throughout the course of a sunspot cycle. Such dipole reversals for the past three sunspot cycles are shown in Figure 3, where the latitude and longitude of the dipole axis are plotted with time. It is found that the dipole axis spends much of its time in the polar regions, and for only about 12–18 months during these cycles it is located equatorward of ±45°.

Zoom In Zoom Out Reset image size

During these reversals, which occur at or near maximal activity intervals, the energy in the dipole never completely disappears. We find that the reduction in the energy (B⁰₁)² in the axial dipole is partially offset by an increase in the energy (B⁻¹₁)² + (B₁¹)² in the equatorial dipole. This results in a reduction of the total energy ∑_m(B^m₁)² in all dipolar modes only by about an order of magnitude from its axial-dipole-dominated value at solar minimum, as shown in Figure 4(a).

Zoom In Zoom Out Reset image size

Figure 3 indicates that, during a reversal when the axial dipolar component is weak, the axis of the equatorial dipolar component wanders in longitude. This seemingly aimless wandering occurs because the longitude of the dipole axis is primarily determined by an interplay among the strongest active regions on the photosphere at the time of observation. As older active regions decay and newer active regions emerge onto the photosphere, the equatorial-dipolar axis responds in kind.

The evolution of the energy contained in the quadrupolar (ℓ = 2) modes exhibit much more variation than the dipole, as shown in Figure 4. As with the equatorial dipole components, all of the quadrupolar modes have more power during greater activity intervals than during quieter periods, as illustrated in the evolution of the various quadrupolar modes plotted in Figure 5. Furthermore, when large amounts of activity occur, it is possible for the total energy ∑_m(B^m₂)² in all quadrupolar modes to be greater than the energy ∑_m(B^m₁)² in the dipolar modes at the photosphere. The ratio between these two groups of modes is shown in Figure 4(c), from which it is evident that during each of the past three sunspot cycles there have been periods of time when the quadrupolar energy exceeded the dipolar energy by as much as a factor of 10. The corona, in turn, reflects the relative strength of a strong quadrupolar configuration of photospheric magnetic fields by creating complex sectors and possibly multiple current sheets that extend into the heliosphere. One example of such complex field geometry is suggested by the potential-field source-surface model of Figure 6, where a quadrupolar configuration having an axis of symmetry lying almost in the equatorial plane is seen to predominate.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As with the quadrupole, the octupolar modes contain more power during periods of high activity and less power during minimum conditions, as illustrated in Figure 7. The exception is the axial octupolar coefficient B⁰₃, plotted in panel (a) of Figure 7, which is nonzero during solar minima and exhibits sign reversals during sunspot maxima in a manner similar to the axial dipole coefficient B⁰₁.

Zoom In Zoom Out Reset image size

The behavior of the various m = 0 modes can be understood by considering their functional symmetry: the Y⁰_ℓ functions are antisymmetric in θ (i.e., antisymmetric across the equator) when the degree ℓ is odd, whereas for even ℓ the Y⁰_ℓ functions are symmetric in θ. The presence of polar caps during solar minimum, a highly antisymmetric configuration, is reflected in the similar evolution of the B⁰₁ and B⁰₃ coefficients, which correspond to the axial dipole (ℓ = 1 and m = 0) and octupole (ℓ = 3 and m = 0) modes. The axial quadrupole (ℓ = 2 and m = 0) mode does not share this behavior because, as a symmetric mode, it is not sensitive to the presence of the polar caps during solar minima.

The dependence of the B⁰_ℓ coefficients on the degree ℓ is illustrated in Figure 8, where the time-averaged energies from the MDI data (spanning solar Cycle 23) as a function of degree ℓ are plotted. Prior to averaging, the spectra were placed in two classes: CRs for which the SSN is relatively large (defined as when SSN > 100) and rotations for which the SSN is relatively small (defined as when SSN < 50), thus capturing the state of the Sun when it is either overtly active or overtly quiet. The figure indicates that the even–odd behavior is more pronounced during quiet periods, and these occur near and during solar minimum when the polar-cap field is significant. During active periods the even–odd trend is still recognizable, but because the polar caps are weak and the active regions are primarily oriented east–west (i.e., in the equatorial plane and thus contributing little power to the axisymmetric modes) the even–odd trend is less pronounced. We will further discuss the behavior of axisymmetric modes in the context of BL dynamo models in Section 4.2.

Zoom In Zoom Out Reset image size

One property of the spherical harmonic functions Y^m_ℓ(θ, ϕ) is that the degree ℓ is equal to the number of node lines (i.e., contours in θ and ϕ, where Y^m_ℓ = 0). In other words, the spatial scale represented by any harmonic mode (i.e., the distance between neighboring node lines) is determined by its spherical harmonic degree ℓ. As a result, the range of ℓ values containing the greatest amount of energy indicates the dominant spatial scales of the magnetic field.

To this end, we have averaged the non-axisymmetric power spectra from each of the data sets both over time and over m, and have displayed the result in Figure 9. As with Figure 8, we have divided the spectra into active and quiet classes depending on SSN. In the figure, it can be seen that the magnetic power spectra form a broad peak with a maximum degree occurring at $\ell _{P_{\rm max}}\approx$ 25, corresponding to a size scale of about 360$^\circ /\ell _{P_{\rm max}}\approx$ 15° in heliographic coordinates. Stated another way, this indicates that much of the magnetic energy can be found (not surprisingly) on the spatial scales of solar active regions or their decay products.

Zoom In Zoom Out Reset image size

Energy spectra determined from WSO charts (not shown) do not show the same broad peak at $\ell _{P_{\rm max}}\approx$ 25 as found in the curves from the MDI-derived data shown in Figure 9. This is an effect of the lower spatial resolution of the WSO magnetograph (which has 180'' pixels and is stepped by 90'' in the east–west direction and 180'' in the north–south direction when constructing a magnetogram) versus MDI (which has a plate scale of 2'' pixel⁻¹ in full-disk mode). The WSO magnetograph, as a result, does not adequately resolve modes higher than about ℓ = 15, creating aliasing effects even at moderate values of ℓ in the energy spectra. Accordingly, as longer time series of data from newer, higher-resolution magnetograph instrumentation are assembled, the high-ℓ behavior of the energy spectra (such as those shown in Figures 8 and 9) may change due to better observations of finer scales of magnetic field.

The projection of the solar surface magnetic fields onto spherical harmonic degrees allows us to delineate the main symmetries of the magnetic field. As noted in Section 2, the harmonic modes can be classified as either axisymmetric (m = 0) or non-axisymmetric (m ≠ 0). Separately, the harmonic modes can be either antisymmetric (odd ℓ + m) or symmetric (even ℓ + m) with respect to the equator (Gubbins & Zhang 1993). Some authors refer to antisymmetric modes as dipolar and symmetric modes as quadrupolar (presumably because the axial dipole and quadrupole modes usually possess the most power in their respective categories), while others synonymously assign these modes to either the primary and secondary family (e.g., McFadden et al. 1991 when characterizing Earth's magnetic field geometry), respectively. In this article, we adopt the primary- and secondary-family nomenclature when describing the equatorial symmetry because this avoids the confusion that may otherwise occur when, for example, it is realized that the equatorial dipole mode (ℓ = 1 and m = 1) formally belongs to the "quadrupolar" family of modes (since ℓ + m is even for this mode).

One important result put forward by the geomagnetic community is that the relative strengths of the primary and secondary families are different during geomagnetic field reversals and excursions. During reversals, the modes associated with the secondary family predominate over primary-family modes, and during excursions this is not the case (Leonhardt & Fabian 2007). We now investigate whether analogous behavior is occurring in the solar setting, by determining which harmonic modes are most correlated with the axial dipole and axial quadrupole.

When applied to two variables, the Spearman rank correlation index ρ ∈ [ − 1, 1] indicates the degree to which two variables are monotonically related. The index ρ is positive when both variables tend to increase and decrease at the same points in time. A rank correlation analysis is more general than a more common Pearson correlation analysis, which specifically measures how well two variables are linearly related, whereas the rank correlation analysis enables a determination of whether the time evolution of two mode amplitudes follow a similar pattern in time without regard to their (unknown) functional dependency.

In Table 1, we list the degrees ℓ and orders m corresponding to the harmonic coefficients B^m_ℓ(t) that have the highest ρ (positive correlation) when compared with the axial dipole and axial quadrupole coefficients B⁰₁(t) and B⁰₂(t) (which peak at different phases of the sunspot cycle). The corresponding harmonic modes comprise the strongest modes in the primary and secondary families, respectively. We find that, among the mode amplitudes that are positively correlated with B⁰₁, two out of three belong to the primary family. Similarly, for B⁰₂, 7 of the 10 most-correlated modes are members of the secondary family.

These correlations indicate a preference in the solar dynamo, at least as inferred from its surface characteristics, for modes belonging to the same family and thus having the same north–south symmetry characteristics to be excited nearly in phase. This preference is demonstrated further in Figures 10 and 11, in which the long-term trends of the time evolution of the first several axisymmetric mode coefficients are shown, after smoothing with a boxcar filter having a width of 1 year. (We focus here on the axisymmetric mode properties because these modes are the only ones considered in most mean-field dynamo models, as discussed further in Section 4.2.) In Figure 10, there is a clear correlation among the first few odd-ℓ and among the first couple of even-ℓ mode coefficients, a trend which is emphasized in Figure 11 in which these same mode coefficients are overplotted. We note that the mode groupings are not precisely in phase, as evidenced, for example, by the lag in ℓ = 3 and especially the ℓ = 5 modes reversing signs with respect to the ℓ = 1 mode. When ℓ ⩾ 6 or so, these trends become much weaker among the axisymmetric modes (although Table 1 indicates that this is not necessarily true for the non-axisymmetric modes), presumably because as smaller and smaller scales are considered the effects of the global organization associated with the 11 year sunspot cycle are less important in structuring the surface magnetic field.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 11 additionally illustrates that the modes of the secondary family attain amplitudes of about 25% of the primary mode amplitudes. Furthermore, the primary and secondary mode families are out of phase: during reversals the primary modes become weak at the same time as the amplitudes of the modes associated with the secondary family become maximal, which was shown previously in Figure 4(c). This same pattern is observed to occur during reversals of the axial dipole field of the geodynamo. As in the geodynamo case, we ascribe the relative amplitudes and phase relation between the primary and secondary families observed during solar dipole reversals as a strong indication that the interplay of the mode families play a key role in the process by which the axial dipole reverses. Hence, any realistic model of the solar dynamo must excite both families of modes to similar amplitude levels, and must exhibit similar coupling between modes belonging to the primary and secondary families.

To illustrate such a dynamical system, following the work of Pétrélis & Fauve (2008), we assume that the axial dipole and quadrupole modes are nonlinearly coupled. We can then define a variable A(t) = B⁰₁(t) + iB₂⁰(t), where we have used the time-varying mode coefficients defined in Equation (1), and write an evolution equation that satisfies the symmetry invariance found in the induction equation, i.e., $\boldsymbol{B} \rightarrow -\boldsymbol{B}$. It then follows that the symmetry A → −A must also be satisfied, and that such an equation to leading order is

$$\begin{equation} \frac{dA}{dt} = \mu _1 A + \mu _2 \bar{A} + \nu _1 A^3 + \nu _2 A^2 \bar{A} + \nu _3 A \bar{A}^2 + \nu _4 \bar{A}^3, \end{equation} \tag{ 5 }$$

where μ_i and ν_i are complex coefficients and $\bar{A}=B_1^0 - i B_2^0$ is the complex conjugate of A, and the quadratic terms have vanished due to symmetry considerations. As discussed in Pétrélis & Fauve (2008) and Pétrélis et al. (2009), such dynamical systems are subject to bifurcations. In particular, they demonstrate that this dynamical system can be characterized by a saddle-node bifurcation when comparing its properties with so-called normal form equations (Guckenheimer & Holmes 1982). In such a bifurcated system, both stable and unstable equilibria (fixed points) exist, as illustrated in the left panel of Figure 12. For instance, if the solution lies at a stable point (for example, where the dipole axis is oriented northward), fluctuations in the system may disturb the equilibrium and push the magnetic axis away from its stable location. If such fluctuations are not strong enough, the evolution of the dynamical system resists the deterministic evolution of the system and the system returns to its original configuration (in the example, resulting in an excursion of the dipole), such as seen in the geomagnetic field (Leonhardt & Fabian 2007). If instead the fluctuations are large enough to push the system past the unstable point, the magnetic field then evolves toward the opposite stable fixed point allowed by $\boldsymbol{B} \rightarrow -\boldsymbol{B}$ (in the example, resulting in a reversal that changes the dipole axis to a southward orientation). Such behavior is also seen in the VKS experiment, from which is observed irregular magnetic activity combined with both excursions and reversals. Reversals result in an asymmetric temporal profile, with the dipole evolving slowly away from its equilibrium followed by a swift flip (cf., Figure 3 of Pétrélis et al. 2009).

Zoom In Zoom Out Reset image size

In Figure 13, we have overplotted the last three sunspot cycle reversals such that the zero crossings of the axial dipole coefficients B⁰₁ for each cycle are aligned. It has been recently shown that the 10 major geomagnetic reversals for which detailed records exist occurring during the past 180 Myr possess a characteristic shape upon suitable normalization (Valet & Fournier 2012). This shape can be described as comprising a precursory event lasting of order 2500 years, a quick reversal not exceeding 1000 years, and a rebound event of order 2500 years. Pétrélis et al. (2009) show that the magnetic field in a simplified VKS laboratory experiment exhibits differing behavior during reversals and excursions. During reversals, the magnetic field has an asymmetric profile that contains a slow decrease in the dipole, followed by a rapid change of polarity and buildup of the opposite polarity, whereas excursions are more symmetric. Additionally, after reversals the magnetic dipole overshoots its eventual value before settling down, whereas during excursions no such overshooting is measured (see Figure 3 of Pétrélis et al. 2009). For the solar cases displayed in Figure 13, we find that only the (green) curve of the reversal of Cycle 22 exhibits an overshoot, whereas the other two cycles do not. Further, the rates at which the solar dynamo approaches and recovers from the reversal appear to be equal, leading to a symmetric profile, in contrast with the VKS results. Therefore, the Sun seems to reverse its magnetic field in a less systematic way than other systems that have shown such behavior.

Zoom In Zoom Out Reset image size

Analyzing such systems from a dynamical systems perspective, when changing the control parameter (here, R_m) past the bifurcation point, the stable and unstable points coalesce and merge and the saddle nodes disappear, as shown in the right panel of Figure 12. This act transforms the system from one containing fixed points to one containing limit cycles with no equilibria (e.g., Guckenheimer & Holmes 1982), yielding an oscillatory solution that manifests itself as cyclic magnetic activity. Typically, large fluctuations are required in order to put the dynamical system above the saddle-node bifurcation threshold.

In the case of the Sun, both the primary and secondary families are excited efficiently and a strong coupling between them is exhibited. The model of Pétrélis & Fauve (2008), in spirit very close to the studies of Knobloch & Landsberg (1996) or Melbourne et al. (2001), may be used to guide our interpretation of the solar data. As illustrated in Figure 11, the axisymmetric dipole and quadrupole are out of phase, such that their coupling may lead to global reversals of the solar poloidal field. To the best of our knowledge, however, the solar dynamo does not exhibit excursions of its magnetic field (unlike the geodynamo) but instead undergoes fairly regular reversals that take about one or two years to transpire (cf., Figure 3). The solar dynamo is thus better approximated by a model in which a limit cycle is present. One may presume that the difference between the geodynamo and the solar dynamo may be a result of the large degree of turbulence present in the solar convection zone, whereas the mantle of the Earth has a more laminar convective flow and thus is below the bifurcation threshold where fixed points are still present.

It may be the case that the solar dynamo is better described by a Hopf bifurcation, in which a limit cycle arises (branches from a fixed point) as the bifurcation parameter is changed. The dynamo instability that occurs as a result of the interaction of magnetic fields and fluid flows (such as αω dynamos typically used to model the Sun, as summarized in Tobias 2002) often arises from a Hopf bifurcation. This allows the system to pass through domains having different properties, such as the aperiodic oscillations that characterize the grand minima and nonuniform sunspot cycle strengths of the solar dynamo (e.g., Spiegel 2009 and references therein). The data analysis shown here does not favor a particular type of bifurcation, but does indicate efficient coupling between the primary and secondary families.

Yet another approach toward investigating magnetic reversals is to develop detailed numerical simulations solving the full set of MHD equations. Such three-dimensional numerical simulations in spherical geometry of Earth's geodynamo (Glatzmaier & Roberts 1995; Li et al. 2002; Nishikawa & Kusano 2008; Olson et al. 2011) or of the solar global dynamo (Brun et al. 2004; Browning et al. 2006; Racine et al. 2011) have looked at the behavior of the polar dipole versus multipolar modes. Even though such models have large numerical resolution and thus possess a large number of modes, all have the property that the dominant polarity of the magnetic field follows the temporal evolution of a few low-order modes, even if in some cases the magnetic energy spectrum peaks at higher angular degree ℓ. These findings suggest that the coupling between the primary and secondary family remains an important factor in characterizing polarity reversals for these simulations and is thought to be linked to a symmetry breaking of the convective flow (Nishikawa & Kusano 2008; Olson et al. 2011). Some studies of the geomagnetic field (e.g., Clement 2004) even advocate for a coupling between two modes of the same primary family, such as the axial dipole and octupole. While in the solar data these modes are well correlated, the coupling between the primary and secondary families of modes seems more likely to be at the origin of the reversal, as demonstrated in Section 3.6.

Mean-field dynamo models are found to capture the essence of the large-scale solar dynamo (Moffatt 1978; Ossendrijver 2003; Dikpati et al. 2004; Charbonneau 2010). At present, the most favored model is the mean-field BL dynamo model (e.g., Dikpati et al. 2004), in which the mean magnetic induction equation is solved using empirical guidance for both the differential rotation and meridional circulation profiles, as well as for parameterizations of the α-effect and poloidal-field source terms. In this section, we use the Stellar Elements (STELEM; T. Emonet & P. Charbonneau 1998, private communication) code (see Appendix A of Jouve & Brun 2007 for more details) to solve the axisymmetric BL dynamo equations, and investigate some of the consequences of the coupling between modes from the primary and secondary families. In the interest of brevity, we refer interested readers to Appendix A.1 for a listing of the governing equations associated with BL dynamos.

In many BL solar dynamo models, the parameters governing the imposed flows and the poloidal-field source terms are chosen based on their solar counterparts. When carefully chosen, these terms favor modes from the primary family and thus antisymmetric field configurations, since this is what the Sun apparently favors much of the time. This is a result of the commonly used latitudinal profiles of the key dynamo ingredients (symmetric large-scale flows and antisymmetric alpha effect) combined with the parity in the BL mean-field dynamo equations, leading to a situation where modes of the primary family remain uncoupled to modes of the secondary family that allows both dynamo families to coexist without much interaction. We consider the symmetry of the BL equations used here in Appendix A.2 (see also Roberts & Stix 1972; Gubbins & Zhang 1993 for broader discussions on this topic).

To demonstrate these characteristics, we now consider a typical solution of the standard BL mean-field dynamo as calculated by STELEM. Figure 14 presents the time evolution of the resulting magnetic field patterns and is thus analogous to the standard solar butterfly diagram. Performing a Legendre transform on the magnetic field reveals the degrees ℓ of the dominant axisymmetric modes. Figure 15 illustrates that the odd ℓ modes from the primary family dominate over the even ones by about 5 orders of magnitude in this model. This differs significantly from what is observed on the Sun, where the amplitude of the quadrupole is measured to be about 25% of the dipole amplitude for most of the time, becoming dominant only during reversals (cf., Figure 4(c)). The behavior of the standard BL model of Figure 15 arises because of the symmetry characteristics of the BL dynamo equations. Because the model was initialized with a dipolar field, no modes from the secondary family are excited in the standard BL model shown in Figure 15 because no coupling exists between the primary and secondary families.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

If instead the calculation were initialized with a quadrupolar field (belonging to the secondary family), we find that the system eventually transitions to a state in which the primary-family modes predominate, as shown in Figure 16. The growth of the primary-family modes is due to the presence of a very weak dipole (likely of numerical origin) at the onset of the simulations. In these models, the BL source term of Equation (A8) quenches the growth of the magnetic field once a certain threshold is passed, and as a result the maximum total amplitude of the magnetic field is capped. The reason why the primary-family modes are preferred stems from the fact that the thresholds for dynamo action (based on the parameter C_s in Equation (A4)) are found empirically to be lower for the dipole than for the quadrupole (C_s ∼ 6.12 versus C_s ∼ 6.25), meaning the dipole-like modes have a higher growth rate than the quadrupole-like modes. In this model, only briefly during the transition phase does the model possess a quadrupole of order 25% of the dipole, as in the Sun.

Zoom In Zoom Out Reset image size

Observations of solar photospheric fields, however, indicate that the Sun excites both families and does not strongly favor members of one family over the other, a situation that has apparently existed over many centuries. Even during the Maunder minimum, evidence suggests that this interval may have been dominated by a hemispherical dynamo with magnetic activity located primarily in the southern hemisphere (Ribes & Nesme-Ribes 1993), which can only be formed by a state in which primary- and secondary-family modes possess nearly equal amplitudes (Tobias 1997; Gallet & Pétrélis 2009). Consequently, the solutions presented in Figures 15 and 16, in which modes from only one family are preferred, are thus not a satisfactory model of the Sun.

As advocated by Roberts & Stix (1972) and Gubbins & Zhang (1993) following their symmetry-based study of the solar dynamo and the induction equation, and more recently by Nishikawa & Kusano (2008) in their geodynamo simulations, a north–south asymmetry of the flow field specified in the BL dynamo, or alternatively an asymmetric poloidal-field source term, may allow the coexistence of both the primary and secondary families. To investigate this effect we have performed two additional BL dynamo calculations, one with a BL source term and one with a meridional circulation that each generate both symmetric and antisymmetric fields (introduced via the parameter of Equations (A9) and (A13) of Appendix A.1), the combination of which yields asymmetric magnetic field patterns.

We have run several dynamo cases, with the antisymmetry arising either in the BL source term or in the meridional flow profile, and with a range of amplitudes for the parameter from 10⁻⁴ to 10⁻¹. All cases were initialized with a dipolar field. We find that when is about 10⁻³, the modes in secondary family grow until they reach about 35% of the dominant dipolar mode, as illustrated in Figures 17 and 18. This result holds true regardless of whether the antisymmetry is introduced in the BL source term or in the meridional circulation profile, with very little difference in the resulting mode amplitudes. As expected, using a smaller results in solutions where the primary-family modes dominate, while using a larger yields a state where the secondary-family modes are comparable to the primary-family modes. Such results may indicate that the Sun need only possess a weak degree of north–south asymmetry in order to behave as it does.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Here, we briefly list the equations governing the axisymmetric mean-field dynamo models calculated in Section 4.2. A more detailed explanation can be found in, e.g., Jouve et al. (2008). Following Moffatt (1978), the mean-field induction equation is

$$\begin{eqnarray} \frac{\partial {\langle \boldsymbol{B} \rangle}} {\partial t} &=& \boldsymbol{\nabla }\boldsymbol{\times } (\langle \boldsymbol{V} \rangle\boldsymbol{\times }\langle\boldsymbol{B} \rangle)\nonumber\\ && + \boldsymbol{\nabla }\boldsymbol{\times } \langle\boldsymbol{v^{\prime }}\boldsymbol{\times }\boldsymbol{b^{\prime }}\rangle - \boldsymbol{\nabla }\boldsymbol{\times }\left(\eta \boldsymbol{\nabla }\boldsymbol{\times }\left< \boldsymbol{B}\right>\right), \end{eqnarray} \tag{ A1 }$$

where the variables $\left<\boldsymbol{B}\right>$ and $\left<\boldsymbol{V}\right>$ refer to the mean parts of the magnetic and velocity fields, and $\boldsymbol{v^{\prime }}$ and $\boldsymbol{b^{\prime }}$ to their respective fluctuating components. The function η is the magnetic diffusivity and is not necessarily a constant. The terms "mean" and "fluctuating" refer to the fact that a separation of scales has been performed, such that the mean quantities are computed by averaging over some appropriate intermediate size scale and the fluctuating quantities are the residuals.

Working in spherical coordinates (r, θ, ϕ) and under the assumption of axisymmetry, we perform a poloidal–toroidal decomposition and write the mean magnetic field B and mean velocity field V (for clarity the angle brackets 〈 and 〉 are omitted going forward) as

$$\begin{eqnarray} \boldsymbol{B}(r,\theta,t)&=\boldsymbol{\nabla }\boldsymbol{\times } (A_{\phi }\hat{\boldsymbol{e}}_{\phi })+B_{\phi } \hat{\boldsymbol{e}}_{\phi } \end{eqnarray} \tag{ A2 }$$

$$\begin{eqnarray} \boldsymbol{V}(r,\theta)&=\boldsymbol{v_{p}} + \Omega \, r\sin \theta \, \hat{\boldsymbol{e}}_{\phi }, \end{eqnarray} \tag{ A3 }$$

where the poloidal stream function A_ϕ(r, θ, t) and toroidal field B_ϕ(r, θ, t) are used to generate B. The velocity field is time-independent and is prescribed by profiles for the meridional circulation $\boldsymbol{v_{p}}(r,\theta)$ and differential rotation Ω(r, θ).

Rewriting the mean induction Equation (A1) in terms of A_ϕ and B_ϕ, we arrive at two coupled partial differential equations for A_ϕ and B_ϕ,

$$\begin{eqnarray} \frac{\partial {A_{\phi }}}{\partial t} &=\frac{\eta }{\eta _{t}} \left(\nabla ^{2}-\frac{1}{\varpi ^{2}}\right)A_{\phi }- R_{e}\frac{\boldsymbol{v_p}}{\varpi }\boldsymbol{\cdot }\boldsymbol{\nabla }(\varpi A_{\phi })+C_{s}S \nonumber\\ \end{eqnarray} \tag{ A4 }$$

$$\begin{eqnarray} \frac{\partial {B_{\phi }}}{\partial t} &=& \frac{\eta }{\eta _{t}} \left(\nabla ^{2}-\frac{1}{\varpi ^{2}}\right)B_{\phi } +\frac{1}{\varpi }\frac{\partial (\varpi B_{\phi })}{\partial r}\frac{\partial (\eta /\eta _{t})}{\partial r}\nonumber\\ && -R_{e}\varpi \boldsymbol{v_p} \boldsymbol{\cdot }\boldsymbol{\nabla }\left(\frac{B_{\phi }}{\varpi }\right)-R_{e}B_{\phi }\boldsymbol{\nabla } \boldsymbol{\cdot }\boldsymbol{v_p}\nonumber\\ && +C_{\Omega }\varpi [\boldsymbol{\nabla }\boldsymbol{\times }(A_{\phi }\hat{\boldsymbol{e}}_{\phi })]\boldsymbol{\cdot }\boldsymbol{\nabla }\Omega, \end{eqnarray} \tag{ A5 }$$

where ϖ = rsin θ. The contribution to the transport term in the mean induction Equation (A1) that arises from the fluctuating fields, namely, the $\boldsymbol{\nabla }\boldsymbol{\times }\left<\boldsymbol{v^{\prime }}\times \boldsymbol{b^{\prime }}\right>$ term, is present in the A_ϕ equation above and in general is assumed to take a specific form in terms of the mean magnetic field (cf., Babcock 1961; Leighton 1969; Wang & Sheeley 1991; Dikpati & Charbonneau 1999; Jouve & Brun 2007). Here, we use a surface BL term S(r, θ, B_ϕ) for this purpose which serves to generate new poloidal field.

Additionally, Equations (A4) and (A5) have been nondimensionalized by using R_☉ as the characteristic length scale and R²_☉/η_t as the characteristic timescale, where η_t = 10¹¹ cm² s⁻¹ is representative of the turbulent magnetic diffusivity in the convective zone. This rescaling leads to the appearance of three dimensionless control parameters C_Ω = Ω₀R²_☉/η_t, C_s = s₀R_☉/η_t, and the Reynolds number R_e = v₀R_☉/η_t, where Ω₀, s₀, and v₀ are, respectively, the rotation rate and the typical amplitude of the surface source term and of the meridional flow.

Equations (A4) and (A5) are solved with the STELEM (T. Emonet & P. Charbonneau 1998, private communication) code (see Appendix A of Jouve & Brun 2007 for more details) in an annular meridional plane with the colatitude θ ∈ [0, π] and the dimensionless radius r ∈ [0.6, 1], i.e., from slightly below the tachocline (r ≈ 0.7) up to the solar surface R_☉. The STELEM code has been thoroughly tested and validated via an international mean-field dynamo benchmarking process involving eight different codes (Jouve et al. 2008). At the latitudinal boundaries at θ = 0 and θ = π, and at the lower radial boundary at r = 0.6, both A_ϕ and B_ϕ vanish. At the upper radial boundary at r = 1, the solution is matched to an external potential field. Usual initial conditions involve setting a confined dipolar field configuration, i.e., A_ϕ is set to (sin θ)/r² in the convective zone and to 0 below the tachocline. To create the simulation shown in Figure 16, the simulation was initialized using a quadrupolar configuration with an A_ϕ of (3cos θ sin θ)/(2r³) in the convection zone. In both cases, the toroidal field is initialized to 0 everywhere.

The rotation profile used in the series of models discussed in this work captures many aspects of the true solar angular velocity profile, such as deduced from helioseismic inversions (Thompson et al. 2003). We thus assume solid-body rotation below r = 0.66 and a differential rotation above this tachocline interface as given by the following rotation profile:

$$\begin{eqnarray} \Omega (r,\theta) &=& \Omega _c+\frac{1}{2}\left[1+{\rm erf}\left(\frac{2(r-r_c)}{d_1}\right)\right] \nonumber\\ &&\times (\Omega _{\rm eq}+a_2\cos ^2\theta +a_4\cos ^4\theta -\Omega _c). \end{eqnarray} \tag{ A6 }$$

The parameters Ω_eq = 1, Ω_c = 0.93944, r_c = 0.7, d₁ = 0.05, a₂ = −0.136076, and a₄ = −0.145713. With this profile for Ω, the radial shear is maximal at the tachocline.

We assume that the diffusivity in the envelope η is dominated by its turbulent contribution, whereas in the stable interior η_c ≪ η_t. We smoothly match the two different constant values with an error function which enables us to continuously transition from η_c to η_t,

$$\begin{equation} \eta (r)=\eta _c+\frac{\left(\eta _t-\eta _c\right)}{2}\left[1+{\rm erf}\left(\frac{r-r_c}{d}\right)\right], \end{equation} \tag{ A7 }$$

with ${\eta _{\rm c}}=10^9 \,\rm cm^2\ \rm s^{-1}$ and d = 0.03.

In BL flux-transport dynamo models, the poloidal field owes its origin to the tilt of magnetic loops emerging at the solar surface. Thus, the source has to be confined to a thin layer just below the surface and since the process is fundamentally non-local, the source term depends on the variation of B_ϕ at the base of the convection zone. We use the following expression (which is a slightly modified version of that used in Jouve & Brun 2007) in order to better confine the activity belt to low latitudes:

$$\begin{eqnarray} S(r,\theta,B_{\phi }) &=& \frac{1}{2}\left[1+{\rm erf}\left(\frac{r-r_2}{d_2}\right)\right]\left[1-{\rm erf}\left(\frac{r-1}{d_2}\right)\right] \nonumber\\ &&\times \left[1+\left({\frac{B_\phi (r_c,\theta,t)}{B_0}}\right)^2\right]^{-1} \nonumber\\ &&\times \cos \theta \,\sin ^3\theta \,B_\phi (r_c,\theta,t), \end{eqnarray} \tag{ A8 }$$

where r₂ = 0.95, d₂ = 0.01, and B₀ = 10⁵ G. In the particular case of an imposed asymmetry between the north and southern hemisphere (Figure 17) we introduce a modified source term, modulated by the parameter , as follows:

$$\begin{eqnarray} S_{\rm asym}(r,\theta,B_\phi) &=& \frac{1}{2}\left[1+{\rm erf}\left(\frac{r-r_2}{d_2}\right)\right]\left[1-{\rm erf}\left(\frac{r-1}{d_2}\right)\right] \nonumber\\ &&\times \left[1+\left({\frac{B_\phi (r_c,\theta,t)}{B_0}}\right)^2\right]^{-1} \nonumber\\ &&\times (\cos \theta + \epsilon \sin \theta) \sin ^3\theta \, B_\phi (r_c,\theta,t). \end{eqnarray} \tag{ A9 }$$

In BL flux-transport dynamo models, meridional circulation is used to link the two sources of the magnetic field, namely, the base of the convection zone (where toroidal field is created via the latitudinal shear) and the solar surface (where poloidal field is introduced via the BL source term). In the series of models discussed in this paper, the meridional circulation is equatorially symmetric, having one large single cell per hemisphere. Flows are directed poleward at the surface and equatorward at depth (as in the Sun), vanishing at the bottom boundary at r = 0.6. The equatorward branch penetrates slightly beneath the tachocline. To model the single cell meridional circulation, we consider a stream function with the following expression (Jouve et al. 2008):

$$\begin{equation} \psi (r,\theta)=-\frac{2(r-r_b)^2}{\pi (1-r_b)}\sin \left(\frac{\pi (r-r_b)}{1-r_b}\right)\cos \theta \sin \theta, \end{equation} \tag{ A10 }$$

which gives, through the relation $\boldsymbol{v_p}=\boldsymbol{\nabla }\boldsymbol{\times }(\psi \hat{\boldsymbol{e}_\phi })$, the following components of the meridional flow,

$$\begin{eqnarray} v_{r}&=-\frac{2(1-r_b)}{\pi r}\frac{(r-r_b)^2}{(1-r_b)^2}\sin \left(\frac{\pi (r-r_b)}{1-r_b}\right)(3\cos ^2\theta -1) \nonumber\\ \end{eqnarray} \tag{ A11 }$$

$$\begin{eqnarray} v_{\theta } &=& \left[\frac{3r-r_b}{1-r_b} \sin \left(\frac{\pi (r-r_b)}{1-r_b}\right)+\frac{r\pi }{1-r_b}\frac{(r-r_b)}{(1-r_b)}\right. \nonumber\\ &&\times\left. \cos \left(\frac{\pi (r-r_b)}{1-r_b}\right)\right]\frac{2(1-r_b)}{\pi r}\frac{(r-r_b)}{(1-r_b)}\cos \theta \sin \theta, \nonumber\\ \end{eqnarray} \tag{ A12 }$$

with r_b = 0.6.

To introduce asymmetry into the model, an alternative to using the asymmetric source term of Equation (A9) is to introduce an asymmetry into the meridional flow profile. Such an asymmetric meridional flow profile can be constructed using the following stream function:

$$\begin{eqnarray} \psi _{\rm asym}(r,\theta) &=& -\frac{2(r-r_b)^2}{\pi (1-r_b)}\sin \left(\frac{\pi (r-r_b)}{1-r_b}\right) \nonumber\\ &&\times (\cos \theta + \epsilon \sin \theta)\sin \theta, \end{eqnarray} \tag{ A13 }$$

which leads to the following components of the meridional flow,

$$\begin{eqnarray} v_{r,{\rm asym}} &=& -\frac{2(1-r_b)}{\pi r}\frac{(r-r_b)^2}{(1-r_b)^2}\sin \left(\frac{\pi (r-r_b)}{1-r_b}\right) \nonumber\\ &&\times (3\epsilon \sin \theta \cos \theta +3\cos ^2\theta -1) \end{eqnarray} \tag{ A14 }$$

$$\begin{eqnarray} v_{\theta,{\rm asym}} &=& \left[\frac{3r-r_b}{1-r_b} \sin \left(\frac{\pi (r-r_b)}{1-r_b}\right)+\frac{r\pi }{1-r_b}\frac{(r-r_b)}{(1-r_b)}\right. \nonumber\\ &&\times\left. \cos \left(\frac{\pi (r-r_b)}{1-r_b}\right)\right]\frac{2(1-r_b)}{\pi r}\frac{(r-r_b)}{(1-r_b)} \nonumber\\ &&\times (\cos \theta + \epsilon \sin \theta)\sin \theta, \end{eqnarray} \tag{ A15 }$$

again with r_b = 0.6.

Following Gubbins & Zhang (1993) it is straightforward to assess symmetry properties of various mathematical operators and equations. We adopt the superscripts ^A and ^S to indicate whether the scalars or vectors are antisymmetric or symmetric across the equator, respectively. For example, products between a scalar and a vector of the form $a \boldsymbol{F}=\boldsymbol{G}$, where a and F are of like symmetry, yield a symmetric result ($a^S\boldsymbol{F}^S\rightarrow \boldsymbol{G}^S$ and $a^A\boldsymbol{F}^A\rightarrow \boldsymbol{G}^S$), whereas products between quantities of differing symmetries are antisymmetric ($a^A\boldsymbol{F}^S\rightarrow \boldsymbol{G}^A$ and $a^S\boldsymbol{F}^A\rightarrow \boldsymbol{G}^A$). For the vector cross product $\boldsymbol{F}\boldsymbol{\times }\boldsymbol{G}=\boldsymbol{H}$, when the two vectors F and G have the same symmetry properties, the result will be antisymmetric ($\boldsymbol{F}^S\boldsymbol{\times }\boldsymbol{G}^S\rightarrow \boldsymbol{H}^A$ and $\boldsymbol{F}^A\boldsymbol{\times }\boldsymbol{G}^A\rightarrow \boldsymbol{H}^A$), while the cross product between two vectors having opposing symmetries will yield a symmetric result ($\boldsymbol{F}^A\boldsymbol{\times }\boldsymbol{G}^S\rightarrow \boldsymbol{H}^S$). Additionally, the curl operator reverses symmetry ($\boldsymbol{\nabla }\boldsymbol{\times }\boldsymbol{G}^A\rightarrow \boldsymbol{H}^S$ and $\boldsymbol{\nabla }\boldsymbol{\times }\boldsymbol{G}^S\rightarrow \boldsymbol{H}^A$), while the Laplacian operator preserves symmetry ($\nabla ^2\boldsymbol{F}^S\rightarrow \boldsymbol{H}^S$ and $\nabla ^2\boldsymbol{F}^A\rightarrow \boldsymbol{H}^A$).

With these properties established, the analysis of the symmetry properties of the magnetic induction equation,

$$\begin{equation} \frac{\partial \boldsymbol{B}}{\partial t}=\boldsymbol{\nabla }\boldsymbol{\times }(\boldsymbol{V}\boldsymbol{\times }\boldsymbol{B}) + \eta \nabla ^2\boldsymbol{B}, \end{equation} \tag{ A16 }$$

follows in a straightforward manner. For cases possessing a symmetric velocity field $\boldsymbol{V}^S$ with respect to the equator, both terms on the right-hand side of Equation (A16) preserve the symmetry of B. Thus, a dynamo having only a symmetric field $\boldsymbol{B}^S$ will remain symmetric over time, since both the transport term and the diffusion term of Equation (A16) generate symmetric field only. Likewise, a dynamo possessing an antisymmetric field $\boldsymbol{B}^A$ will preserve its antisymmetry over time. Because Equation (A16) is linear in B, it follows that a magnetic field possessing mixed symmetry in the midst of a symmetric velocity field can be considered to be operating two independent, noninteracting dynamos: one that is symmetric and one that is antisymmetric.

However, in cases with an antisymmetric velocity field $\boldsymbol{V}^A$, the transport term on the right-hand side of Equation (A16) provides a mechanism by which the symmetric and antisymmetric modes of B can couple. This coupling arises because an initially symmetric field $\boldsymbol{B}^S$ will first generate an antisymmetric field that in turn leads to the presence of both antisymmetric and symmetric fields, according to the right-hand side of Equation (A16). Analogously, initializing with a purely antisymmetric field $\boldsymbol{B}^A$ will generate fields of mixed symmetry over time.

This analysis procedure can further be applied to the mean-field induction equation. For example, an analysis of Equation (A1) with an α–ω dynamo (i.e., where $\left<\boldsymbol{v}^{\prime }\boldsymbol{\times }\boldsymbol{b}^{\prime }\right>$ is set to $\alpha \left<\boldsymbol{B}\right>$),

$$\begin{equation} \frac{\partial {\left<\boldsymbol{B}\right>}}{\partial t} = \boldsymbol{\nabla }\boldsymbol{\times }\left(\left<\boldsymbol{V}\right> \boldsymbol{\times }\left<\boldsymbol{B}\right> + \alpha \left<\boldsymbol{B}\right> \right) - \eta \nabla ^2 \left(\left<\boldsymbol{B}\right> \right) \end{equation} \tag{ A17 }$$

indicates that, for an assumed symmetric mean velocity field $\left<\boldsymbol{V}\right>^S$ and an antisymmetric alpha effect α^A (which is the natural outcome of helical turbulence in a rotating fluid), such a mean-field dynamo will preserve the symmetry (or antisymmetry) of the initial fields. Hence, as pointed out by Roberts & Stix (1972) and McFadden et al. (1991), a symmetric mean velocity field and an antisymmetric alpha effect do not couple magnetic field modes belonging to different families. Alternatively, if instead an antisymmetric mean flow $\left<\boldsymbol{V}\right>^A$ or a symmetric alpha effect α^S are considered, this now enables a coupling between symmetric and antisymmetric mean fields.

In a similar vein, the BL Equations (A4) and (A5), which are determined by performing the poloidal–toroidal decomposition on Equation (A17), can also be analyzed for symmetry. It is important to note that, by Equation (A2), the poloidal stream function 〈A_ϕ〉 has a symmetry opposite to that of the mean magnetic field $\langle \boldsymbol{B} \rangle$ it generates (and thus also to the corresponding toroidal field 〈B_ϕ〉). We established above that the diffusion term in the mean-field Equation (A17) preserves the symmetry of $\langle \boldsymbol{B}\rangle$, and so it follows that the corresponding diffusion terms in the BL dynamo Equations (A4) and (A5) will serve to preserve the symmetries A_ϕ and B_ϕ. Likewise, because the large-scale transport term $\boldsymbol{\nabla }\boldsymbol{\times } (\langle \boldsymbol{V} \rangle \boldsymbol{\times } \langle\boldsymbol{B}\rangle)$ term in Equation (A17) preserves the symmetry of $\left<\boldsymbol{B}\right>$ whenever $\left<\boldsymbol{V}\right>$ is symmetric, it follows that the analogous terms in the Equations (A4) and (A5) also preserve the symmetry of the system as long as the imposed velocity field is symmetric. For the BL dynamo considered here, an antisymmetric poloidal velocity stream function ψ^A, as in Equation (A10), yields a symmetric meridional flow profile $\boldsymbol{v_p}^S$, since $\boldsymbol{v_p}=\boldsymbol{\nabla }\boldsymbol{\times }(\psi \hat{\boldsymbol{e}_\phi })$, which in turn gives a symmetric mean velocity $\left<\boldsymbol{V}\right>^S$ from Equation (A3). Therefore, the imposed velocity field as defined by Equations (A3) and (A10) will preserve the symmetries of A_ϕ and B_ϕ. Lastly, the source term S as defined by Equation (A8) also preserves the symmetry of A_ϕ, since it is comprised of a series of symmetric coefficients multiplied by cos θ B_ϕ. Therefore, an antisymmetric toroidal field implies a symmetric source term that in turn serves to preserve the symmetry of A_ϕ (and thus $\left<\boldsymbol{B}\right>$), and the same is true when the toroidal field is symmetric.

For these reasons, the dynamo whose characteristics are illustrated in Figure 15, which was initialized with a dipolar field (which is antisymmetric), preserves its antisymmetry with time since all of the terms in Equations (A4) and (A5) preserve the initial symmetry. Indeed, the amplitude of the secondary-family modes (which are symmetric) remain low in this model, as shown in Figures 15(b) and (d). In the dynamos whose characteristics are displayed in Figures 17 or 18, this effect is responsible for the growth of symmetric mean fields, even though both models were initialized with the same antisymmetric mean magnetic field.

It is therefore a direct outcome of symmetry considerations that in standard mean-field dynamo models, either one or the other families of magnetic fields is excited. In the experiments discussed earlier in Section 4.2, we controlled the degree to which the symmetries were mixed via the parameter in Equation (A9) and in Equations (A14) and (A15), which led to the dynamos illustrated in Figures 17 or 18, respectively. In both cases, was chosen to yield a dynamo where the end state contained secondary-family amplitudes of about 25%, as is observed on the Sun; other choices of will lead to end states with different ratios.