1 Introduction

Looking at a full disk magnetogram (a map showing spatially the line of sight flux density of the magnetic field) of the solar photosphere one sees that the most prominent large scale pattern of magnetic flux concentrations on the solar surface are the bipolar active regions (see Figure 1).

Figure 1:
figure 1

A full disk magnetogram from the Kitt Peak Solar Observatory showing the line of sight magnetic flux density on the photosphere of the Sun on May 11, 2000. White (Black) color indicates a field of positive (negative) polarity.

When observed in white light (see Figure 2), an active region usually contains sunspots and is sometimes called a sunspot group.

Figure 2:
figure 2

A continuum intensity image of the Sun taken by the MDI instrument on board the SOHO satellite on the same day as Figure 1. It shows the sunspots that are in some of the active regions in Figure 1.

Active regions are so named because they are centers of various forms of solar activity (such as solar flares) and sites of X-ray emitting coronal loops (see Figure 3).

Figure 3:
figure 3

A full disk soft X-ray image of the solar coronal taken on the same day as Figure 1 from the soft X-ray telescope on board the Yohkoh satellite. Active regions appear as sites of bright X-ray emitting loops.

Despite the turbulent nature of solar convection which is visible from the granulation pattern on the photosphere, the large scale bipolar active regions show remarkable order and organization as can be seen in Figure 1. The active regions are roughly confined into two latitudinal belts which are located nearly symmetrically on the two hemispheres. Over the course of each 11-year solar cycle, the active region belts march progressively from mid-latitude of roughly 35° toward the equator on both hemispheres (Maunder, 1922). The polarity orientations of the bipolar active regions are found to obey the well-known Hale polarity law (Hale et al., 1919; Hale and Nicholson, 1925) outlined as follows. The line connecting the centers of the two magnetic polarity areas of each bipolar active region is usually nearly east-west oriented. Within each 11-year solar cycle, the leading polarities (leading in the direction of solar rotation) of nearly all active regions on one hemisphere are the same and are opposite to those on the other hemisphere (see Figure 1), and the polarity order reverses on both hemispheres with the beginning of the next cycle. The magnetic fields at the solar north and south poles are also found to reverse sign every 11 years near sunspot maximum (i.e. near the middle of a solar cycle). Therefore, the complete magnetic cycle, which corresponds to the interval between successive appearances at mid-latitudes of active regions with the same polarity arrangement, is in fact 22 years.

Besides their highly organized behavior during each solar cycle, active regions are found to possess some interesting asymmetries between their leading and following polarities. Observations show that the axis connecting the leading and the following polarities of each active region is nearly east-west oriented (or toroidal) but on average shows a small tilt relative to the east-west direction with the leading polarity of the region being slightly closer to the equator than the following (see Figure 1). This small mean tilt angle is found to increase approximately linearly with the latitude of the active region (Wang and Sheeley Jr, 1989, 1991; Howard, 1991a,b; Fisher et al., 1995; Kosovichev and Stenflo, 2008). This observation of active region tilts is originally summarized in Hale et al. (1919) and is generally referred to as Joy’s law. Note that Joy’s law describes the statistical mean behavior of the active region tilts. The tilt angles of individual active regions also show a large scatter about the mean (Wang and Sheeley Jr, 1989, 1991; Howard, 1991a,b; Fisher et al., 1995; Kosovichev and Stenflo, 2008). Another intriguing asymmetry is found in the morphology of the leading and the following polarities of an active region. The flux of the leading polarity tends to be concentrated in large well-formed sunspots, whereas the flux of the following polarity tends to be more dispersed and to have a fragmented appearance (see Bray and Loughhead, 1979). Observations also show that the magnetic inversion lines (the neutral lines separating the fluxes of the two opposite polarities) in bipolar active regions are statistically nearer to the main following polarity spot than to the main leading spot (van Driel-Gesztelyi and Petrovay, 1990; Petrovay et al., 1990). Furthermore for young growing active regions, there is an asymmetry in the east-west proper motions of the two polarities, with the leading polarity spots moving prograde more rapidly than the retrograde motion of the following polarity spots (see Chou and Wang, 1987; van Driel-Gesztelyi and Petrovay, 1990; Petrovay et al., 1990).

More recently, vector magnetic field observations of active regions on the photosphere have shown that active region magnetic fields have a small but statistically significant mean twist that is left-handed in the northern hemisphere and right-handed in the southern hemisphere (see Pevtsov et al., 1995, 2001). The twist is measured in terms of the quantity α ≡ 〈J z/B z〉, the ratio of the vertical electric current over the vertical magnetic field averaged over an active region. The measured α for individual solar active regions show considerable scatter, but there is clearly a statistically significant trend for negative α (left-handed field line twist) in the northern hemisphere and positive α (right-handed field line twist) in the southern hemisphere. In addition, soft X-ray observations of solar active regions sometimes show hot plasma of S or inverse-S shapes called “sigmoids” with the northern hemisphere preferentially showing inverse-S shapes and the southern hemisphere forward-S shapes (Rust and Kumar, 1996; Pevtsov et al., 2001, see Figure 4 for an example).

Figure 4:
figure 4

A soft X-ray image of the solar coronal on May 27, 1999, taken by the Yohkoh soft X-ray telescope. The arrows point to two “sigmoids” at similar longitudes north and south of the equator showing an inverse-S and a forward-S shape respectively.

This hemispheric preference of the sign of the active region field line twist and the direction of X-ray sigmoids do not change with the solar cycle (see Pevtsov et al., 2001). New high resolution vector magnetic field observations from the Hinode space mission show further evidence for the emergence of twisted active region magnetic flux in association with the formation of active region filaments (see review by Lites, 2009).

The cyclic large scale magnetic field of the Sun with a period of 22 years is believed to be sustained by a dynamo mechanism (see e.g. review by Charbonneau, 2005). The Hale polarity law of solar active regions indicates the presence of a large scale subsurface toroidal magnetic field generated by the solar cycle dynamo mechanism. The picture of how and where the large scale solar dynamo operates has undergone substantial revision due in part to new knowledge from helioseismology regarding the solar internal rotation profile (see Deluca and Gilman, 1991; Gilman, 2000). Evidence now points to the tachocline, the thin shear layer at the base of the solar convection zone, where solar rotation changes from the latitudinal differential rotation of the solar convective envelope to the nearly solid-body rotation of the radiative interior, as the site for the generation and amplification of the large scale toroidal magnetic field from a weak poloidal magnetic field (see Charbonneau and MacGregor, 1997; Dikpati and Charbonneau, 1999; Dikpati and Gilman, 2001). Furthermore, with its stable (weakly) subadiabatic stratification, the thin overshoot region in the upper part of the tachocline layer (Gilman, 2000) allows storage of strong toroidal magnetic fields against their magnetic buoyancy for time scales comparable to the solar cycle period (Parker, 1975, 1979; van Ballegooijen, 1982; Moreno-Insertis et al., 1992; Fan and Fisher, 1996; Moreno-Insertis et al., 2002; Rempel, 2003). Thus with toroidal magnetic fields being generated and stored in the tachocline layer at the base of the solar convection zone, these fields need to traverse the entire convection zone before they can emerge at the photosphere to form the observed solar active regions.

High resolution observations have shown that magnetic fields on the solar photosphere are in a fibril state, i.e. contain discrete flux tubes of high field strength (B ≳ 103 G in equipartition with the thermal pressure) having a hierarchy of cross-sectional sizes that range from sunspots of active regions down to below the limit of observational resolution (see Zwaan, 1987; Stenflo, 1989; Domínguez Cerdeña, et al., 2003; Khomenko et al., 2003; Socas-Navarro and Sánchez Almeida, 2003; Orozco Suárez et al., 2007). It is thus likely that the subsurface magnetic fields in the solar convection zone are also concentrated into discrete flux tubes. One mechanism that can concentrate magnetic flux in a turbulent conducting fluid, such as the solar convection zone, into high field strength flux tubes is the process known as “flux expulsion”, i.e. magnetic fields are expelled from the interior of convecting cells into the boundaries. This process has been studied by MHD simulations of the interaction between convection and magnetic fields (see Galloway and Weiss, 1981; Nordlund et al., 1992). In particular, the 3D simulations of magnetic fields in convecting flows by Nordlund et al. (1992) show the formation of strong discrete flux tubes in the vicinity of strong downdrafts. In addition, Parker (put forth an interesting argument that supports the fibril form of magnetic fields in the solar convection zone. He points out that although the magnetic energy is increased by the compression from a continuum field into the fibril state, the total energy of the convection zone (thermal + gravitational + magnetic) is reduced by the fibril state of the magnetic field by avoiding the magnetic inhibition of convective overturning. Assuming an idealized polytropic atmosphere, he was able to derive the filling factor of the magnetic fields that corresponds to the minimum total energy state of the atmosphere. By applying an appropriate polytropic index for the solar convection zone, he computed the filling factor which yielded fibril magnetic fields of about 1–5 kG, roughly in agreement with the observed fibril fields at the solar surface.

Since both observational evidence and theoretical arguments support the fibril picture of solar magnetic fields, the idealized concept of isolated magnetic flux tubes surrounded by “field-free” plasmas has been developed and widely used in modeling magnetic fields in the solar convection zone (see Parker, 1979; Spruit, 1981; Vishniac, 1995a,b), even though in reality there is no field-free region in the quite-sun photosphere (e.g Orozco Suárez et al., 2007). The manner in which individual bipolar active regions emerge at the photosphere (see Zwaan, 1987) and the well-defined order of the active regions as described by the Hale polarity rule suggest that they correspond to coherent and discrete flux tubes rising through the solar convection zone and reaching the photosphere in a reasonably cohesive fashion, not severely distorted by convection. It is this process — the formation of buoyant flux tubes from the toroidal magnetic field stored in the overshoot region and their dynamic rise through the convection zone to form solar active regions — that is the central focus of this review.

It should be noted that bipolar magnetic regions emerge on the photosphere with a wide range of size scales that span at least 5 orders of magnitude (from below 1018 Mx to 1023 Mx) in absolute flux content, ranging from the large, sunspot-containing active regions to small, ephemeral regions (ERs) that appear in the quiet Sun (e.g. Harvey, 1993; Hagenaar et al., 2008). The well organized pattern and cycle dependence as described by the butterfly diagram, the Hale polarity rule, and Joy’s law exhibited by active regions are progressively less well obeyed by the smaller bipoles. Small ERs emerge in both the closed-field, mixed polarity quite-sun regions as well as the more unipolar coronal hole regions (e.g. Hagenaar et al., 2008). The nature and origin of ERs are not certain. The ER flux may originate close to the surface, produced by a “local dynamo” due to small-scale convective motions near the surface (e.g. Cattaneo et al., 2003; Bercik et al., 2005). Alternatively, ERs may correspond to flux sheared off from emerging or decaying active region flux tubes. The study by Hagenaar et al. (2008) using MDI magnetogram sequences have shown an interesting dependence of the ER emergence rate on the local flux imbalance, with lower emergence rate in regions of larger flux imbalance. This functional dependence is found to be the same for both the closed-field quiet-sun regions and the more unipolar coronal holes. Such a dependence of the ER emergence rate may, however, result from either of the above two scenarios for the origin of ERs (Hagenaar et al., 2008).

High resolution vector magnetic field observations from the Hinode satellite have also revealed new, unprecedented details of the photospheric magnetic field at the solar polar region (Tsuneta et al., 2008). It is found that the polar magnetic field is characterized by unipolar vertical kilogauss patches with super-equipartition field strength, and ubiquitous weaker transient horizontal fields. The origin of these unipolar strong flux tubes is not clear but they may simply be the surviving fragments of the following polarity of the decaying active regions being transported to the polar region through the combined actions of diffusion by supergranular motion and advection by meridional flows (see e.g. Wang et al., 1989).

Although magnetic fields are generated on all scales in the solar convection zone (e.g Schüssler, 2002), in this review, we only focus on the emergence process of active region scale flux tubes, which are generally thought to originate from a deep seated solar cycle dynamo operating at the base of the solar convection zone (see e.g. review by Charbonneau, 2005). The remainder of the review will be organized as follows.

  • Section 2 gives a brief overview of some of the simplifying models and computational approaches that have been applied to studying the dynamic evolution of magnetic flux tubes in the solar convection zone. In particular, the thin flux tube model is proven to be a very useful tool for understanding the global dynamics of emerging active region flux tubes in the solar convective envelope and (as discussed in the later sections) has produced results that explain the origin of several basic observed properties of solar active regions.

  • Section 3 discusses the storage and equilibrium properties of large scale toroidal magnetic fields in the stable overshoot region below the solar convection zone.

  • Section 4 focuses on the buoyancy instabilities associated with the equilibrium toroidal magnetic fields and the formation of buoyant flux tubes from the base of the solar convection zone.

  • Section 5 reviews results on the dynamic evolution of emerging flux tubes in the solar convection zone.

    • Section 5.1 discusses major findings from various thin flux tube simulations of emerging flux loops.

    • Section 5.2 discusses some of the recent studies using local helio-seismology techniques to look for thermodynamic perturbations and plasma flow signatures that may be produced by subsurface emerging flux.

    • Section 5.3 discusses the observed hemispheric trend of the twist of the magnetic field in solar active regions and the models that explain its origin.

    • Section 5.4 reviews results from direct MHD simulations with regard to the minimum twist necessary for tube cohesion.

    • Section 5.5 discusses further constraint on the twist of subsurface emerging tubes due to Joy’s law of active region tilts and the development of asymmetric twist between the leading and following sides of the emerging tube, based on results from a set of 3D spherical-shell anelastic MHD simulations.

    • Section 5.6 discusses the kink evolution of highly twisted emerging tubes.

    • Section 5.7 reviews the influence of 3D stratified convection on the evolution of buoyant flux tubes.

  • Section 6 discusses results from 3D MHD simulations on the asymmetric transport of magnetic flux (or turbulent pumping of magnetic fields) by stratified convection penetrating into a stable overshoot layer.

  • Section 7 discusses an alternative mechanism of magnetic flux amplification by converting the potential energy associated with the stratification of the convection zone into magnetic energy.

  • Section 8 gives a brief overview of our current understanding of active region flux emergence into the atmosphere and the post-emergence evolution of the subsurface fields.

  • Section 9 summarizes and discusses the basic conclusions.

2 Models and Computational Approaches

2.1 The thin flux tube model

The well-defined order of the solar active regions (see description of the observational properties in Section 1) suggests that their precursors at the base of the solar convection zone should have a field strength that is at least B eq, where B eq is the field strength that is in equipartition with the kinetic energy density of the convective motions: \(B_{{\rm{eq}}}^2/8\pi = \rho v_{\rm{c}}^2\) If we use the results from the mixing length models of the solar convection zone for the convective flow speed v c, then we find that in the deep convection zone B eq is on the order of 104 G. In the past two decades, direct 3D numerical simulations have led to a new picture for solar convection that is non-local, driven by the concentrated downflow plumes formed by radiative cooling at the surface layer, and with extreme asymmetry between the upward and downward flows (see reviews by Spruit et al., 1990; Spruit, 1997). Hence it should be noted that the B eq derived based on the local mixing length description of solar convection may not really reflect the intensity of the convective flows in the deep solar convection zone. With this caution in mind, we nevertheless refer to B eq ∼ 104 G as the field strength in equipartition with convection in this review.

Assuming that in the deep solar convection zone the magnetic field strength for flux tubes responsible for active region formation is at least 104 G, and given that the amount of flux observed in solar active regions ranges from ∼ 1020 Mx to 1022 Mx (see Zwaan, 1987), one then finds that the cross-sectional sizes of the flux tubes are small in comparison to other spatial scales of variation, e.g. the pressure scale height. For an isolated magnetic flux tube that is thin in the sense that its cross-sectional radius a is negligible compared to both the scale height of the ambient unmagnetized fluid and any scales of variation along the tube, the dynamics of the flux tube may be simplified with the thin flux tube approximation (see Spruit, 1981; Longcope and Klapper, 1997) which corresponds to the lowest order in an expansion of MHD in powers of a/L, where L represents any of the large length scales of variation. Under the thin flux tube approximation, all physical quantities of the tube, such as position, velocity, field strength, pressure, density, etc. are assumed to be averages over the tube cross-section and they vary spatially only along the tube. Furthermore, because of the much shorter sound crossing time over the tube diameter compared to the other relevant dynamic time scales, an instantaneous pressure balance is assumed between the tube and the ambient unmagnetized fluid:

$$p + {{{B^2}} \over {8\pi }} = {p_{\rm{e}}}$$
(1)

where p is the tube internal gas pressure, B is the tube field strength, and p e is the pressure of the external fluid. Applying the above assumptions to the ideal MHD momentum equation, Spruit (1981) derived the equation of motion of a thin untwisted magnetic flux tube embedded in a field-free fluid. Taking into account the differential rotation of the Sun, Ω e(r) = Ωe(r), the equation of motion for the thin flux tube in a rotating reference frame of angular velocity Ω = Ωis (Ferriz-Mas and Schüssler, 1993; Caligari et al., 1995)

$$\matrix{{\rho {{d{\bf{v}}} \over {dt}} = 2\rho ({\bf{v}} \times \Omega) + \rho ({\Omega ^2} - \Omega _{\rm{e}}^2)\varpi \hat \varpi + (\rho - {\rho _{\rm{e}}}){{\bf{g}}_{{\rm{eff}}}}} \hfill \cr {\quad \quad \quad + {\bf{\hat l}}{\partial \over {\partial s}}\left({{{{B^2}} \over {8\pi }}} \right) + {{{B^2}} \over {4\pi }}{\bf{k}} - {C_{\rm{D}}}{{{\rho _{\rm{e}}}\vert{{({{\bf{v}}_{{\rm{rel}}}})}_ \bot }\vert{{({{\bf{v}}_{{\rm{rel}}}})}_ \bot }} \over {\pi {{(\Phi/B)}^{1/2}}}},} \hfill \cr}$$
(2)

where

$${{\bf{g}}_{{\rm{eff}}}} = {\bf{g}} + \Omega _{\rm{e}}^2\varpi \hat \varpi,$$
(3)
$${({{\bf{v}}_{{\rm{rel}}}})_ \bot } = {[{\bf{v}} - ({{\bf{\Omega }}_{\rm{e}}} - {\bf{\Omega }}) \times {\bf{r}}]_ \bot }.$$
(4)

In the above, r, v, B, p, ρ, denote the position vector, velocity, magnetic field strength, plasma pressure and density of a Lagrangian tube element respectively, each of which is a function of time t and the arc-length s measured along the tube, ρ e(r) denotes the external density at the position r of the tube element, is the unit vector pointing in the direction of the solar rotation axis, \(\hat{\varpi}\) denotes the unit vector perpendicular to and pointing away from the rotation axis at the location of the tube element and ϖ denotes the distance to the rotation axis, \(\hat{\bf l}\equiv\partial{\bf r}/\partial s\) is the unit vector tangential to the flux tube, k 2r/∂s 2 is the tubes curvature vector, the subscript ⊥ denotes the vector component perpendicular to the local tube axis, g is the gravitational acceleration, and C D is the aerodynamic drag coefficient which is believed to be of order unity. The drag term (the last term on the right hand side of the equation of motion (2)) is added to approximate the opposing force experienced by the flux tube as it moves relative to the ambient fluid. The term is derived based on the case of incompressible flows past a rigid cylinder under high Reynolds number conditions, in which a turbulent wake develops behind the cylinder, creating a pressure difference between the up- and down-stream sides and hence a drag force on the cylinder (see Batchelor, 1967).

If one considers only the solid body rotation of the Sun, then the Equations (2), (3), and (4) can be simplified by letting Ω e = Ω. Calculations using the thin flux tube model (see Section 5.1) have shown that the effect of the Coriolis force 2ρ(v × Ω) acting on emerging flux loops can lead to east-west asymmetries in the loops that explain several well-known properties of solar active regions.

Note that in the equation of motion (2), the effect of the “enhanced inertia” caused by the back-reaction of the fluid to the relative motion of the flux tube is completely ignored. This effect has sometimes been incorporated by treating the inertia for the different components of Equation (2) differently, with the term ρ(dv/dt) on the left-hand-side of the perpendicular component of the equation being replaced by (ρ + ρ e)(dv/dt) (see Spruit, 1981). This simple treatment is problematic for curved tubes and the proper ways to treat the back-reaction of the fluid are controversial in the literature (Cheng, 1992; et al., 1994; Moreno-Insertis et al., 1996; Osin et al., 1999). Since the enhanced inertial effect is only significant during the impulsive acceleration phases of the tube motion, which occur rarely in the thin flux tube calculations of emerging flux tubes, and the results obtained do not depend significantly on this effect, many later calculations have taken the approach of simply ignoring it (see Caligari et al., 1995, 1998; Fan and Fisher, 1996).

Equations (1) and (2) are to be complemented by the following equations to completely describe the dynamic evolution of a thin untwisted magnetic flux tube:

$${d \over {dt}}\left({{B \over \rho }} \right) = {B \over \rho }\left[ {{{\partial ({\bf{v}} \cdot {\bf{\hat l}})} \over {\partial s}} - {\bf{v}} \cdot {\bf{k}}} \right],$$
(5)
$${1 \over \rho }{{d\rho } \over {dt}} = {1 \over {\gamma p}}{{dp} \over {dt}} - {{{\nabla _{{\rm{ad}}}}} \over p}{{dQ} \over {dt}},$$
(6)
$$p = {{\rho RT} \over \mu },$$
(7)

where ∇ad ≡ ( ln T/∂ ln p)s. Equation (5) describes the evolution of the tube magnetic field and is derived from the ideal MHD induction equation (Spruit, 1981). Equation (6) is the energy equation for the thin flux tube (Fan and Fisher, 1996), in which dQ/dt corresponds to the volumetric heating rate of the flux tube by non-adiabatic effects, e.g. by radiative diffusion (Section 3.2). Equation (7) is simply the equation of state for an ideal gas. Thus the five Equations (1), (2), (5), (6), and (7) completely determine the evolution of the five dependent variables v (t, s), B (t, s), p (t, s), p (t, s), and T (t, s) for each Lagrangian tube element of the thin flux tube.

Spruit’s original formulation for the dynamics of a thin isolated magnetic flux tube as described above assumes that the tube consists of untwisted flux \({\bf B}=B\hat{\bf l}\). Longcope and Klapper (1997) extend the above model to include the description of a weak twist of the flux tube, assuming that the field lines twist about the axis at a rate q whose magnitude is 2π/L w, where L w is the distance along the tube axis over which the field lines wind by one full rotation and |qa| ≪ 1. Thus in addition to the axial component of the field B, there is also an azimuthal field component in each tube cross-section, which to lowest order in qa is given by B θ = qr B, where r denotes the distance to the tube axis. An extra degree of freedom for the motion of the tube element — the spin of the tube cross-section about the axis — is also introduced, whose rate is denoted by w (angle per unit time). By considering the kinematics of a twisted ribbon with one edge corresponding to the tube axis and the other edge corresponding to a twisted field line of the tube, Longcope and Klapper (1997) derived an equation that describes the evolution of the twist q in response to the motion of the tube:

$${{dq} \over {dt}} = - {{d\ln \delta s} \over {dt}}q + {{\partial \omega } \over {\partial s}} + ({\bf{\hat l}} \times {\bf{k}}) \cdot {{d{\bf{\hat l}}} \over {dt}},$$
(8)

where δs denotes the length of a Lagrangian tube element. The first term on the right-hand-side describes the effect of stretching on q: Stretching the tube reduces the rate of twist q. The second term is simply the change of q resulting from the gradient of the spin along the tube. The last term is related to the conservation of total magnetic helicity which, for the thin flux tube structure, can be decomposed into a twist component corresponding to the twist of the field lines about the axis, and a writhe component corresponding to the “helicalness” of the axis (see discussion in Longcope and Klapper, 1997). It describes how the writhing motion of the tube axis can induce twist of the opposite sense in the tube.

Furthermore, by integrating the stresses over the surface of a tube segment, Longcope and Klapper (1997) evaluated the forces experienced by the tube segment. They found that for a weakly twisted (|qa| ≪ 1) thin tube (|a∂ s| ≪ 1, where s denotes the inverse of the length scale of variation along the tube), the equation of motion of the tube axis differs very little from that for an untwisted tube — the leading order term in the difference is O[qa 2 s] (see also Ferriz-Mas and Schüssler, 1990). Thus the equation of motion (2) applies also to a weakly twisted thin flux tube. By further evaluating the torques exerted on a tube segment, Longcope and Klapper (1997) also derived an equation for the evolution of the spin w:

$${{d\omega } \over {dt}} = - {2 \over a}{{da} \over {dt}}\omega + v_{\rm{a}}^2{{\partial q} \over {\partial s}},$$
(9)

where \(v_{\rm a}=B/\sqrt{4\pi\rho}\) is the Alfvén speed. The first term on the right hand side simply describes the decrease of spin due to the expansion of the tube cross-section as a result of the tendency to conserve angular momentum. The second term, in combination with the second term on the right hand side of Equation (8), describes the propagation of torsional Alfvén waves along the tube.

The two new Equations (8) and (9) — derived by Longcope and Klapper (1997) — together with the earlier Equations (1), (2), (5), (6), and (7) provide a description for the dynamics of a weakly twisted thin flux tube. Note that the two new equations are decoupled from and do not have any feedback on the solutions for the dependent variables described by the earlier equations. One can first solve for the motion of the tube axis using Equations (1), (2), (5), (6), and (7), and then apply the resulting motion of the tube axis to Equations (8) and (9) to determine the evolution of the twist of the tube. If the tube is initially twisted, then the twist q can propagate and redistribute along the tube as a result of stretching (1st term on the right-hand-side of Equation (8)) and the torsional Alfvén waves (2nd term on the right-hand-side of Equation (8)). Twist can also be generated due to writhing motion of the tube axis (last term on the right-hand-side of Equation (8)), as required by the conservation of total helicity.

2.2 MHD simulations

The thin flux tube (TFT) model described above is physically intuitive and computationally tractable. It provides a description of the dynamic motion of the tube axis in a three-dimensional space, taking into account large scale effects such as the curvature of the convective envelope and the Coriolis force due to solar rotation. The Lagrangian treatment of each tube segment in the TFT model allows for preserving perfectly the frozen-in condition of the tube plasma. Thus there is no magnetic diffusion in the TFT model. However, the TFT model ignores variations within each tube cross-section. It is only applicable when the flux tube radius is thin (Section 2.1) and the tube remains a cohesive object (Section 5.4). Clearly, to complete the picture, direct MHD calculations that resolve the tube cross-section and its interaction with the surrounding fluid are needed. On the other hand, direct MHD simulations that discretize the spatial domain are subject to numerical diffusion. The need to adequately resolve the flux tube — so that numerical diffusion does not have a significant impact on the dynamical processes of interest (e.g. the variation of magnetic buoyancy) — severely limits the spatial extent of the domain that can be modeled. So far the MHD simulations cannot address the kinds of large scale dynamical effects that have been studied by the TFT model (Section 5.1). Thus the TFT model and the resolved MHD simulations complement each other.

For the bulk of the solar convection zone, the fluid stratification is very close to being adiabatic with δ ≪ 1, where δ ≡ ∇ − ∇ad is the non-dimensional superadiabaticity with ∇ = d ln T/d ln p and ∇ad = (d ln T/d ln p)ad denoting the actual and the adiabatic logarithmic temperature gradient of the fluid respectively, and the convective flow speed v c is expected to be much smaller than the sound speed c s: v c/c sδ 1/2 ≪ 1 (see Schwarzschild, 1958; Lantz, 1991). Furthermore, the plasma β defined as the ratio of the thermal pressure to the magnetic pressure (βp/(B 2/8π)) is expected to be very high (β ≫ 1) in the deep convection zone. For example for flux tubes with field strengths of order 105 G, which is significantly super-equipartition compared to the kinetic energy density of convection, the plasma β is of order 105. Under these conditions, a very useful computational approach for modeling subsonic magnetohydrodynamic processes in a pressure dominated plasma is the well-known anelastic approximation (see Gough, 1969; Gilman and Glatzmaier, 1981; Glatzmaier, 1984; Lantz and Fan, 1999). The main feature of the anelastic approximation is that it filters out the sound waves so that the time step of numerical integration is not limited by the stringent acoustic time scale which is much smaller than the relevant dynamic time scales of interest as determined by the flow velocity and the Alfvén speed.

Listed below is the set of anelastic MHD equations (see Gilman and Glatzmaier, 1981; Lantz and Fan, 1999, for details of the derivations):

$$\nabla \cdot ({\rho _0}{\bf{v}}) = 0,$$
(10)
$${\rho _0}\left[ {{{\partial {\bf{v}}} \over {\partial t}} + ({\bf{v}}\cdot\nabla){\bf{v}}} \right] = - \nabla {p_1} + {\rho _1}{\bf{g}} + {1 \over {4\pi }}(\nabla \times {\bf{B}}) \times {\bf{B}} + \nabla \cdot{\bf{\Pi }},$$
(11)
$${\rho _0}{T_0}\left[ {{{\partial {s_1}} \over {\partial t}} + ({\bf{v}}\cdot\nabla)({s_0} + {s_1})} \right] = \nabla \cdot(K{\rho _0}{T_0}\nabla {s_1}) + {1 \over {4\pi }}\eta \vert\nabla \times {\bf{B}}{\vert^2} + ({\bf{\Pi }}\cdot\nabla)\cdot{\bf{v}},$$
(12)
$$\nabla \cdot {\bf{B}} = 0,$$
(13)
$${{\partial {\bf{B}}} \over {\partial t}} = \nabla \times ({\bf{v}} \times {\bf{B}}) - \nabla \times (\eta \nabla \times {\bf{B}}),$$
(14)
$${{{\rho _1}} \over {{\rho _0}}} = {{{p_1}} \over {{p_0}}} - {{{T_1}} \over {{T_0}}},$$
(15)
$${{{s_1}} \over {{c_p}}} = {{{T_1}} \over {{T_0}}} - {{\gamma - 1} \over \gamma }{{{p_1}} \over {{p_0}}},$$
(16)

where s 0(z), p 0(z), ρ 0(z), and T 0(z) correspond to a time-independent, background reference state of hydrostatic equilibrium and nearly adiabatic stratification, and velocity v, magnetic field B, thermodynamic fluctuations s 1, p 1, ρ 1, and T 1 are the dependent variables to be solved that describe the changes from the reference state. The quantity Π is the viscous stress tensor given by

$${\Pi _{ij}} \equiv \mu \left({{{\partial {v_i}} \over {\partial {x_j}}} + {{\partial {v_j}} \over {\partial {x_i}}} - {2 \over 3}(\nabla \cdot{\bf{v}}){\delta _{ij}}} \right),$$

and μ, K and η denote the dynamic viscosity, and thermal and magnetic diffusivity, respectively. The anelastic MHD equations (10)(16) are derived based on a scaled-variable expansion of the fully compressible MHD equations in powers of δ and β −1, which are both assumed to be quantities ≪ 1. To first order in δ, the continuity equation (10) reduces to the statement that the divergence of the mass flux equals to zero. As a result sound waves are filtered out, and pressure is assumed to adjust instantaneously in the fluid as if the sound speed was infinite. Although the time derivative of density no longer appears in the continuity equation, density ρ 1 does vary in space and time and the fluid is compressible but on the dynamic time scales (as determined by the flow speed and the Alfvén speed) not on the acoustic time scale, thus allowing convection and magnetic buoyancy to be modeled in the highly stratified solar convection zone. Fan (2001a) has shown that the anelastic formulation gives an accurate description of the magnetic buoyancy instabilities under the conditions of high plasma β and nearly adiabatic stratification.

Fully compressible MHD simulations have also been applied to study the dynamic evolution of a magnetic field in the deep solar convection zone using non-solar but reasonably large β values such as β ∼ 10 to 1000. In several cases comparisons have been made between fully compressible simulations using large plasma β and the corresponding anelastic MHD simulations, and good agreement was found between the results (see Fan et al., 1998a; Rempel, 2002). Near the top of the solar convection zone, neither the TFT model nor the anelastic approximation are applicable because the active region flux tubes are no longer thin (Moreno-Insertis, 1992) and the velocity field is no longer subsonic. Fully compressible MHD simulations are necessary for modeling flux emergence near the surface (Section 8).

3 Equilibrium Conditions of Toroidal Magnetic Fields Stored at the Base of the Solar Convection Zone

3.1 The mechanical equilibria for an isolated toroidal flux tube or an extended magnetic layer

The Hale’s polarity rule of solar active regions indicates a subsurface magnetic field that is highly organized, of predominantly toroidal direction, and with sufficiently strong field strength (super-equipartition compared to the kinetic energy density of convection) such that it is not subjected to strong deformation by convective motions. It is argued that the weakly subadiabatically stratified overshoot layer at the base of the solar convection zone is the most likely site for the storage of such a strong coherent toroidal magnetic field against buoyant loss for time scales comparable to the solar cycle period (see Parker, 1979; van Ballegooijen, 1982).

It is not clear if the toroidal magnetic field is in the state of isolated flux tubes or stored in the form of a more diffuse magnetic layer. Moreno-Insertis et al. (1992) have considered the mechanical equilibrium of isolated toroidal magnetic flux tubes (flux rings) in a subadiabatic layer using the thin flux tube approximation (Section 2.1). The forces experienced by an isolated toroidal flux ring at the base of the convection zone is illustrated in Figure 5(a).

Figure 5:
figure 5

Schematic illustrations based on Schüssler and Rempel (2002) of the various forces involved with the mechanical equilibria of an isolated toroidal flux ring (a) and a magnetic layer (b) at the base of the solar convection zone. In the case of an isolated toroidal ring (see the black dot in (a) indicating the location of the tube cross-section), the buoyancy force has a component parallel to the rotation axis, which cannot be balanced by any other forces. Thus mechanical equilibrium requires that the buoyancy force vanishes and the magnetic curvature force is balanced by the Coriolis force resulting from a prograde toroidal flow in the flux ring. For a magnetic layer (as indicated by the shaded region in (b)), on the other hand, a latitudinal pressure gradient can be built up, so that an equilibrium may also exist where a non-vanishing buoyancy force, the magnetic curvature force and the pressure gradient are in balance with vanishing Coriolis force (vanishing longitudinal flow).

The condition of total pressure balance (1) and the presence of a magnetic pressure inside the flux tube require a lower gas pressure inside the flux tube compared to the outside. Thus either a lower density or a lower temperature (or a combination of the two) inside the flux tube is needed to achieve the lower gas pressure required for pressure balance. If the flux tube is in thermal equilibrium with the surrounding, then the density inside needs to be lower and the flux tube is buoyant. The buoyancy force associated with a magnetic flux tube in thermal equilibrium with its surrounding is often called the magnetic buoyancy (Parker, 1975). It can be seen in Figure 5(a) that a radially directed buoyancy force has a component that is parallel to the rotation axis, which cannot be balanced by any other forces associated with the toroidal flux ring. Thus for the toroidal flux ring to be in mechanical equilibrium, the tube needs to be in a neutrally buoyant state with vanishing buoyancy force, and with the magnetic curvature force pointing towards the rotation axis being balanced by a Coriolis force produced by a faster rotational speed of the flux ring (see Figure 5(a)). Such a neutrally buoyant flux ring (with equal density between inside and outside) then requires a lower internal temperature than the surrounding plasma to satisfy the total pressure balance. If one starts with a toroidal flux ring that is initially in thermal equilibrium with the surrounding and rotates at the same ambient angular velocity, then the flux ring will move radially outward due to its buoyancy and latitudinally poleward due to the unbalanced poleward component of the tension force. As a result of its motion, the flux ring will lose buoyancy due to the subadiabatic stratification and attain a larger internal rotation rate with respect to the ambient field-free plasma due to the conservation of angular momentum, evolving towards a mechanical equilibrium configuration. The flux ring will undergo superposed buoyancy and inertial oscillations around this mechanical equilibrium state. It is found that the oscillations can be contained within the stably stratified overshoot layer and also within a latitudinal range of Δθ ≲ 20° to be consistent with the active region belt, if the field strength of the toroidal flux ring B ≲ 105 G and the subadiabaticity of the overshoot layer is sufficiently strong with δ ≡ ∇ − ∇ad ≲ −10−5, where ∇ ≡ d ln T/d ln P is the logarithmic temperature gradient and ∇ad is ∇ for an adiabatically stratified atmosphere. Flux rings with significantly larger field strength cannot be kept within the low latitude zones of the overshoot region.

Rempel et al. (2000) considered the mechanical equilibrium of a layer of an axisymmetric toroidal magnetic field of 105 G in a subadiabatically stratified region near the bottom of the solar convection zone in full spherical geometry. In this case, as illustrated in Figure 5(b), a latitudinal pressure gradient can be built up, allowing for force balance between a non-vanishing buoyancy force, the magnetic curvature force, and the pressure gradient without requiring a prograde toroidal flow. Thus a wider range of equilibria can exist. Rempel et al. (2000) found that under the condition of a strong subadiabatic stratification such as the radiative interior with δ ∼ −0.1, the magnetic layer tends to establish a mechanical equilibrium where a latitudinal pressure gradient is built up to balance the poleward component of the magnetic tension, and where the net radial component of the buoyancy and magnetic tension forces is efficiently balanced by the strong sub-adiabaticity. The magnetic layer reaches this equilibrium solution in a time scale short compared to the time required for a prograde toroidal flow to set up for the Coriolis force to be significant. For this type of equilibrium where a latitudinal pressure gradient is playing a dominant role in balancing the poleward component of the magnetic curvature force, there is significant relative density perturbation (≫ 1/β) in the magnetic layer compared to the background stratification. On the other hand, under the condition of a very weak subadiabatic stratification such as that in the overshoot layer near the bottom of the convection zone with δ ∼ −10−5, the magnetic layer tends to evolve towards a mechanical equilibrium which resembles that of an isolated toroidal flux ring, where the relative density perturbation is small (≪ 1/β), and the magnetic curvature force is balanced by the Coriolis force induced by a prograde toroidal flow in the magnetic layer. Thus regardless of whether the field is in the state of an extended magnetic layer or isolated flux tubes, a 105 G toroidal magnetic field stored in the weakly subadiabatically stratified overshoot region is preferably in a mechanical equilibrium with small relative density perturbation and with a prograde toroidal flow whose Coriolis force balances the magnetic tension. The prograde toroidal flow necessary for the equilibrium of the 105 G toroidal field is about 200 ms−1, which is approximately 10% of the mean rotation rate of the Sun. Thus one may expect significant changes in the differential rotation in the overshoot region during the solar cycle as the toroidal field is being amplified (Rempel et al., 2000). Detecting these toroidal flows and their temporal variation in the overshoot layer via helioseismic techniques is a means by which we can probe and measure the toroidal magnetic field generated by the solar cycle dynamo.

3.2 Effect of radiative heating

Storage of a strong super-equipartition field of 105 G at the base of the solar convection zone requires a state of mechanical equilibrium since convective motion is not strong enough to counteract the magnetic stress (Section 5.7). For isolated flux tubes stored in the weakly subadiabatic overshoot layer, the mechanical equilibrium corresponds to a neutrally buoyant state with a lower internal temperature (Section 3.1). Therefore flux tubes will be heated by radiative diffusion due to the mean temperature difference between the tube and the surrounding field-free plasma (see Parker, 1979; van Ballegooijen, 1982). Moreover, it is not adequate to just consider this zeroth order contribution due to the mean temperature difference in evaluating the radiative heat exchange between the flux tube and its surroundings. Due to the convective heat transport, the temperature gradient in the overshoot region and the lower convection zone is very close to being adiabatic, deviating significantly from that of a radiative equilibrium, and hence there is a non-zero divergence of radiative heat flux (see Spruit, 1974; van Ballegooijen, 1982). Thus an isolated magnetic flux tube with internally suppressed convective transport should also experience a net heating due to this non-zero divergence of radiative heat flux, provided that the radiative diffusion is approximately unaffected within the flux tube (Fan and Fisher, 1996; Moreno-Insertis et al., 2002; Rempel, 2003). In the limit of a thin flux tube, the rate of radiative heating (per unit volume) experienced by the tube is estimated to be (Fan and Fisher, 1996)

$${{dQ} \over {dt}} = - \nabla \cdot ({{\bf{F}}_{{\rm{rad}}}}) - \kappa (x_1^2/{a^2})(\overline T - {\overline T _{\rm{e}}})$$

where F rad is the unperturbed radiative energy flux, κ is the unperturbed radiative conductivity, x 1 is the first zero of the Bessel function J 0(x), a is the tube radius, \(\overline{T}\) is the mean temperature of the flux tube, and \(\overline{T}_{\rm e}\) is the corresponding unperturbed temperature at the location of the tube. Under the conditions prevailing near the base of the solar convection zone and for flux tubes that are responsible for active region formation, the first term due to the non-vanishing divergence of the radiative heat flux is found in general to dominate the second term. In the overshoot region, it can be shown that for these flux tubes the time scale for the heating to significantly increase their buoyancy from an initial neutrally buoyant state is long compared to the dynamic time scale characterized by the Brunt-Väisälä frequency. Thus the radiative heating is found to cause a quasi-static rise of the toroidal flux tubes, during which the tubes remain close to being neutrally buoyant. The upward drift velocity is estimated to be ∼ 10−3|δ|−1 cm s−1 which does not depend sensitively on the field strength of the flux tube (Fan and Fisher, 1996; Rempel, 2003). This implies that maintaining toroidal flux tubes in the overshoot region for a period comparable to the solar cycle time scale requires a strong subadiabaticity of δ < −10−4, which is significantly more subadiabatic than the values obtained by most of the overshoot models based on the non-local mixing length theory (see van Ballegooijen, 1982; Schmitt et al., 1984; Skaley and Stix, 1991).

On the other hand if the spatial filling factor of the toroidal flux tubes is large, or if the toroidal magnetic field is stored in the form of an extended magnetic layer, then the suppression of convective motion by the magnetic field is expected to alter the overall temperature stratification in the overshoot region. Rempel (2003) performed a 1D thermal diffusion calculation to model the change of the mean temperature stratification in the overshoot region when convective heat transport is being significantly suppressed. It is found that a reduction of the convective heat conductivity by a factor of 100 leads to the establishment of a new thermal equilibrium of significantly more stable temperature stratification with δ ∼ −10−4 in a time scale of a few months. Thus as the toroidal magnetic field is being amplified by the solar dynamo process, it may improve the conditions for its own storage by reducing the convective energy transport and increasing the subadiabaticity in the overshoot region.

4 Destabilization of a Toroidal Magnetic Field and Formation of Buoyant Flux Tubes

In the previous section, we have reviewed the equilibrium properties of a strong (∼ 105 G) toroidal magnetic field stored at the base of the solar convection zone. In this section we focus on the stability of the equilibria and the mechanisms by which the magnetic field can escape in the form of discrete buoyant flux tubes.

4.1 The buoyancy instability of isolated toroidal magnetic flux tubes

By linearizing the thin flux tube dynamic equations (1), (2), (5), (6), and (7), the stability of neutrally buoyant toroidal magnetic flux tubes to isentropic perturbations have been studied (see Spruit and van Ballegooijen, 1982a,b; Ferriz-Mas and Schüssler, 1993, 1995).

In the simplified case of a horizontal neutrally buoyant flux tube in a plane parallel atmosphere, ignoring the effects of curvature and solar rotation, the necessary and sufficient condition for instability is (Spruit and van Ballegooijen, 1982a,b)

$${k^2}H_p^2 < {{\beta/2} \over {1 + \beta }}(1/\gamma + \beta \delta),$$
(17)

where k is the wavenumber along the tube of the undulatory perturbation, H p is the local pressure scale height, βp/(B 2/8π) is the ratio of the plasma pressure divided by the magnetic pressure of the flux tube, δ = ∇ − ∇ad is the superadiabaticity, and γ is the ratio of the specific heats. If all values of k are allowed, then the condition for the presence of instability is

$$\beta \delta > - 1/\gamma.$$
(18)

Note that k → 0 is a singular limit. For perturbations with k = 0 which do not involve bending the field lines, the condition for instability becomes (Spruit and van Ballegooijen, 1982a)

$$\beta \delta > {2 \over \gamma }\left({{1 \over \gamma } - {1 \over 2}} \right)\sim 0.12$$
(19)

which is a significantly more stringent condition than (18), even more stringent than the convective instability for a field-free fluid (δ > 0). Thus the undulatory instability (with k ≠ 0) is of a very different nature and is easier to develop than the instability associated with uniform up-and-down motions of the entire flux tube. The undulatory instability can develop even in a convectively stable stratification with δ < 0 as long as the field strength of the flux tube is sufficiently strong (i.e. β is of sufficiently small amplitude) such that |βδ| is smaller than 1/γ. In the regime of −1 < βδ< (2γ)(1/γ − 1/2) where only the undulatory modes with k ≠ 0 are unstable, a longitudinal flow from the crests to the troughs of the undulation is essential for driving the instability. Since the flux tube has a lower internal temperature and hence a smaller pressure scale height inside, upon bending the tube, matter will flow from the crests to the troughs to establish hydrostatic equilibrium along the field. This increases the buoyancy of the crests and destabilizes the tube (Spruit and van Ballegooijen, 1982a).

Including the curvature effect of spherical geometry, but still ignoring solar rotation, Spruit and van Ballegooijen (1982a,b) have also studied the special case of a toroidal flux ring in mechanical equilibrium within the equatorial plane. Since the Coriolis force due to solar rotation is ignored, the flux ring in the equatorial plane needs to be slightly buoyant to balance the inward tension force. For latitudinal motions out of the equatorial plane, the axisymmetric component is unstable, which corresponds to the poleward slip of the tube as a whole. But this instability can be suppressed when the Coriolis force is included (Ferriz-Mas and Schüssler, 1993). For motions within the equatorial plane, the conditions for instabilities are (Spruit and van Ballegooijen, 1982a,b)

$$\matrix{{{1 \over 2}\beta \delta > ({m^2} - 3 - s){f^2} + 2f/\gamma - 1(2\gamma)} \hfill \& {(m \geq 1),} \hfill \cr {{1 \over 2}\beta \delta > {f^2}(1 - s) - 2f/\gamma + {1 \over \gamma }\left({{1 \over \gamma } - {1 \over 2}} \right)} \hfill \& {(m = 0)} \hfill \cr}$$
(20)

where fH p/r 0 is the ratio of the pressure scale height over the radius of the bottom of the solar convection zone, m (having integer values 0, 1, …) denotes the azimuthal order of the undulatory mode of the closed toroidal flux ring, i.e. the wavenumber k = m/r 0, s is a parameter that describes the variation of the gravitational acceleration: gr s. Near the base of the solar convection zone, f ∼ 0.1, s ∼ −2. Thus conditions (20) show that it is possible for m = 0, 1, 2, 3, 4 modes to become unstable in the weakly subadiabatic overshoot region, and that the instabilities of m = 1, 2, 3 modes require less stringent conditions than the instability of m = 0 mode. Since Equation (20) is derived for the singular case of an equilibrium toroidal ring in the equatorial plane, its applicability is very limited.

The general problem of the linear stability of a thin toroidal flux ring in mechanical equilibrium in a differentially rotating spherical convection zone at arbitrary latitudes has been studied in detail by Ferriz-Mas and Schüssler (1993, 1995). For general non-axisymmetric perturbations, a sixth-order dispersion relation is obtained from the linearized thin flux tube equations. It is not possible to obtain analytical stability criteria. The dispersion relation is solved numerically to find instability and the growth rates of the unstable modes. The regions of instability in the (B 0, λ 0) plane (with B 0 being the magnetic field strength of the flux ring and λ 0 being the equilibrium latitude), under the conditions representative of the overshoot layer at the base of the solar convection zone are shown in Figure 6 (from Caligari et al., 1995).

Figure 6:
figure 6

Upper panel: Regions of unstable toroidal flux tubes in the (B 0, λ 0)-plane (with B 0 being the magnetic field strength of the flux tubes and λ 0 being the equilibrium latitude). The subadiabaticity at the location of the toroidal flux tubes is assumed to be δ ≡ ∇ − ∇ad = −2.6 × 10−6. The white area corresponds to a stable region while the shaded regions indicate instability. The degree of shading signifies the azimuthal wavenumber of the most unstable mode. The contours correspond to lines of constant growth time of the instability. Thicker lines are drawn for growth times of 100 days and 300 days. Lower panel: Same as the upper panel except that the subadiabaticity at the location of the toroidal tubes is δ ≡ ∇ − ∇ad = −1.9 × 10−7. From Caligari et al. (1995).

The basic parameters that determine the stability of an equilibrium toroidal flux ring are its field strength and the subadiabaticity of the external stratification. In the case δ ≡ ∇ − ∇ad = −2.6 × 10−6 (upper panel of Figure 6), unstable modes with reasonably short growth times (less than about a year) only begin to appear at sunspot latitudes for B 0 ≳ 1.2 × 105 G. These unstable modes are of m = 1 and 2. In case of a weaker subadiabaticity, δ ≡ ∇ − ∇adad ≡ −1.9 × 10−7 (lower panel of Figure 6), reasonably fast growing modes (growth time less than a year) begin to appear at sunspot latitudes for B 0 ≳ 5 × 104 G, and the most unstable modes are of m = 1 and 2. These results suggest that toroidal magnetic fields stored in the overshoot layer at the base of the solar convection zone do not become unstable until their field strength becomes significantly greater than the equipartition value of 104 G.

Thin flux tube simulations of the non-linear growth of the non-axisymmetric instabilities of initially toroidal flux tubes and the emergence of Ω-shaped flux loops through the solar convective envelope will be discussed in Section 5.1.

4.2 Breakup of an equilibrium magnetic layer and formation of buoyant flux tubes

It is possible that the toroidal magnetic field stored at the base of the convection zone is in the form of an extended magnetic layer, instead of individual magnetic flux tubes for which the thin flux tube approximation can be applied. The classic problem of the buoyancy instability of a horizontal magnetic field \({\bf B}=B(z)\hat{\bf x}\) in a plane-parallel, gravitationally stratified atmosphere with a constant gravity −g, pressure p(z), and density ρ(z), in hydrostatic equilibrium,

$${d \over {dz}}\left({p + {{{B^2}} \over {8\pi }}} \right) = - \rho g,$$
(21)

has been studied by many authors in a broad range of astrophysics contexts including

  • magnetic fields in stellar convection zones (see Newcomb, 1961; Parker, 1979; Hughes and Cattaneo, 1987),

  • magnetic flux emergence into the solar atmosphere (see Shibata et al., 1989),

  • stability of prominence support by a magnetic field (see Zweibel and Bruhwiler, 1992),

  • and the instability of the interstellar gas and magnetic field (see Parker, 1966).

The linear stability analysis of the above equilibrium horizontal magnetic layer (Newcomb, 1961) showed that the necessary and sufficient condition for the onset of the general 3D instability with non-zero wavenumbers (k x ≠ 0, k y ≠ 0) in both horizontal directions parallel and perpendicular to the magnetic field is that

$${{d\rho } \over {dz}} > - {{{\rho ^2}g} \over {\gamma p}},$$
(22)

is satisfied somewhere in the stratified fluid. On the other hand the necessary and sufficient condition for instability of the purely interchange modes (with k x = 0 and k y ≠ 0) is that

$${{d\rho } \over {dz}} > - {{{\rho ^2}g} \over {\gamma p + {B^2}/4\pi }}.$$
(23)

is satisfied somewhere in the fluid — a more stringent condition than (22). Note in Equations (22) and (23), p and ρ are the plasma pressure and density in the presence of the magnetic field. Hence the effect of the magnetic field on the instability criteria is implicitly included. As shown by Thomas and Nye (1975) and Acheson (1979), the instability conditions (22) and (23) can be alternatively written as

$${{v_{\rm{a}}^2} \over {c_{\rm{s}}^2}}{{d\ln B} \over {dz}} < - {1 \over {{c_p}}}{{ds} \over {dz}}$$
(24)

for instability of general 3D undulatory modes and

$${{v_{\rm{a}}^2} \over {c_{\rm{s}}^2}}{d \over {dz}}\left[ {\ln \left({{B \over \rho }} \right)} \right] < - {1 \over {{c_p}}}{{ds} \over {dz}}$$
(25)

for instability of purely 2D interchange modes, where v a is the Alfvén speed, c s is the sound speed, c p is the specific heat under constant pressure, and ds/dz is the actual entropy gradient in the presence of the magnetic field. The development of these buoyancy instabilities is driven by the gravitational potential energy that is made available by the magnetic pressure support. For example, the magnetic pressure gradient can “puff-up” the density stratification in the atmosphere, making it decrease less steeply with height (causing condition (22) to be met), or even making it top heavy. This raises the gravitational potential energy and makes the atmosphere unstable. In another situation, the presence of the magnetic pressure can support a layer of cooler plasma with locally reduced temperature embedded in an otherwise stably stratified fluid. This can also cause the instability condition (22) to be met locally in the magnetic layer. In this case the pressure scale height within the cooler magnetic layer is smaller, and upon bending the field lines, plasma will flow from the crests to the troughs to establish hydrostatic equilibrium, thereby releasing gravitational potential energy and driving the instability. This situation is very similar to the buoyancy instability associated with the neutrally buoyant magnetic flux tubes discussed in Section 4.1.

The above discussion on the buoyancy instabilities considers ideal adiabatic perturbations. It should be noted that the role of finite diffusion is not always stabilizing. In the solar interior, it is expected that ηK and νK, where η, ν, and K denote the magnetic diffusivity, the kinematic viscosity, and the thermal diffusivity respectively. Under these circumstances, it is shown that thermal diffusion can be destabilizing (see Gilman, 1970; Acheson, 1979; Schmitt and Rosner, 1983). The diffusive effects are shown to alter the stability criteria of Equations (24) and (25) by reducing the term ds/dz by a factor of η/K (see Acheson, 1979). In other words, efficient heat exchange can significantly “erode away” the stabilizing effect of a subadiabatic stratification. This process is an example of the double-diffusive instabilities.

Direct multi-dimensional MHD simulations have been carried out to study the break-up of a horizontal magnetic layer by the non-linear evolution of the buoyancy instabilities and the formation of buoyant magnetic flux tubes (see Cattaneo and Hughes, 1988; Cattaneo et al., 1990; Matthews et al., 1995; Wissink et al., 2000; Fan, 2001a).

Cattaneo and Hughes (1988), Matthews et al. (1995), and Wissink et al. (2000) have carried out a series of 2D and 3D compressible MHD simulations where they considered an initial horizontal magnetic layer that supports a top-heavy density gradient, i.e. an equilibrium with a lower density magnetic layer supporting a denser plasma on top of it. It is found that for this equilibrium configuration, the most unstable modes are the Rayleigh-Taylor type 2D interchange modes. Two-dimensional simulations of the non-linear growth of the interchange modes (Cattaneo and Hughes, 1988) found that the formation of buoyant flux tubes is accompanied by the development of strong vortices whose interactions rapidly destroy the coherence of the flux tubes. In the non-linear regime, the evolution is dominated by vortex interactions which act to prevent the rise of the buoyant magnetic field. Matthews et al. (1995) and Wissink et al. (2000) extend the simulations of Cattaneo and Hughes (1988) to 3D allowing variations in the direction of the initial magnetic field. They discovered that the flux tubes formed by the initial growth of the 2D interchange modes subsequently become unstable to a 3D undulatory motion in the non-linear regime due to the interaction between neighboring counter-rotating vortex tubes, and consequently the flux tubes become arched. Matthews et al. (1995) and Wissink et al. (2000) pointed out that this secondary undulatory instability found in the simulations is of similar nature as the undulatory instability of a pair of counter-rotating (non-magnetic) line vortices investigated by Crow (1970). Wissink et al. (2000) further considered the effect of the Coriolis force due to solar rotation using a local f-plane approximation, and found that the principle effect of the Coriolis force is to suppress the instability. Further 2D simulations have also been carried out by Cattaneo et al. (1990) where they introduced a variation of the magnetic field direction with height into the previously unidirectional magnetic layer of Cattaneo and Hughes (1988). The growth of the interchange instability of such a sheared magnetic layer results in the formation of twisted, buoyant flux tubes which are able to inhibit the development of vortex tubes and rise cohesively.

On the other hand, Fan (2001a) has considered a different initial equilibrium state for a horizontal unidirectional magnetic layer, where the density stratification remains unchanged from that of an adiabatically stratified polytrope, but the temperature and the gas pressure are lowered in the magnetic layer to satisfy the hydrostatic condition. For such a neutrally buoyant state with no density change inside the magnetic layer, the 2D interchange instability is completely suppressed and only 3D undulatory modes (with non-zero wavenumbers in the field direction) are unstable. The strong toroidal magnetic field stored in the weakly subadiabatic overshoot region below the bottom of the convection zone is likely to be close to such a neutrally buoyant mechanical equilibrium state (see Section 3.1). Anelastic MHD simulations (Fan, 2001a) of the growth of the 3D undulatory instability of this horizontal magnetic layer show formation of significantly arched magnetic flux tubes (see Figure 7) whose apices become increasingly buoyant as a result of the diverging flow of plasma from the apices to the troughs.

Figure 7:
figure 7

mpg-Movie (249 KB) Still from a movie showing The formation of arched flux tubes as a result of the non-linear growth of the undulatory buoyancy instability of a neutrally buoyant equilibrium magnetic layer perturbed by a localized velocity field. From Fan (2001a). The images show the volume rendering of the absolute magnetic field strength |B|. Only one half of the wave length of the undulating flux tubes is shown, and the left and right columns of images show, respectively, the 3D evolution as viewed from two different angles (For video see appendix).

The decrease of the field strength B at the apex of the arched flux tube as a function of height is found to follow approximately the relation \(B/\sqrt{\rho}={\rm constant}\), or, the Alfvén speed being constant, which is a significantly slower decrease of B with height compared to that for the rise of a horizontal flux tube without any field line stretching, for which case B/ρ should remain constant. The variation of the apex field strength with height following \(B/\sqrt{\rho}={\rm constant}\) found in the 3D MHD simulations of the arched flux tubes is in good agreement with the results of the thin flux tube models of emerging Ω-loops (see Moreno-Insertis, 1992) during their rise through the lower half of the solar convective envelope where the stratification is very close to being adiabatic as is assumed in the 3D simulations.

Kersalé et al. (2007) studied the nonlinear 3D evolution of the magnetic buoyancy instability resulting from a smoothly stratified horizontal magnetic field, and with the instability continually driven via the boundary conditions. They considered the case where the prescribed magnetic pressure gradient is such that the equilibrium is unstable to the 3D modes but stable to 2D interchange modes. One important distinction of this work compared to many of the previous studies is that the instability is continually driven through imposing a fixed magnetic pressure gradient at the top and bottom boundaries (Figure 8) which are stress-free and impermeable.

Figure 8:
figure 8

Horizontal average of the magnetic field B x as a function depth for the initial state (dotted line), at a later time when the instability saturates (dashed line), and in the final steady state (solid line). The magnetic pressure gradient is maintained at the top and bottom boundaries during the non-linear evolution of the magnetic buoyancy instability. From Kersalé et al. (2007). Figure reproduced by permission of the AAS.

The initial growth of the instabilities from a random perturbation results in the formation of arched flux tubes. In the non-linear stage, the system is found to establish a modulated periodic state where discrete flux tube concentrations with field strength significantly stronger than the initial mean field form periodically as modulated traveling waves (see Figures 9 and 10). The development of isolated flux tube concentrations results from convergent downflows continually driven by the instability (Figure 10). This result provides an interesting mechanism for the formation of strong active region flux tubes from dynamo generated large scale field at the base of the convection zone.

Figure 9:
figure 9

Evolution of the kinetic energy density. The system eventually establishes a modulated periodic state with two disparate time scales. From Kersalé et al. (2007). Figure reproduced by permission of the AAS.

Figure 10:
figure 10

Space-time plots for fixed values of x and z of the magnetic energy density (left), the transverse horizontal (y) velocity (middle), and the vertical velocity (right), with the horizontal axis being the y-axis and vertical axis denoting the time. From Kersalé et al. (2007). Figure reproduced by permission of the AAS.

4.3 Buoyancy breakup of a shear-generated magnetic layer

Instead of prescribing an unstable equilibrium of an initial magnetic flux tube or layer, Vasil and Brummell (2008) carried out a series of 3D MHD simulations of the generation of a strong layer of horizontal magnetic field by the action of a vertical shear on a weak vertical field in a subadiabatically stratified atmosphere, and examine the subsequent breakup of the resulting magnetic configuration via magnetic buoyancy instabilities (see Figure 11).

Figure 11:
figure 11

From Vasil and Brummell (2008). A 3D MHD simulation of the build up and subsequent buoyancy break up of a layer of horizontal magnetic field forced by a vertical shear on an initially weak vertical field in a subadiabatically stratified atmosphere. The sequence of images show the volume renderings of the magnetic field strength. Figure reproduced by permission of the AAS.

The aim of these simulations is to examine under what conditions the radial shear of differential rotation operating in the thin solar tachocline layer can amplify a strong enough large scale toroidal magnetic field that undergoes magnetic buoyancy instabilities and develops buoyantly rising structures. The numerical simulations together with a subsequent analytical study (Vasil and Brummell, 2009) show that magnetic buoyancy instabilities can indeed develop in the shear-generated magnetic layer (Figure 11) if the forcing that drives the shear flow is sufficiently large. The needed forcing is such that, in the absence of the magnetic field, it imposes a hydrodynamically unstable shear. It is found that the imposed shear needs to have a Richardson number R i being less than 1, where R i measures the relative importance of the stabilizing effect of the stratification over the strength of the shear to overturn the fluid (Vasil and Brummell, 2009; Silvers et al., 2009). This result is not surprising because in order for the magnetic layer to be buoyantly unstable, the imposed shear flow needs to transfer enough energy to the magnetic field for it to overcome the stable background stratification (Silvers et al., 2009). It is not clear whether such strong forcing of the shear exists in the solar tachocline. For the observed shear in the solar tachocline, the Richardson number is estimated to be much greater than 1, R i ∼ 103 ∼ 105 (Gough, 2007). However, the observed shear in the tachocline may not correspond to the forcing shear, but is the end steady-state reached when the forcing is balanced by the built-up magnetic stress and turbulent transport.

Silvers et al. (2009) further extend the studies of Vasil and Brummell (2008) and Vasil and Brummell (2009) by considering the fact that the ratio of the magnetic diffusivity (η) over the thermal diffusivity (κ) in the solar tachocline is very small: ξ = η/κ ≪ 1. Under such conditions the double-diffusive magnetic buoyancy instabilities can develop at a much less steep magnetic pressure gradient for the magnetic layer compared to that required for magnetic buoyancy instabilities under the assumption of adiabatic evolution (see end of Section 4.2). The stabilizing effect of the subadiabatic stratification is significantly reduced by the thermal diffusion. Simulations by Silvers et al. (2009) verify that double-diffusive magnetic buoyancy instabilities indeed can develop for a magnetic layer generated by a weak forcing shear that is hydrodynamically stable (with R i < 2.96).

5 Dynamic Evolution of Emerging Flux Tubes in the Solar Convection Zone

5.1 Results from thin flux tube simulations of emerging loops

Beginning with the seminal work of Moreno-Insertis (1986) and Choudhuri and Gilman (1987), a large body of numerical simulations solving the thin flux tube dynamic equations (1), (2), (5), (6), and (7) — or various simplified versions of them — have been carried out to model the evolution of emerging magnetic flux tubes in the solar convective envelope (see Choudhuri, 1989; D’Silva and Choudhuri, 1993; Fan et al., 1993, 1994; Schüssler et al., 1994; Caligari et al., 1995; Fan and Fisher, 1996; Caligari et al., 1998; Fan and Gong, 2000). The results of these numerical calculations have contributed greatly to our understanding of the basic properties of solar active regions and provided constraints on the field strengths of the toroidal magnetic fields at the base of the solar convection zone.

Most of the earlier calculations (see Choudhuri and Gilman, 1987; Choudhuri, 1989; D’Silva and Choudhuri, 1993; Fan et al., 1993, 1994) considered initially buoyant toroidal flux tubes by assuming that they are in temperature equilibrium with the external plasma. Various types of initial undulatory displacements are imposed on the buoyant tube so that portions of the tube will remain anchored within the stably stratified overshoot layer and other portions of the tube are displaced into the unstable convection zone which subsequently develop into emerging Ω-shaped loops.

Later calculations (see Schüssler et al., 1994; Caligari et al., 1995, 1998; Fan and Gong, 2000) considered more physically self-consistent initial conditions where the initial toroidal flux ring is in the state of mechanical equilibrium. In this state the buoyancy force is zero (neutrally buoyant) and the magnetic curvature force is balanced by the Coriolis force resulting from a prograde toroidal motion of the tube plasma. It is argued that this mechanical equilibrium state is the preferred state for the long-term storage of a toroidal magnetic field in the stably stratified overshoot region (Section 3.1). In these simulations, the development of the emerging Ω-loops is obtained naturally by the non-linear, adiabatic growth of the undulatory buoyancy instability associated with the initial equilibrium toroidal flux rings (Section 4.1). As a result there is far less degree of freedom in specifying the initial perturbations. The eruption pattern needs not be prescribed in an ad hoc fashion but is self-consistently determined by the growth of the instability once the initial field strength, latitude, and the subadiabaticity at the depth of the tube are given. For example Caligari et al. (1995) modeled emerging loops developed due to the undulatory buoyancy instability of initial toroidal flux tubes located at different depths near the base of their model solar convection zone which includes a consistently calculated overshoot layer according to the non-local mixing-length treatment. They choose values of initial field strengths and latitudes that lie along the contours of constant instability growth times of 100 days and 300 days in the instability diagrams (see Figure 6), given the subadiabaticity at the depth of the initial tubes. The tubes are then perturbed with a small undulatory displacement which consists of a random superposition of Fourier modes with azimuthal order ranging from m = 1 through m = 5, and the resulting eruption pattern is determined naturally by the growth of the instability.

On the other hand, non-adiabatic effects may also be important in the destabilization process. It has been discussed in Section 3.2 that isolated magnetic flux tubes with internally suppressed convective transport experience a net heating due to the non-zero divergence of radiative heat flux in the weakly subadiabatically stratified overshoot region and also in the lower solar convection zone. The radiative heating causes a quasi-static upward drift of the toroidal flux tube with a drift velocity ∼ 10−3|δ|−1 cm s−1. Thus the time scale for a toroidal flux tube to drift out of the stable overshoot region may not be long compared to the growth time of its undulatory buoyancy instability. For example if the subadiabaticity δ is ∼ −10−6, the time scale for the flux tube to drift across the depth of the overshoot region is about 20 days, smaller than the growth times (∼ 100–300 days) of the most unstable modes for tubes of a ∼ 105 G field as shown in Figure 6. Therefore radiative heating may play an important role in destabilizing the toroidal flux tubes. The quasi-static upward drift due to radiative heating can speed-up the development of emerging Ω-loops (especially for weaker flux tubes) by bringing the tube out of the inner part of the overshoot region of stronger subadiabaticity, where the tube is stable or the instability growth is very slow, to the outer overshoot region of weaker subadiabaticity or even into the convection zone, where the growth of the undulatory buoyancy instability occurs at a much shorter time scale.

A possible scenario in which the effect of radiative heating helps to induce the formation of Ω-shaped emerging loops has been investigated by Fan and Fisher (1996). In this scenario, the initial neutrally buoyant toroidal flux tube is not exactly uniform, and lies at non-uniform depths with some portions of the tube lying at slightly shallower depths in the overshoot region. Radiative heating and quasi-static upward drift of this non-uniform flux tube bring the upward protruding portions of the tube first into the unstably stratified convection zone. These portions can become buoyantly unstable (if the growth of buoyancy overcomes the growth of tension) and rise dynamically as emerging loops. In this case the non-uniform flux tube remains close to a mechanical equilibrium state during the initial quasi-static rise through the overshoot region. The emerging loop develops gradually as a result of radiative heating and the subsequent buoyancy instability of the outer portion of the tube entering the convection zone.

In the following subsections we review the major findings and conclusions that have been drawn from the various thin flux tube simulations of emerging flux loops.

5.1.1 Latitude of flux emergence

Axisymmetric simulations of the buoyant rise of toroidal flux rings in a rotating solar convective envelope by Choudhuri and Gilman (1987) first demonstrate the significant influence of the Coriolis force on the rising trajectories. The basic effect is that the Coriolis force acting on the radial outward motion of the flux tube (or the tendency for the rising tube to conserve angular momentum) drives a retrograde motion of the tube plasma. This retrograde motion then induces a Coriolis force directed towards the Sun’s rotation axis which acts to deflect the trajectory of the rising tube poleward. The amount of poleward deflection by the Coriolis force depends on the initial field strength of the emerging tube, being larger for flux tubes with weaker initial field. For flux tubes with an equipartition field strength of B ∼ 104 G, the effect of the Coriolis force is so dominating that it deflects the rising tubes to emerge at latitudes poleward of the sunspot zones even though the flux tubes start out from low latitudes at the bottom of the convective envelope. In order for the rising trajectory of the flux ring to be close to radial so that the emerging latitudes are within the observed sunspot latitudes, the field strength of the toroidal flux ring at the bottom of the solar convection zone needs to be ∼ 105 G. This basic result is confirmed by later simulations of non-axisymmetric, Ω-shaped emerging loops rising through the solar convective envelope (see Choudhuri, 1989; D’Silva and Choudhuri, 1993; Fan et al., 1993; Schüssler et al., 1994; Caligari et al., 1995, 1998; Fan and Fisher, 1996).

Simulations by Caligari et al. (1995) of emerging loops developed self-consistently due to the undulatory buoyancy instability show that, for tubes with initial field strength ≳ 105 G, the trajectories of the emerging loops are primarily radial with poleward deflection no greater than 3°. For tubes with initial field strength exceeding 4 ×104 G, the poleward deflection of the emerging loops remain reasonably small (no greater than about 6°). However, for a tube with equipartition field strength of 104 G, the rising trajectory of the emerging loop is deflected poleward by about 20°. Such an amount of poleward deflection is too great to explain the observed low latitudes of active region emergence. Furthermore, it is found that with such a weak initial field the field strength of the emerging loop falls below equipartition with convection throughout most of the convection zone. Such emerging loops are expected to be subjected to strong deformation by turbulent convection and may not be consistent with the observed well defined order of solar active regions.

Fan and Fisher (1996) modeled emerging loops that develop as a result of radiative heating of non-uniform flux tubes in the overshoot region. The results on the poleward deflection of the emerging loops as a function of the initial field strength are very similar to that found in Caligari et al. (1995). Figure 12 shows the latitude of loop emergence as a function of the initial latitude at the base of the solar convection zone. It can be seen that tubes of 105 G emerge essentially radially with very small poleward deflection (< 3°), and for tubes with B ≳ 3 × 104 G, the poleward deflections remain reasonably small so that the emerging latitudes are within the observed sunspot zones.

Figure 12:
figure 12

Latitude of loop emergence as a function of the initial latitude at the base of the solar convection zone, for tubes with initial field strengths B = 30, 60, and 100 kG and fluxes Φ = 1021 and 1022 Mx. From Fan and Fisher (1996).

5.1.2 Active region tilts

A well-known property of the solar active regions is the so called Joy’s law of active region tilts. The averaged orientation of bipolar active regions on the solar surface is not exactly toroidal but is slightly tilted away from the east-west direction, with the leading polarity (the polarity leading in the direction of rotation) being slightly closer to the equator than the following polarity. The mean tilt angle is a function of latitude, being approximately ∝ sin(latitude) (Wang and Sheeley Jr, 1989, 1991; Howard, 1991a,b; Fisher et al., 1995; Kosovichev and Stenflo, 2008).

Using thin flux tube simulations of the non-axisymmetric eruption of buoyant Ω-loops in a rotating solar convective envelope, D’Silva and Choudhuri (1993) were the first to show that the active region tilts as described by Joy’s law can be explained by Coriolis forces acting on the flux loops. As the emerging loop rises, there is a relative expanding motion of the mass elements at the summit of the loop. The Coriolis force induced by this diverging, expanding motion at the summit is to tilt the summit clockwise (counter-clockwise) for loops in the northern (southern) hemisphere as viewed from the top, so that the leading side from the summit is tilted equatorward relative to the following side. Since the component of the Coriolis force that drives this tilting has a sin(latitude) dependence, the resulting tilt angle at the apex is approximately ∝ sin(latitude).

Caligari et al. (1995) studied tilt angles of emerging loops developed self-consistently due to the undulatory buoyancy instability of flux tubes located at the bottom as well as just above the top of their model overshoot region, with selected values of initial field strengths and latitudes lying along contours of constant instability growth times (100 days and 300 days). The resulting tilt angles at the apex of the emerging loops (see Figure 13) produced by these sets of unstable tubes (whose field strengths are within the range of 4 × 104 G to 1.5 × 105 G) show good agreement with the observed tilt angles for sunspot groups measured by Howard (1991b).

Figure 13:
figure 13

Tilt angles at the apex of the emerging flux loops as a function of the emergence latitudes. The squares and the asterisks denote loops originating from initial toroidal tubes located at different depths with different local subadiabaticity (squares: δ ≡ ∇ − ∇ad = −2.6 × 10−6 and field strength ranges between 105 G and 1.5 × 105 G; asterisks: δ ≡ ∇ − ∇ad = −1.9 − 10−7 and field strength ranges between 4 × 104 G and 6 × 104 G). The shaded region indicates the range of the observed tilt angles of sunspot groups measured by Howard (1991b). From Caligari et al. (1995).

They also found that loops formed from toroidal flux tubes with an equipartition field strength of 104 G develop a tilt angle (−9°) of the wrong sign at the loop apex.

Similar results of loop tilt angles (see Figure 14) are found by Fan and Fisher (1996) who considered formation of emerging loops by gradual radiative heating of non-uniform toroidal flux tubes initially in mechanical equilibrium in the overshoot region.

Figure 14:
figure 14

Tilt angles at the apex of the emerging loops versus emerging latitudes resulting from thin flux tube simulations of Fan and Fisher (1996). The numbers used as data points indicate the corresponding initial field strength values in units of 10 kG (‘X’s represent 100 kG). Also plotted (solid line) is the least squares fit: tilt angle = 15.7° × sin(latitude), obtained in Fisher et al. (1995) by fitting to the measured tilt angles of 24701 sunspot groups observed at Mt. Wilson, same data set studied in Howard (1991b).

It is found that emerging loops with initial field strength of the range 4 × 104 G to 105 G show tilt angles that are in good agreement with the observed tilt angles of sunspot groups. Tilts of the wrong sign or direction begin to appear for tubes with initial field strength ≲ 3 × 104 G.

Thin flux tube simulations also show that the tilt angles of the emerging loops tend to be smaller for tubes with smaller total flux or radius (see D’Silva and Choudhuri, 1993; Fan and Fisher, 1996). This is because smaller tubes are more influenced by the drag force, whose direct effect is to oppose the tilting motion of the loops. This prediction on the dependence of the tilt angle on active region flux is confirmed by Fisher et al. (1995), who studied the relation between the tilt angle and the mean separation distance between the leader and follower spots of a sunspot group. Wang and Sheeley Jr (1989) and Howard (1992) have shown that the mean polarity separation distance of an active region is a good proxy for the total magnetic flux in the region. Fisher et al. (1995) found that for a fixed latitude, the tilt angle of sunspot groups decreases with decreasing polarity separation, and hence decreasing total flux, consistent with the results from the thin flux tube calculations. This agreement adds support to the explanation that the Coriolis force acting on rising flux loops is the main cause of active region tilts and argues against the suggestion by Babcock (1961) that the active region tilt angles simply reflect the orientation of the underlying toroidal magnetic field stretched out by the latitudinal differential rotation.

Using Mount Wilson sunspot group data, Fisher et al. (1995) further studied the dispersion or scatter of spot-group tilts away from the mean tilt behavior described by Joy’s law. First they found that the magnitude of the tilt dispersion is significantly greater than the level expected from measurement errors, suggesting a solar origin of the tilt angle scatter. Furthermore, they found that the root-mean-square tilt scatter decreases with increasing polarity separation (or total flux) and does not vary with latitude (in contrast to the latitudinal dependence of the mean tilts). This result is consistent with the picture that scatter of active region tilts away from the mean Joy’s law behavior results from buffeting of emerging loops by convective motions during their rise through the solar convection zone (Longcope and Fisher, 1996).

Non-linear simulations of the two-dimensional MHD tachocline (Cally et al., 2003) show that bands of toroidal magnetic fields in the solar tachocline may become tipped relative to the azimuthal direction by an amount that is within +/− 10° at sunspot latitudes due to the non-linear evolution of the 2D global joint instability of differential rotation and toroidal magnetic fields. This tipping may either enhance or reduce the observed tilt in bipolar active regions depending on from which part of the tipped band the emerging loops develop. Thus the basic consequence of the possible tipping of the toroidal magnetic fields in the tachocline is to contribute to the spread of the tilts of bipolar active regions.

More recently, using a series of 96 minute cadence magnetograms from SOHO MDI and analyzing 715 bipolar magnetic regions which emerged within 30° from the central meridian and outside already existing active regions, Kosovichev and Stenflo (2008) investigated how the active region tilt angle evolves during flux emergence and how it correlates with other properties of the emerging region. The study shows that at the beginning of emergence the tilt angles are random, and the mean tilt angle is about zero (see Figure 15(a)).

Figure 15:
figure 15

The distribution of the tilt angle as a function of sine latitude: (a) at the beginning of flux emergence, (b) at the middle of the emergence period, and (c) at the end of emergence. From Kosovichev and Stenflo (2008). Figure reproduced with permission of the AAS.

However by the middle of the emergence period (flux growth period), the tilt angles clearly show a systematic mean as a function of latitude that follows Joy’s law (Figure 15(b)). At the end of the emergence period, the Joy’s law dependence has become more pronounced as the scatter from the systematic mean tilt decreases (Figure 15(c)). The above result that the systematic mean tilt following Joy’s law is established during the flux emergence period (flux growth period) suggests that the tilt of the emerging flux tube has developed in the interior before reaching the surface. This is consistent with the above model of rising flux tubes where the tilt angle is caused by the effect of the Coriolis force during the rise. However, Kosovichev and Stenflo (2008) also found that the tilt angle does not show a systematic dependence on the flux of the active region, in contradiction to the expectation of the rising flux tube model (e.g. Fisher et al., 1995). Furthermore, Kosovichev and Stenflo (2008) found that there is no tendency for the systematic active region mean tilt to relax towards the east-west direction after the emergence has ceased and the driving Coriolis force has vanished, at which time the tension of the flux tube is expected to act to restore the original toroidal orientation of the tube at the base of the solar convection zone. The latter result may be understood if the active region magnetic fields on the photosphere become dynamically disconnected from the interior flux tubes soon after emergence (Section 8.3). On the other hand, as suggested in Kosovichev and Stenflo (2008), it may be that Joy’s law of solar active regions reflects not the Coriolis effect of the rising flux tubes but the spiral orientation of the nearly toroidal magnetic field lines in the interior generated by the latitudinal differential rotation (Babcock, 1961).

5.1.3 Morphological asymmetries of active regions

An intriguing property of solar active regions is the asymmetry in morphology between the leading and following polarities. The leading polarity of an active region tends to be in the form of large sunspots, whereas the following polarity tends to appear more dispersed and fragmented; moreover, the leading spots tend to be longer lived than the following. Fan et al. (1993) offered an explanation for the origin of this asymmetry. In their thin flux tube simulations of the non-axisymmetric eruption of buoyant Ω-loops through a rotating model solar convective envelope, they found that an asymmetry in the magnetic field strength develops between the leading and following legs of an emerging loop, with the field strength of the leading leg being about 2 times that of the following leg. The field strength asymmetry develops because the Coriolis force, or the tendency for the tube plasma to conserve angular momentum, drives a counter-rotating flow of plasma along the emerging loop, which, in conjunction with the diverging flow of plasma from the apex to the troughs, gives rise to an effective asymmetric stretching of the two legs of the loop with a greater stretching and hence a stronger field strength in the leading leg. Fan et al. (1993) argued that the stronger field along the leading leg of the emerging loop makes it less subject to deformation by the turbulent convection and therefore explains the more coherent and less fragmented appearance of the leading polarity of an active region.

However, the calculations by Caligari et al. (1995) and Caligari et al. (1998) show that the field strength asymmetry found in Fan et al. (1993) depends on the choices of the initial conditions. They found that if one uses the more self-consistent mechanical equilibrium state (instead of the buoyant, temperature equilibrium state used in Fan et al. (1993)) for the initial toroidal tube, for which there exists an initial prograde toroidal motion, the subsequent differential stretching of the emerging loop is quite different. It is found that a consistently stronger field along the leading leg of the emerging loop only occurs for cases with relatively weak initial field strengths (∼ 104 G). For stronger fields (∼ 105 G), which seem to fit the observed properties of solar active regions better in all the other aspects, e.g. the latitude of flux emergence and the tilt angles, the field strength asymmetry becomes very small and even reverses its sense in the upper half of the solar convection zone.

Further work by Fan and Fisher (1996) examined the field strength asymmetry for emerging flux loops that form as a result of gradual radiative heating of non-uniform flux tubes initially in mechanical equilibrium. They find (see Figure 16) significantly stronger fields in the leading leg compared to the following in the lower solar convection zone for all initial field strengths, but nearly equal field strengths of the two legs in the upper convection zone for stronger initial fields (B ≳ 6 × 104 G).

Figure 16:
figure 16

Plots of the magnetic field strength as a function of depth along the emerging loops calculated from the thin flux tube model of Fan and Fisher (1996) showing the asymmetry in field strength between the leading leg (solid curve) and the following leg (dash-dotted curve) of each loop. Panels (a), (b), and (c) correspond to the cases with initial toroidal field strengths of 3 × 104 G, 6 × 104 G, and 105 G respectively. The flux Φ = 1022 Mx and the initial latitude θ = 5° are the same for the three cases shown.

Clearly, further investigations, and perhaps fully resolved three-dimensional MHD simulations of emerging flux loops, are needed to understand the origin of the observed morphological asymmetries of the two polarities of solar active regions.

5.1.4 Geometrical asymmetry of emerging loops and the asymmetric proper motions of active regions

Another asymmetry in the emerging loop generated by the effect of the Coriolis force is the asymmetry in the east-west inclinations of the two sides of the loop. This asymmetry is first shown in the thin flux tube calculations of Moreno-Insertis et al. (1994) and Caligari et al. (1995) who modeled emerging loops that develop self-consistently as a result of the buoyancy instability of toroidal magnetic flux tubes initially in mechanical equilibrium. Moreno-Insertis et al. (1994) and Caligari et al. (1995) found that as the emerging loop rises, the Coriolis force, or the tendency for the tube to conserve angular momentum, drives a counter-rotating motion of the tube plasma, which causes the summit of the loop to move retrograde relative to the valleys, resulting in an asymmetry in the inclinations of the two legs of the loop with the leading leg being inclined more horizontally with respect to the surface than the following leg. This asymmetry in inclination can be clearly seen in Figure 17, which shows a view from the north pole of an asymmetric emerging loop obtained from a simulation by Caligari et al. (1995).

Figure 17:
figure 17

A view from the north pole of the configuration of an emerging loop obtained from a thin flux tube simulation of a buoyantly unstable initial toroidal flux tube. The initial field strength is 1.2 × 105 G, and the initial latitude is 15°. Note the strong asymmetry in the east-west inclination of the two sides of the emerging loop. From Caligari et al. (1995).

Caligari et al. (1998) and Fan and Fisher (1996) show that the geometrical asymmetry is a robust result for models of emerging loops that originate from initial flux tubes in mechanical equilibrium. Loops with initial field strength ranging from 3 × 104 to 105 G consistently show this asymmetry.

The observational consequences of this geometric asymmetry are discussed in Moreno-Insertis et al. (1994) and Caligari et al. (1995). The emergence of such an eastward inclined loop is expected to produce apparent asymmetric east-west proper motions of the two polarities of the emerging region, with a more rapid motion of the leading polarity spots away from the emerging region compared to the motion of the following polarity spots. Such asymmetric proper motions are observed in young active regions and sunspot groups (see Chou and Wang, 1987; van Driel-Gesztelyi and Petrovay, 1990; Petrovay et al., 1990). Furthermore, the asymmetry in the inclination of the emerging loop may also explain the observation that the magnetic inversion line in bipolar regions is statistically nearer to the main following spot than to the main proceeding one (van Driel-Gesztelyi and Petrovay, 1990; Petrovay et al., 1990).

5.1.5 Other properties at the apex of the emerging loop

Besides the basic asymmetric properties of the emerging loops described in the previous sections, the thin flux tube model also provides useful information with regard to the evolution of the field strength, rise speed, etc. of the rising tube under the perfect flux frozen-in condition (without being subject to numerical diffusion). Figure 18 shows the evolution of a set of quantities at the apex of a rising flux tube as it traverse the convection zone, based on a thin flux tube simulation of an emerging Ω-loop developed due to the Parker instability of an initial 105 G toroidal flux ring in mechanical equilibrium at the base of the convection zone (Fan and Gong, 2000, corresponding to the case shown in the top panel of Figure 1 in that paper). It can be seen that the rise velocity v r remains ≲ 200 ms−1 and the Alfvén speed v a remains nearly constant (being much greater than both the rise and the convective flow speed) in the bulk of the convection zone, until the top few tens of Mm of the convection zone, where v r accelerates steeply and v a decreases rapidly due to the steep super-adiabaticity in this top layer. At a depth of roughly 20 Mm, the radius of the tube a exceeds the local pressure scale height H p and v r also exceeds v a. At this point, the thin flux tube approximation breaks down and the tube is likely to be severely distorted and fragmented. Nevertheless, if one continues to use the v r from this point on as an estimate, one finds that the tube will rise through the last 20 Mm depth of the convection zone in only about 7 hours. These results provide some basic information for local helioseismology (see e.g. review by Gizon and Birch, 2005) to estimate possible helioseismic signatures (e.g. wave travel time changes due to thermodynamic perturbations and plasma flows) for detecting subsurface emerging active region flux tubes.

Figure 18:
figure 18

The evolution at the tube apex of (top-left panel) the Alfvén velocity v a, rise velocity v r, convective velocity v conv from a model solar convection zone of Christensen-Dalsgaard (Christensen-Dalsgaard et al., 1993), the azimuthal velocity v φ, (top-right panel) the magnetic field strength B, (bottom-left panel) time elapsed since the onset of the Parker instability, and (bottom-right panel) the tube radius a and the local pressure scale height H p, as a function of depth, resulting from a thin flux tube simulation of an emerging Ω-shaped tube described in Fan and Gong (2000, corresponding to the case shown in the top panel of Figure 1 in that paper). From Fan (2009a).

5.2 Helioseismic probing of subsurface emerging flux

There have been several local helioseismology studies looking for wave speed perturbations and plasma flow signatures that may be produced by subsurface rising flux tubes in emerging active regions (e.g. Kosovichev and Duvall Jr, 2008; Komm et al., 2009). So far there have been no definitive detections of any significant perturbations or signatures associated with the emerging flux prior to the appearance of the flux at the visible surface.

Time-distance helioseismology analysis of a large emerging active region (AR 10488) observed by SoHO MDI show that the onset of wave speed perturbations within the top 20 Mm layer correlates well with the growth of the magnetic flux at the surface, with no significant perturbations prior to the growth of magnetic flux at the surface (Kosovichev and Duvall Jr, 2008). This suggests that the magnetic flux emerges very rapidly through the top 20 Mm layer, consistent with the estimate based on the thin flux tube model given above, such that there is no significant time difference (on the scale ≳ a few hours, which is the temporal resolution for the measured wave speed perturbation) between the onset of the sound-speed perturbation associated with the emerging active region and the appearance of the photospheric magnetic flux. On the other hand, time-distance helioseismology measurement of subsurface flows in the same region of flux emergence show that there is possibly an onset of diverging horizontal flow (and the associated upflow) at the depth range of 1–6 Mm, a few hours prior to the growth of the photospheric magnetic flux. However the signal fluctuates. There also appears to be a localized shear flow forming at the depth of 2 Mm below the photosphere at the location of the first magnetic field signal a few hours before the appearance of the magnetic flux. These interesting signatures of the horizontal flows need to be further studied for more emerging active regions to examine their significance.

By analyzing about five years of GONG high-resolution Doppler data with ring-diagram analysis, Komm et al. (2009) have studied the temporal variation of subsurface horizontal flows of 788 active regions and 978 quiet regions in the depth range of 0–16 Mm, during their disk passage within 60° CMD. A subsurface vertical velocity component is also derived from the divergence of the measured horizontal flows using the requirement of mass conservation, and assuming the flows are subsonic. The regions are sorted based on the variation of their unsigned flux during their disk passage into five subsets of equal size ranging from emerging-flux to decaying-flux subsets. It is found that the average vertical flows of the emerging-flux subset are systematically shifted toward upflows (or diverging horizontal flows) compared to the grand average values of the complete data set, whereas the average flows of the decaying-flux subset show much more pronounced downflows (or converging flows). A study of the averaged cross-correlation between the temporal variation of the unsigned flux and the vertical velocity for the emerging-flux subset suggests that flows in the near-surface and the deeper layers might change about one day before flux emerges at the surface. Thus the change in the divergence of the subsurface horizontal flows might be a precursor of the flux changes. Clearly further observational studies are needed to determine whether significant subsurface horizontal flow signatures associated with emerging flux tubes can be detected prior to the appearance of active region flux at the surface.

5.3 Hemispheric trend of the twist in solar active regions

Vector magnetic field observations of active regions on the photosphere have revealed that on average the solar active regions have a small but statistically significant mean twist that is left-handed in the northern hemisphere and right-handed in the southern hemisphere (see Pevtsov et al., 1995, 2001, 2003). What is being measured is the quantity α ≡ 〈J z/B z〉, the ratio of the vertical electric current over the vertical magnetic field averaged over the active region. When plotted as a function of latitude, the measured α for individual solar active regions show considerable scatter, but there is clearly a statistically significant trend for negative (positive) α in the northern (southern) hemisphere (see Figure 19).

Figure 19:
figure 19

The figure shows the latitudinal profile of α best (see Pevtsov et al., 1995, for the exact way of determining α best) for (a) 203 active regions in cycle 22 (Longcope et al., 1998), and (b) 263 active regions in cycle 23. Error bars (when present) correspond to 1 standard deviation of the mean α best from multiple magnetograms of the same active region. Points without error bars correspond to active regions represented by a single magnetogram. The solid line shows a least-squares best-fit linear function. From Pevtsov et al. (2001).

A linear least squares fit to the data of α as a function of latitude (Figure 19a) found that α = −2.7 × 10−10θ deg m−1, where θ deg is latitude in degrees, and that the r.m.s. scatter of α from the linear fit is Δα = 1.28 ×10−8 m−1 (Longcope et al., 1998). The observed systematic α in solar active regions may reflect a systematic field line twist in the subsurface emerging flux tubes.

If the measured α values are a direct consequence of the emergence of twisted magnetic flux tubes from the interior, then it would imply subsurface emerging tubes with a field line twist of q = α/2 (Longcope and Klapper, 1997), where q denotes the angular rate of field line rotation about the axis over a unit axial distance along the tube. Several subsurface mechanisms for producing twist in emerging flux tubes have been proposed (see e.g. the review by Petrovay et al., 2006). The twist may be due to the current helicity in the dynamo generated toroidal magnetic field, from which buoyant flux tubes form at the base of the convection zone (Gilman and Charbonneau, 1999), or it may be acquired during the rise of the flux tubes through the solar convection zone (Longcope et al., 1998; Choudhuri, 2003; Choudhuri et al., 2004; Chatterjee et al., 2006).

Longcope et al. (1998) explain the origin of the observed twist in emerging active region flux tubes as a result of buffeting by the helical turbulence in the solar convection zone during the rise of the tubes. Applying the dynamic model of a weakly twisted thin flux tube (Section 2.1), Longcope et al. (1998) modeled the rise of a nearly straight, initially untwisted tube, buffeted by a random velocity field representative of the turbulent convection in the solar convection zone, which has a nonzero kinetic helicity due to the effect of solar rotation. The kinetic helicity causes helical distortion of the tube axis, which in turn leads to a net twist of the field lines about the axis in the opposite sense within the tube as a consequence of the conservation of magnetic helicity. This process is termed the Σ-effect by Longcope et al. (1998). Quantitative model calculations of Longcope et al. (1998) show that the Σ-effect can explain the hemispheric sign, magnitude, latitude variation, and the r.m.s. dispersion of the observed α of solar active regions. Furthermore, because the aerodynamic drag force acting on flux tubes by the convective flows has an a −2 dependence, where a is the tube radius, this model predicts that the resulting twist generated should also have an a −2 dependence, suggesting a systematic trend for greater twist in smaller flux tubes.

As discussed in Section 5.1.2, Ω-shaped emerging loops are themselves acted upon by the Coriolis force, developing a “tilt” of the loop. This helical deformation of the tube axis will then also induce a twist of the field lines of the opposite sense within the tube as a consequence of conservation of magnetic helicity. Calculations based on the weakly twisted thin flux tube model show that this twist generated by the large scale tilting of the emerging Ω-loop resulting from the Coriolis force has the right hemispheric sign and latitude dependence, but is of too small a magnitude to account for the observed twist in solar active regions (Longcope and Klapper, 1997; Fan and Gong, 2000).

Another interesting and natural explanation for the origin of twist in emerging flux tubes is the accretion of the background mean poloidal field onto the rising flux tube as it traverse through the solar convection zone (Choudhuri, 2003; Choudhuri et al., 2004; Chatterjee et al., 2006). In a Babcock-Leighton type dynamo, the dispersal of solar active regions with a slight mean tilt angle at the surface generates a mean poloidal magnetic field. The mean tilt angle of solar active regions is produced by the Coriolis force acting on the rising flux tubes (see Section 5.1.2). In the northern hemisphere, when a toroidal flux tube rises into a poloidal field that has been created due to the tilt of the same type of flux tubes emerged earlier, the poloidal field gets wrapped around the flux tube will produce a left-handed twist for the tube. This is illustrated in Figure 20.

Figure 20:
figure 20

This figure illustrates that in the northern hemisphere, when a toroidal flux tube (whose cross-section is the hashed area with a magnetic field going into the paper) rising into a region of poloidal magnetic field (in the clockwise direction) generated by the Babcock-Leighton type α-effect of earlier emerging flux tubes of the same type, the poloidal field gets wrapped around the cross-section of the toroidal tube and reconnects behind it, creating an emerging flux tube with left-handed twist. In this figure, the north-pole is to the left, equator to the right, and the dashed line indicating the solar surface. Note the α-effect for the Babcock-Leighton type solar dynamo model mentioned above is not to be confused with the α value measured in solar active region discussed in this section. From Choudhuri (2003).

Using a circulation-dominated (or flux-transport) Babcock-Leighton type mean-field dynamo model, Choudhuri et al. (2004) did a rough estimate of the twist acquired by an emerging flux tube rising through the solar convection zone. Figure 21 shows the resulting butterfly diagram indicating the sign of α of the emerging regions as a function of latitude and time.

Figure 21:
figure 21

Simulated butterfly diagram of active region emergence based on a circulation-dominated mean-field dynamo model with Babcock-Leighton α-effect. The sign of the twist of the emerging active region flux tube is determined by considering poloidal flux accretion during its rise through the convection zone. Right handed twist (left handed twist) is indicated by plus signs (circles). From Choudhuri et al. (2004).

It is found that at the beginning of a solar cycle, there is a short duration where the sign of α is opposite to the preferred sign for the hemisphere. This is because of the phase relation between the toroidal and poloidal magnetic fields produced by this types of solar dynamo models. At the beginning of a cycle, the mean poloidal magnetic field in the convection zone is still dominated by that generated by the emerging flux tubes of the previous cycle, and toroidal flux tubes of the new cycle emerging into this poloidal field gives rise to a right-handed (left-handed) twist of the tube in the northern (southern) hemisphere. However, for the rest of the cycle starting from the solar maximum, the poloidal magnetic field changes sign and the twist for the emerging tubes becomes consistent with the hemispheric preference. For the whole cycle, it is found that about 67% of the emerging regions have a sign of α consistent with the hemispheric rule. The rough estimate also shows that the magnitude of the α values produced by poloidal flux accretion is consistent with the observed values, and that there is an α −2 dependence on the radius a of the emerging tube, i.e. smaller sunspots should have greater α values. This prediction is also made by the Σ-effect mechanism.

Given the frozen-in condition of the magnetic field, it is expected that the accreted poloidal flux be confined in a sheath at the outer periphery of the rising tube, and that in order to produce a twist within the tube, some form of turbulent diffusion needs to be invoked (Chatterjee et al., 2006). By solving the induction equation in a co-moving Lagrangian frame following the rising flux tube and using several simplifying assumptions, Chatterjee et al. (2006) modeled the evolution of the magnetic field in the rising tube cross-section as a result of poloidal flux accretion and penetration due to a field strength dependent turbulent diffusivity. They found that with plausible choices of assumptions and parameter values an value comparable to the observations is obtained.

When a buoyant magnetic flux tube formed at the base of the solar convection zone from the dynamo generated (pre-dominantly) toroidal magnetic field, it should already obtain an initial twist due to the weak poloidal mean field contained in the magnetic layer. This initial twist will then be further augmented or altered due to poloidal flux accretion and also due to the Σ-effect as the tube rises through the solar convection zone. MHD simulations of the formation and rise of buoyant magnetic flux tubes directly incorporating the mean field profiles from dynamo models (for both the fields at the base and in the bulk of the convection zone) as the initial state is necessary to quantify the initial twist and contribution from poloidal field accretion. Such simulations should be done for dynamo mean fields at different phases of the cycle to access the cycle variation of the twist in the emerging flux tubes.

Observationally, looking for any systematic variations of the the α value (or twist) of solar active regions with the solar cycle phase is helpful for identifying the main mechanisms for the origin of the twist. So far results in this area have been inconclusive (Bao et al., 2000; Pevtsov et al., 2001). On the other hand, observational studies of the correlation between active region twist (as measured by α) and tilt angles have revealed interesting results (Holder et al., 2004; Tian et al., 2005; Nandy, 2006). The Σ-effect predicts that the twist being generated in the tube is uncorrelated to the local tilt of the tube at the apex (Longcope et al., 1998). However, due to the Coriolis force, active region Ω-tubes acquire a mean tilt that has a well defined latitudinal dependence as described by the Joy’s law. The mean twist generated by the Σ-effect also has a latitudinal dependence that is consistent with the observed hemispheric rule. Thus, due to the mutual dependence of their mean values on latitude there should be a correlation between the tilt angle and twist of solar active regions and the correlation is expected to be positive if one assigns negative (positive) sign to a clockwise (counter-clockwise) tilt. However, Holder et al. (2004) found a statistically significant negative correlation between the twist and tilt for the 368 bipolar active regions studied, opposite to that expected from the mutual dependence on latitude of mean twist and mean tilt of active regions. Removing the effect of the mutual dependence of the mean tilt and mean twist on latitude (by either determining the correlation at a fixed latitude, or by subtracting off the fitted mean tilt and mean twist at the corresponding latitude), they found that the negative correlation is enhanced. Furthermore, it is found that the negative correlation is mainly contributed by those active regions (174 out of 368 regions) that deviate significantly from Joy’s law (by > 6σ), while regions that obey Joy’s law to within 6σ show no significant correlation between their twist and tilt. A separate study by Tian et al. (2005) found that a sample of 104 complex δ-configuration active regions, more than half of which have tilts that are opposite to the direction prescribed by Joy’s law, show a significant negative correlation between their twist and tilt (after correcting for their definition of the sign of tilt which is opposite to the definition used in Holder et al. (2004)). The results of Holder et al. (2004) and Tian et al. (2005) both indicate that there is a significant population (about one half in the case of Holder et al. (2004)) of solar active regions whose twist/tilt properties cannot be explained by the Σ-effect together with the effect of the Coriolis force alone. These active regions are consistent with the situation where the buoyant flux tube form at the base of the solar convection zone has acquired an initial twist, such that as it rises upward due to buoyancy into an Ω-tube, it develops a writhe that is of the same sense as the initial twist of the tube. In the extreme case, the twist can be so large that the flux tube becomes kink unstable (see Section 5.6). The resulting tilt at the apex due to the writhe has the negative correlation with the twist as described in the above observations.

5.4 On the minimum twist needed for maintaining cohesion of rising flux tubes in the solar convection zone

As described in Section 5.1, simulations based on the thin flux tube approximation have revealed many interesting results with regard to the global dynamics of active region emerging flux loops in the solar convective envelope, which provide explanations for several basic observed properties of solar active regions. However, one major question ignored by the thin flux tube model is how a flux tube remains a discrete and cohesive object as it moves in the solar convection zone. The manner in which solar active regions emerge on the photosphere suggests that they are coherent flux bundles rising through the solar convection zone and reaching the photosphere in a reasonably cohesive fashion. To address this question, 3D MHD models that fully resolve the rising flux tubes are needed.

As a natural first step, 2D MHD simulations have been carried out to model buoyantly rising, infinitely long horizontal magnetic flux tubes in a stratified layer representing the solar convection zone, focusing on the dynamic evolution of the tube cross-section. The first of such calculations was done in fact much earlier by Schüssler (1979) and later, simulations of higher numerical resolutions have been performed (Moreno-Insertis and Emonet, 1996; Longcope et al., 1996; Fan et al., 1998a; Emonet and Moreno-Insertis, 1998). The basic result from these 2D models of buoyant horizontal flux tubes is that due to the vorticity generation by the buoyancy gradient across the flux tube cross-section, if the tube is untwisted, it quickly splits into a pair of vortex tubes of opposite circulations, which move apart horizontally and cease to rise. If on the other hand, the flux tube is sufficiently twisted such that the magnetic tension of the azimuthal field can effectively suppress the vorticity generation by the buoyancy force, then most of the flux in the initial tube is found to rise in the form of a rigid body whose rise velocity follows the prediction by the thin flux tube approximation. The result described above is illustrated in Figure 22 which shows a comparison of the evolution of the tube cross-section between the case where the buoyant horizontal tube is untwisted (upper panels) and a case where the twist of the tube is just above the minimum value needed for the tube to rise cohesively (lower panels).

Figure 22:
figure 22

Upper panel: Evolution of a buoyant horizontal flux tube with purely longitudinal magnetic field. Lower panel: Buoyant rise of a twisted horizontal flux tube with twist that is just above the minimum value given by Equation (26). The color indicates the longitudinal field strength and the arrows describe the velocity field. From Fan et al. (1998a). (For a corresponding movie showing the evolution of the tube for the untwisted and the twisted cases refer to Figure 23.)

Figure 23:
figure 23

mpg-Movie (195 KB) Still from a movie showing The evolution of a rising flux tube. From Fan et al. (1998a). For a detailed description see Figure 22 (For video see appendix).

This minimum twist needed for tube cohesion can be estimated by considering a balance between the magnetic tension force from the azimuthal field and the magnetic buoyancy force. For a flux tube near thermal equilibrium whose buoyancy |Δρ/ρ| ∼ 1/β, where Δρρρ e denotes the density difference between the inside and the outside of the tube and βp/(B 2/8π) denotes the ratio of the gas pressure over the magnetic pressure, such an estimate (Moreno-Insertis and Emonet, 1996) yields the condition that the pitch angle Ψ of the tube field lines on average needs to reach a value of order

$$\tan \Psi \equiv {{{B_\phi }} \over {{B_z}}} \mathbin{\lower.3ex\hbox{$\buildrel>\over{\smash{\scriptstyle\sim}\vphantom{_x}}$}} {\left({{a \over {{H_p}}}} \right)^{1/2}}.$$
(26)

In Equation (26), B z and B φ denote the axial and azimuthal field of the horizontal tube respectively, a is the characteristic radius of the tube, and H p is the local pressure scale height. The above result on the minimum twist can also be expressed in terms of the rate of field line rotation about the axis per unit length along the tube q. For a uniformly twisted flux tube, B φ = qrB z, where r is the radial distance to the tube axis. Then q needs to reach a value of order (Longcope et al., 1999)

$$q \mathbin{\lower.3ex\hbox{$\buildrel>\over{\smash{\scriptstyle\sim}\vphantom{_x}}$}} {\left({{1 \over {{H_p}a}}} \right)^{1/2}}$$
(27)

for the flux tube to maintain cohesion during its rise. Note that the conditions given by Equations (26) and (27) and also the 2D simulations described in this section all assume buoyant flux tubes with initial buoyancy |Δρ/ρ| ∼ 1/β. For tubes with lower level of buoyancy, the necessary twist is smaller with tan Ψ and q both ∝ |Δρ/ρ|1/2 (see Emonet and Moreno-Insertis, 1998).

Longcope et al. (1999) pointed out that the amount of twist given by Equation (27) is about an order of magnitude too big compared to the twist deduced from vector magnetic field observations of solar active regions on the photosphere. They assumed that the averaged αJ z/B z (the ratio of the vertical electric current over the vertical magnetic field) measured in an active region on the photosphere directly reflects the twist in the subsurface emerging tube, i.e. q = α/2 (Longcope and Klapper, 1997; Longcope et al., 1998). If this is true then it seems that the measured twists in solar active regions directly contradict the condition for the cohesive rise of a horizontal flux tube with buoyancy as large as |Δρ/ρ| ∼ 1/β.

More recently, 3D simulations of Ω-shaped arched flux tubes have been carried out (Abbett et al., 2000, 2001; Fan, 2001a). Fan (2001a) performed 3D simulations of arched flux tubes which form from an initially neutrally buoyant horizontal magnetic layer as a result of its undulatory buoyancy instability (see Section 4.2 and Figure 7). It is found that without any initial twist the flux tubes that form rise through a distance of about one density scale height included in the simulation domain without breaking up. This significantly improved cohesion of the 3D arched flux tubes compared to the previous 2D models of buoyant horizontal tubes is not only due to the additional tension force made available by the 3D nature of the arched flux tubes, but also due largely to the absence of an initial buoyancy and a slower initial rise (Fan, 2001a). With a neutrally buoyant initial state, both the buoyancy force and the magnetic tension force grow self-consistently from zero as the flux tube arches. The vorticity source term produced by the growing magnetic tension as a result of bending and braiding the field lines is found to be able to effectively counteract the vorticity generation by the growing buoyancy force in the apex cross-section, preventing it from breaking up into two vortex rolls. The 2D models (Moreno-Insertis and Emonet, 1996; Longcope et al., 1996; Fan et al., 1998a; Emonet and Moreno-Insertis, 1998) on the other hand considered an initially buoyant flux tube for which there is an impulsive initial generation of vorticity by the buoyancy force. A significant initial twist is thus required to suppress this initial vorticity generation. Therefore the absence of an initial vorticity generation by buoyancy, and the subsequent magnetic tension force resulting from bending and braiding the field lines allow the arched tube with no net twist in Fan (2001a) to rise over a significantly greater distance without disruption.

Abbett et al. (2000) performed 3D simulations where an initial horizontal flux tube is prescribed with a non-uniform buoyancy distribution along the tube such that it rises into an Ω-shaped loop. As discussed above, due to the prescribed buoyancy in the initial horizontal tube, there is an impulsive initial generation of vorticity by the buoyancy force which breaks up the apex of the rising Ω-loop if there is no initial twist. However the separation of the two vortex fragments at the apex is reduced due to the three-dimensional effect (Abbett et al., 2000). By further including the effect of solar rotation using a local f-plane approximation, Abbett et al. (2001) found that the influence of the Coriolis force significantly suppresses the degree of fragmentation at the apex of the Ω-loop (Figure 24).

Figure 24:
figure 24

The rise of a buoyant Ω-loop with an initial field strength B = 105 G in a rotating model solar convection zone at a local latitude of 15° (from Abbett et al. (2001)). The Ω-loop rises cohesively even though it is untwisted. The loop develops an asymmetric shape with the leading side (leading in the direction of rotation) having a shallower angle relative to the horizontal direction compared to the following side.

They also found that the Coriolis force causes the emerging loop to become asymmetric about the apex, with the leading side (leading in the direction of rotation) having a shallower angle with respect to the horizontal direction compared to the following (Figure 24), consistent with the geometric asymmetry found in the thin flux tube calculations (Section 5.1.4).

Another interesting possibility is suggested by the 3D simulations of Dorch and Nordlund (1998), who showed that a random or chaotic twist with an amplitude similar to that given by Equation (26) or (27) in the flux tube can ensure that the tube rises cohesively. Such a random twist may not be detected in the photosphere measurement of active region twists which is determined by taking some forms of average of the quantity α = J z/B z over the active region.

5.5 A further constraint on the twist of subsurface emerging tubes: results from rotating spherical-shell simulations

Fan (2008) has carried out a set of 3D anelastic MHD simulations of the buoyant rise of active region scale flux tubes in a “quiescent” model solar convective envelope in a rotating spherical shell geometry (see Figure 25 and the associated video).

Figure 25:
figure 25

mpg-Movie (431 KB) Still from a movie showing The evolution of a weakly twisted, buoyantly rising Ω-tube, resulting from a simulation described in Fan (2008, see the LNT run in that paper). From Fan (2008). Figure and movie reproduced with permission of the AAS (For video see appendix).

These simulations have considered twisted, buoyant toroidal flux tubes at the base of the solar convection zone with an initial field strength of 105 G, being ∼ 10 times the equipartition field strength, and thus have neglected the effect of convection. The main finding from these simulations is that the twist of the tube induces a tilt at the apex of the rising Ω-tube that is opposite to the direction of the observed mean tilt of solar active regions, if the sign of the twist follows the observed hemispheric preference. It is found that in order for the tilt driven by the Coriolis force to dominate, such that the emerging Ω-tube shows a tilt consistent with Joy’s law of active region mean tilt, the initial twist rate of the flux tube needs to be smaller than about a half of that required for the tube to rise cohesively. Under such conditions, the buoyant flux tube is found to undergo severe flux loss during its rise, with less than 50% of the initial flux remaining in the final Ω-tube that rises to the surface (see Figures 26a and 26b).

Figure 26:
figure 26

(a) 3D volume rendering of the magnetic field strength of a weakly twisted, rising Ω-tube, whose apex is approaching the top boundary, resulting from a simulation described in Fan (2008, see the LNT run in that paper). (For a corresponding movie see Figure 25.) (b) A cross section of B near the top boundary at r = 0.937R ; (c) selected field lines threading through the coherent apex cross-section of the Ω-tube.

Furthermore, it is found that the Coriolis force drives a retrograde flow along the apex portion, resulting in a relatively greater stretching of the field lines and hence stronger field strength in the leading leg of the tube. With a greater field strength, the leading leg is more buoyant with a greater rise velocity, and remains more cohesive compared to the following leg (see Figure 26a). Figure 26c shows selected field lines threading through the coherent apex cross-section of the final Ω-tube, resulting from the simulation of a weakly twisted buoyant tube described in Fan (2008, see the LNT run in that paper). It can be seen that the field lines in the leading side are winding about each other smoothly in a coherent fashion, while the field lines in the following side are significantly more frayed. By evaluating the local twist rate given by αJ · B/B 2, where J is the electric current density, along each of the selected field lines as a function of depth, Fan (2009a) found that field lines in the leading leg show more coherent values of α, whereas the field lines in the following leg show significantly larger fluctuations and mixed signs of local twist (see Figure 27).

Figure 27:
figure 27

Dots show values of αJ · B/B 2 computed along each of the selected field lines of the final Ω-tube shown in Figure 26(c) as a function of depth for the following side (left panel) and the leading side (right panel). The field-line averaged mean α is shown as the solid curve. From Fan (2009a).

Although the mean value of α averaged over the field lines is not systematically greater along the leading leg compared to the following (Fan, 2009a), the greater buoyancy and hence higher rise velocity of the leading leg can give rise to a greater upward helicity flux in the leading polarity comparing to the following as a result of the emergence of the Ω-tube (Fan et al., 2009). Furthermore, based on a simplified model of active region flux emergence into the corona by Longcope and Welsch (2000), Fan et al. (2009) show that a stronger field strength in the leading tube also produces a faster rotation of the leading polarity sunspot driven by torsional Alfvén waves along the flux tube. This also contributes to a greater helicity injection rate in the leading polarity of an emerging active region. Observational study of bipolar emerging active regions by Tian and Alexander (2009) have found that the helicity injection rate is about 3–10 times greater in the (compact) leading polarity than the (fragmented) following polarity.

Jouve and Brun (2007) have also carried out anelastic MHD simulations in a rotating spherical shell geometry to study the buoyant rise of an axisymmetric toroidal flux ring in an isentropically stratified (non-convecting) envelope. They have considered a even greater initial field strength of 1.8 × 105 G for the initial toroidal flux ring. As was discussed in Fan (2008), the poleward deflection of the rise trajectory of the tube due to the Coriolis force is far more severe for an axisymmetric toroidal ring (where the whole ring is moving away from the rotation axis of the Sun) than for a localized 3D Ω-shaped tube (see Section 3.1 in Fan, 2008). Thus an initial field strength of ≳ 1.8 × 105 G is needed for an axisymmetric toroidal ring to rise nearly radially (Jouve and Brun, 2007). The simulations of Jouve and Brun (2007) also recovered the previous results from Cartesian simulations that if the flux tube is not twisted, it splits into two counter rotating vortices before reaching the top of the envelope.

5.6 The rise of kink unstable magnetic flux tubes and the origin of delta-sunspots

Most of the solar active regions are observed to have very small twists, with an averaged value of α ≡ 〈J z/B z〈 (the averaged ratio of the vertical electric current over the vertical magnetic field) measured to be on the order of 0.01 Mm−1 (see Pevtsov et al., 1995, 2001). However there is a small but important subset of active regions, called the δ-sunspots, which are observed to be highly twisted with α reaching a few times 0.1 Mm−1 (Leka et al., 1996), and to have unusual polarity orientations that are sometimes reversed from Hale’s polarity rule (see Zirin and Tanaka, 1973; Zirin, 1988; Tanaka, 1991). These δ-sunspots are compact structures where umbrae of opposite polarity are contained within a common penumbra. They are found to be the most flare productive active regions (see Zirin, 1988). Through careful analysis of the evolution of flare-active δ-sunspot groups, Tanaka (1991) proposed a model of an emerging twisted flux rope with kinked or knotted geometry to explain the observed evolution of these regions.

Motivated by the observations of δ-sunspots, MHD calculations of the evolution of highly twisted, kink unstable magnetic flux tubes in the solar convection zone have been carried out (Linton et al., 1996, 1998, 1999; Fan et al., 1998b, 1999). For an infinitely long twisted cylindrical flux tube with axial field B z(r), azimuthal field B θ(r) = q (r) rB z(r), and plasma pressure p (r) in hydrostatic equilibrium where \(dp/dr=-(B_{\theta}^{2}/4\pi r)-d(B_{z}^{2}+B_{\theta}^{2})/dr\), a sufficient condition for the flux tube to be kink unstable is (see Freidberg, 1987)

$${r \over 4}{\left({{{q^{\prime}} \over q}} \right)^2} + {{8\pi p^{\prime}} \over {B_z^2}} < 0$$
(28)

to be true somewhere in the flux tube. In Equation (28) the superscript’ denotes the derivative with respect to r. This is known as Suydam’s criterion. Note that condition (28) is sufficient but not necessary for the onset of the kink instability and hence there can be cases which are kink unstable but do not satisfy condition (28). One such example are the force-free twisted flux tubes which are shown to be always kink unstable without line-tying (i.e. infinitely long) (Anzer, 1968), but for which p′ = 0. Force-free fields are the preferred state for coronal magnetic fields under low plasma-β conditions and are not a likely state for magnetic fields in the high-β plasma of the solar interior. Linton et al. (1996) considered the linear kink instability of uniformly twisted cylindrical flux tubes with q = B θ/rB z being constant, confined in a high β plasma. They found that the equilibrium is kink unstable if q exceeds a critical value q cr = a −1, where a −2 is the coefficient for the r 2 term in the Taylor series expansion of the equilibrium axial magnetic field B z about the tube axis: B z(r) = B 0(1 − a −2r 2 + ⋯). This result is consistent with Suydam’s criterion. They further argued that an emerging, twisted magnetic flux loop will tend to have a nearly uniform q along its length since the rise speed through most of the solar convection zone is sub-Alfvénic and torsional forces propagating at the Alfvén speed will equilibrate quasi-statically. Meanwhile expansion of the tube radius at the apex as it rises will result in a decrease in the critical twist q cr = a −1 necessary for the instability. This implies that as a twisted flux tube rises through the solar interior, a tube that is initially stable to kinking may become unstable as it rises, and that the apex of the flux loop will become kink unstable first because of the expanded tube cross-section there (Parker, 1979; Linton et al., 1996).

The non-linear evolution of the kink instability of twisted magnetic flux tubes in a high-β plasma has been investigated by 3D MHD simulations (Linton et al., 1998, 1999; Fan et al., 1998b, 1999). Fan et al. (1998b, 1999) modeled the rise of a kink unstable flux tube through an adiabatically stratified model solar convection zone.

In the case where the initial twist of the tube is significantly supercritical such that the e-folding growth time of the most unstable kink mode is smaller than the rise time scale, they found sharp bending of the flux tube as a result of the non-linear evolution of the kink instability. During the onset of the kink instability, the magnetic energy decreases while the magnetic helicity is approximately conserved. The writhing of the flux tube also significantly increases the axial field strength and hence enhances the buoyancy of the flux tube. The flux tube rises and arches upward at the portion where the kink concentrates, with a rotation of the tube orientation at the apex that exceeds 90° (see Figure 28). The emergence of this kinked flux tube can give rise to a compact magnetic bipole with polarity order inverted from the Hale polarity rule (Figure 29) as often seen in δ-sunspots. The conservation of magnetic helicity requires that the writhing of the tube due to the kink instability is of the same sense as the twist of the field lines. Hence for a kinked emerging tube, the rotation or tilt of the magnetic bipole from the east-west polarity orientation defined by the Hale’s polarity rule should be related to the twist of the tube. The rotation or tilt should be clockwise (counterclockwise) for right-hand-twisted (left-hand-twisted) flux tubes. This tilt-twist relation can be used as a means to test the model of kinked flux tubes as the origin of δ-sunspots (Tanaka, 1991; Leka et al., 1994, 1996; López Fuentes et al., 2003). Observations have found with both consistent and opposing cases (Leka et al., 1996; López Fuentes et al., 2003). A recent study (Tian et al., 2005) which includes a large sample (104) of complex δ-configuration active regions shows that 65−67% of these δ-regions have the same sign of twist and writhe, supporting the model of kinked flux tubes.

Figure 28:
figure 28

mpg-Movie (314 KB) Still from a movie showing The rise of a kink unstable magnetic flux tube through an adiabatically stratified model solar convection zone (result from a simulation in Fan et al. (1999) with an initial right-handed twist that is 4 times the critical level for the onset of the kink instability). In this case, the initial twist of the tube is significantly supercritical so that the e-folding growth time of the most unstable kink mode is smaller than the rise time scale. The flux tube is perturbed with multiple unstable modes. The flux tube becomes kinked and arches upward at the center where the kink concentrates, with a rotation of the tube orientation at the apex that exceeds 90° (For video see appendix).

Figure 29:
figure 29

A horizontal cross-section near the top of the upward arching kinked loop shown in the last panel of Figure 28. The contours denote the vertical magnetic field B z with solid (dotted) contours representing positive (negative) B z. The arrows show the horizontal magnetic field. One finds a compact bipolar region with sheared transverse field at the polarity inversion line. The apparent polarity orientation (i.e. the direction of the line drawn from the peak of the positive pole to the peak of the negative pole) is rotated clockwise by about 145° from the +x direction (the east-west direction) of the initial horizontal flux tube.

5.7 Buoyant flux tubes in a 3D stratified convective velocity field

5.7.1 General considerations

To understand how active region flux tubes emerge through the solar convection zone, it is certainly important to understand how 3D convective flows in the solar convection zone affect the rise and the structure of buoyant flux tubes. The well-defined order of solar active regions as described by the Hale polarity rule suggests that the emerging flux tubes are not subject to strong deformation by the turbulent convection. One can thus consider the following two simplified order-of-magnitude estimates. First, the magnetic buoyancy of the flux tube should probably dominate the downward hydrodynamic force from the convective downflows:

$${{{B^2}} \over {8\pi {H_p}}} > {C_{\rm{D}}}{{\rho v_{\rm{c}}^2} \over {\pi a}} \Rightarrow B > \left({{{2{C_{\rm{D}}}} \over \pi }} \right){\left({{{{H_p}} \over a}} \right)^{1/2}}{B_{{\rm{eq}}}},$$
(29)

where B eq ≡ (4πρ)1/2v c is the field strength at which the magnetic energy density is in equipartition with the kinetic energy density of the convective downdrafts, v c is the flow speed of the downdrafts, H p is the local pressure scale height, a is the tube radius, and C D is the aerodynamic drag coefficient which is of order unity. In Equation (29) we have used the aerodynamic drag force as an estimate for the magnitude of the hydrodynamic forces. The estimate (29) leads to the condition that the field strength of the flux tube needs to be significantly higher than the equipartition field strength by a factor of \(\sqrt{H_{p}/a}\). For flux tubes responsible for active region formation, \(\sqrt{H_{p}/a} > 3\) near the bottom of the solar convection zone. Thus we call B ≳ 3B eq the “magnetic buoyancy dominated regime”, and expect B < 3B eq to be the regime where convective downdrafts become dominant. A second consideration is that the magnetic tension force resulting from bending the flux tube by the convective flows should also dominate the hydrodynamic force due to convection:

$${{{B^2}} \over {8\pi {l_{\rm{c}}}}} > {C_{\rm{D}}}{{\rho v_{\rm{c}}^2} \over {\pi a}} \Rightarrow B > \left({{{{C_{\rm{D}}}} \over \pi }} \right){\left({{{{l_{\rm{c}}}} \over a}} \right)^{1/2}}{B_{{\rm{eq}}}},$$
(30)

where l c denotes the size scales of the convective flows. This leads to a condition for the magnetic field strength very similar to Equation (29) if we consider the largest convective cell scale for l c to be comparable to the pressure scale height. Similar results have also been found in Cline (2003).

5.7.2 Simulations in a local Cartesian geometry without rotation

Fan et al. (2003) carried out direct 3D MHD simulations of the evolution of a buoyant magnetic flux tube in a stratified convective velocity field. The basic result is illustrated in Figure 30.

Figure 30:
figure 30

The evolution of a uniformly buoyant magnetic flux tube in a stratified convective velocity field from the simulations of Fan et al. (2003). Top-left image: A snapshot of the vertical velocity of the 3D convective velocity field in a superadiabatically stratified fluid. The density ratio between the bottom and the top of the domain is 20. Top-right image: The velocity field (arrows) and the tube axial field strength (color image) in the vertical plane that contains the axis of the uniformly buoyant horizontal flux tube inserted into the convecting box. Lower panel: The evolution of the buoyant flux tube with B = B eq(left column) and with B = 10B eq (right column). The color indicates the absolute field strength of the flux tube scaled to the initial tube field strength at the axis. (For a corresponding movie see Figure 31.)

Figure 31:
figure 31

mpg-Movie (328 KB) Still from a movie showing The evolution of a uniformly buoyant magnetic flux tube. From Fan et al. (2003). For a detailed description see Figure 30 (For video see appendix).

They first computed a 3D convective velocity field in a superadiabatically stratified fluid, until the convection reaches a statistical steady state. The resulting velocity field (see top-left image in Figure 30) shows the typical features of overturning convection in a stratified fluid as found in many previous investigations. The surface layer displays a cellular pattern with patches of upflow region surrounded by narrow downflow lanes. In the bulk of the convecting domain, the downflows are concentrated into narrow filamentary plumes, some of which extend all the way across the domain, while the upflows are significantly broader and are of smaller velocity amplitude in comparison to the downdrafts. A uniformly buoyant, twisted horizontal magnetic flux tube having an entropy that is equal to the entropy at the base of the domain is inserted into the convecting box (see top-right image in Figure 30). In the case where the field strength of the tube is in equipartition to the kinetic energy density of the strongest downdraft (left column in the bottom panel of Figure 30), the magnetic buoyancy for this flux tube is weaker than the hydrodynamic force resulting from the convective downflows and the evolution of the tube depends sensitively on the local condition of the convective flows. Despite being buoyant, the portions of the tube in the paths of downdrafts are pushed downward and pinned down to the bottom, while the rise speed of sections within upflow regions is significantly boosted. The pinned-down flux is then further distorted and transported laterally by the horizontal diverging flow at the bottom. On the other hand in the case where the tube field strength is 10 times the equipartition value (right column in the bottom panel of Figure 30), the horizontal flux tube rises under its uniform buoyancy, nearly unaffected by the convection. In this case the horizontal flux tube is sufficiently twisted so that it does not break up into two vortex tubes.

In case B = B eq, it is found that the random north-south tilting of the flux tube caused by convection is of the amplitude ∼ 30°, which is greater than the r.m.s. scatter of the active region tilts away from Joy’s law for large active regions (10°), but is not beyond the r.m.s. tilt scatter for small active regions (30°) (see Fisher et al., 1995). This indicates that the distortion of flux tubes of equipartition field strength by the convective flows during their buoyant rise through the solar convection zone is probably too large to be consistent with the observational constraint of tilt dispersion for large solar active regions. Furthermore it should be noted that the realistic convective flows in the solar convection zone is probably far more turbulent than that computed in Fan et al. (2003), containing flows of scales significantly smaller than the cross-section of the flux tube. Hence the flux tube distortion found in the simulations of Fan et al. (2003) is most likely a lower limit.

The distances between the major downflow plumes are also an important factor in determining the fate of the buoyant flux tubes of equipartition field strength. In Fan et al. (2003), the distance between neighboring downflow plumes can be as large as about 5H p, hence allowing the portion of the buoyant tube between the plumes to rise up to the top of the domain. If on the other hand, the distances between the downflow plumes are below ∼ 2H p, then the tension force for the tube between the pinned-down points will exceed the magnetic buoyancy force (∼ B 2/8πH p) and the entire flux tube will be prevented from emerging to the surface.

5.7.3 Global rotating spherical-shell simulations

Jouve and Brun (2009) have carried out the first set of global anelastic MHD simulations of the buoyant rise of an initially toroidal flux ring in a rotating, fully convective spherical shell, possessing self-consistently generated mean flows such as meridional circulations and differential rotation, representative of the conditions of the solar convective envelope (see e.g. review by Miesch, 2005). They inserted into the fully developed convecting envelope a buoyant toroidal flux ring with different initial field strengths, twist rates, and initial latitudes, to study how the flux tube rises in the presence of convection and the associated mean flows, and how the dynamic evolution depends on the above initial parameters.

It is found that the magnetic field strength corresponding to the value that is in equipartition with the kinetic energy of the strongest downflows is rather high, B eq ≈ 6.1 × 104 G. The initial field strength B of the toroidal flux ring considered in the simulations are all significantly greater than B eq, being 2.5B eq, 5B eq, and 10B eq. Thus, except for the case with B = 2.5B eq, all of the other cases simulated are essentially in the magnetic buoyancy dominated regime (with B > 3B eq as discussed in Section 5.7.1). As a result, the simulations recovered many of the findings obtained from previous simulations in the absence of convective flows. These include the dependence of the poleward deflection of the tube on the initial tube field strength (e.g. Choudhuri and Gilman, 1987; Fan, 2008), the critical dependence on the initial twist for the cohesion of the buoyantly rising flux tube (e.g. Emonet and Moreno-Insertis, 1998; Abbett et al., 2000), and the dependence of the tilt angle of the emerging tube on the initial twist (Section 5.5 and Fan, 2008). Due to the relatively high magnetic diffusivity in the code, flux tubes with a very large initial field strength (ranging from 1.5 × 105 G to 6 ×105 G) and a large radius, corresponding to a total flux on the order of a few times 1023 Mx, significantly greater than the typical active region fluxes, are considered, such that the rise times of the flux tubes are ≲ the diffusive time scale of about 14.5 days. Because most of the cases considered are essentially in the magnetic buoyancy dominated regime, the rising toroidal flux tube only develops rather moderate undulations by the influence of the convective flows (see Figure 32), and Ω tubes with undulations extending the depth of the convection zone are not found.

Figure 32:
figure 32

Cut at the latitude of 30° of the radial velocity (color) and of the magnetic energy (line contours) for three different simulations of the rise of a buoyant toroidal flux ring with different initial field strengths: 2.5B eq (top panel), 5B eq (middle panel), and 10B eq (bottom panel). From Jouve and Brun (2009). Figure reproduced with permission of the AAS.

It is also found that flux tubes introduced at lower latitudes (e.g. at 15°) have difficulty reaching the top of the domain (even with a strong initial field strength of 5B eq ≈ 3 × 105 G), and the authors attributed the cause of this to the differential rotation. For the weakest field strength case (with B = 2.5B eq = 1.5 × 105 G), it is found that portions of the toroidal ring are pinned down by the convective downdrafts, and eventually the tube loses its buoyancy due to magnetic diffusion and is unable to rise to the top (see top panel of Figure 32).

Clearly simulations with a reduced magnetic diffusion are necessary to model the evolution of rising flux tubes in more realistic parameter regimes. Specifically, it is important to model cases with a weaker initial field strength (104 G ≲ B ≲ 105 G) and thus smaller magnetic buoyancy to study whether large-undulation Ω-shaped emerging tubes with properties consistent with solar active regions can develop (see more discussions in Section 9.1). Note also that the twist of the initial toroidal flux rings introduced in the northern hemisphere (see Figure 1 of Jouve and Brun, 2009) in the simulations is right-handed or positive, which is opposite to the observed preferred sign of twist (left-handed) for active regions in the northern hemisphere (e.g. Pevtsov et al., 2001). With a magnitude of the twist that is just above the critical twist needed for the buoyant flux tube to rise cohesively, the orientation of the tilt angle of the final emerging tube is found to be largely determined by the sign of the twist rather than by the Coriolis force (Fan, 2008). Thus the sign of the tilt angle of the final emerging tube in the simulations may reverse (i.e. becomes opposite to the direction of active region tilts) if toroidal flux rings with a negative or left-handed initial twist were introduced in the northern hemisphere.

6 Turbulent Pumping of a Magnetic Field in the Solar Convection Zone

As discussed in Section 5.7, for buoyant flux tubes with significantly super-equipartition field strength B ≳ (H p/a)1/2B eq, where B eq is in equipartition with the kinetic energy density of the convective downflows, the magnetic buoyancy of the tubes dominates the hydrodynamic force from the convective downflows and the flux tubes can rise unaffected by convection. On the other hand if the field strength of the flux tubes is comparable to or smaller than the equipartition value B eq, the magnetic buoyancy is weaker than the hydrodynamic force from the convective downflows and the evolution of the tubes becomes largely controlled by the convective flows. In this regime of convection dominated evolution, due to the strong asymmetry between up- and downflows characteristic of stratified convection, it is found that magnetic flux is preferentially transported downward against its magnetic buoyancy out of the turbulent convection zone into the stably stratified overshoot region below. This process of “turbulent pumping” of a magnetic field has been demonstrated by several high resolution 3D compressible MHD simulations (see Tobias et al., 1998, 2001; Dorch and Nordlund, 2001).

Tobias et al. (2001) carried out a series of 3D MHD simulations to investigate the turbulent pumping of a magnetic field by stratified convection penetrating into a stably stratified layer. A thin slab of a unidirectional horizontal magnetic field is introduced into the middle of an unstably stratified convecting layer, which has a stable overshoot layer attached below. It is found that the fast, isolated downflow plumes efficiently pump flux from the convecting layer into the stable layer on a convective time scale. Tobias et al. (2001) quantify this flux transport by tracking the amount of flux in the unstable layer and that in the stable layer, normalized by the total flux. Within a convective turnover time, the flux in the stable layer is found to increase from the initial value of 0 to a steady value that is greater than 50% of the total flux (reaching 80% in many cases), i.e. more than half of the total flux is settled into the stable overshoot region in the final steady state. Moreover the stable overshoot layer is shown to be an effective site for the storage of a toroidal magnetic field. If the initial horizontal magnetic slab is put in the stable overshoot layer, the penetrative convection is found to be effective in pinning down the majority of the flux against its magnetic buoyancy, preventing it from escaping into the convecting layer. It should be noted however that the fully compressible simulations by Tobias et al. (2001) assumed a polytropic index of m = 1 for the unstable layer which corresponds to a value for the non-dimensional superadiabaticity δ of 0.1. Thus the convective flow speed v c for the strong downdrafts is on the order of \(\sqrt{\delta}\sim 0.3\) times the sound speed, i.e. not very subsonic. On the other hand, the range of initial magnetic field strength considered in the simulations corresponds to a plasma β ranging from 2 × 103 to 2 × 107. Thus comparing to the kinetic energy density of the convective flows, the strongest initial field considered is of the order 0.1 B eq, i.e. below equipartition, although a field with strength up to B eq may be generated as a result of amplification by the strong downflows during the evolution. Therefore the above results obtained with regard to the efficient turbulent pumping of a magnetic field out of the convection zone into the stable overshoot region against its magnetic buoyancy apply only to fields weaker than or at most comparable to B eq.

Abbett et al. (2004) found that turbulent pumping is very weak and ineffective in an MHD convection model without a stable overshoot layer at the bottom. Considering both the results of Tobias et al. (2001) and Abbett et al. (2004) it appears that the presence of the stable overshoot layer below the convection zone is an essential ingredient for effective turbulent pumping.

The turbulent pumping of magnetic flux with field strength BB eq out of the convection zone into the stable overshoot region demonstrated by the high resolution 3D MHD simulations has profound implications for the working of the interface mean field dynamo models (Parker, 1993; Charbonneau and MacGregor, 1997) as discussed by Tobias et al. (2001). The interface dynamo models require efficient transport of the large scale poloidal field generated by the α-effect of the cyclonic convection out of the convection zone into the stably stratified tachocline region where strong rotational shear generates and amplifies the large scale toroidal magnetic field. Turbulent pumping is shown to enhance both the transport of magnetic flux into the stable shear layer and the storage of the toroidal magnetic field there. It further implies that the transport of magnetic flux by turbulent pumping should not be simply treated as an enhanced isotropic turbulent diffusivity in the convection zone as is typically assumed in the mean field dynamo models. Tobias et al. (2001) noted that a better treatment would be to add an extra advective term to the mean field equation characterizing the effect of turbulent pumping, which would correspond to including the antisymmetric part of the α-tensor known as the γ-effect (see Moffatt, 1978).

7 Amplification of a Toroidal Magnetic Field by Conversion of Potential Energy

Thin flux tube models of emerging flux loops through the solar convective envelope (Section 5.1) have inferred a strong super-equipartition field strength of order 105 G for the toroidal magnetic field at the base of the solar convection zone. Generation of such a strong field is dynamically difficult since the magnetic energy density of a 105 G field is about 10–100 times the kinetic energy density of the differential rotation (Parker, 1994; Rempel and Schüssler, 2001). An alternative mechanism for amplifying the toroidal magnetic field has been proposed which converts the potential energy associated with the stratification of the convection zone into magnetic energy. It is shown that upflow of high entropy plasma towards the inflated and “exploding” top of a rising Ω-loop developed from an initial toroidal field of equipartition field strength (∼ 104 G) can significantly intensify the submerged part of the field by extracting plasma out of it (Parker, 1994; Moreno-Insertis et al., 1995; Rempel and Schüssler, 2001). This process is a barometric effect and is caused by the entropy gradient in the solar convection zone maintained by the energy transport.

In the thin flux tube simulations of rising Ω-loops in the solar convection zone, it is found that flux loops with a low initial field strength of ∼ 104 G do not reach the upper half of the solar convection zone before the apexes of the loops loose pressure confinement and effectively “explode” (Moreno-Insertis et al., 1995). This loss of pressure confinement at the top of the emerging loop occurs as plasma inside the tube establishes hydrostatic equilibrium along the tube which happens if the emerging loop rises sufficiently slowly (Moreno-Insertis et al., 1995). The loop rises adiabatically carrying the high entropy plasma from the base of the solar convection zone while the entropy outside the flux tube decreases with height in the superadiabatically stratified convection zone. Hydrostatic equilibrium therefore dictates that the plasma pressure inside the flux tube decreases with height slower than the outside and becomes equal to the external pressure at a certain height where the magnetic field can no longer be confined. This “explosion” height for the emerging loop is found to be a function of the initial tube field strength at the base of the solar convection zone (see Figure 33).

Figure 33:
figure 33

Depth from the surface where the apex of an emerging flux loop with varying initial field strength rising from the bottom of the convection zone looses pressure confinement or “explodes” as a result of tube plasma establishing hydrostatic equilibrium along the tube. The explosion height is computed by considering an isentropic thin flux loop with hydrostatic equilibrium along the field lines (see Moreno-Insertis et al., 1995) in a model solar convection zone of Christensen-Dalsgaard (Christensen-Dalsgaard et al., 1993).

For loops with an initial field strength ∼ 104 G, the explosion height is at about the middle of the solar convection zone. When the loop apex approaches the explosion height, it expands drastically and the buoyancy of the high entropy material in the tube is expected to drive an outflow which extracts plasma out of the lower part of the flux tube at the base of the solar convection zone. This process has been demonstrated by Rempel and Schüssler (2001), who performed MHD simulations of exploding magnetic flux sheets in two-dimensional Cartesian geometry.

The simulations of Rempel and Schüssler (2001) start with a magnetic sheet with a higher value of entropy placed at the bottom of an adiabatically stratified layer (constant entropy layer). This setup avoids the complication of involving convective flows in the simulations while keeping the essential effect of the entropy decrease in the solar convection zone by assuming a constant entropy difference between the flux sheet and the isentropic layer. The central portion of the flux sheet is perturbed upward which subsequently forms a rising loop as a result of its magnetic buoyancy (see Figure 34).

Figure 34:
figure 34

Evolution of magnetic field strength (gray scale: darker gray denotes stronger field) and velocity field (arrows) during the flux loop explosion. The horizontal part of the field is amplified by a factor of 3. From Rempel and Schüssler (2001).

The apex of the rising loop explodes into a cloud of weak magnetic field as it crosses the predicted explosion height (middle panel) and high entropy plasma, driven by buoyancy, continues to flow out of the “stumps”, draining mass from the lower horizontal part of the flux sheet (bottom panel). The field strength of the horizontal part of the flux sheet is visibly intensified. It is found that the final field strength the horizontal part of the field can reach is roughly the value for which the explosion height is close to the top of the stratification. For the solar convection zone, this value corresponds to ∼ 105 G (see Figure 33). This implies that the large field strength of order 105 G may be achieved by the process of flux “explosion”, which draws energy from the potential energy associated with the stratification of the solar convection zone. A coherent picture maybe as follows. Differential rotation in the tachocline shear layer generates and amplifies the toroidal magnetic field to an equipartition value of about 104 G, at which point magnetic buoyancy becomes important dynamically, and magnetic buoyancy instability, radiative heating or other convective perturbations drive the formation of buoyant flux loops rising into the solar convection zone. These rising loops explode in the middle of the convection zone and fail to rise all the way to the surface. These “failed” eruptions “pump” out the plasma from the toroidal field at the base of the solar convection zone (Parker, 1994) and amplify the toroidal field until it reaches strongly super-equipartition field strength of order 105 G, whose eruptions then lead to the emergence of solar active regions at the surface.

8 Flux Emergence at the Surface and Post-Emergence Evolution of Subsurface Fields

8.1 Evolution in the top layer of the solar convection zone and the photosphere

A complete 3D MHD model of emerging active region flux tubes, extending from the base of the solar convection zone up into the visible solar atmosphere, is not yet possible. Thin flux tube calculations of rising flux tubes in the solar convection zone usually extend from the base of the solar convection zone up to roughly 20–30 Mm below the surface, at which depths the validity of the thin flux tube approximation breaks down since the diameter of the flux tube begins to exceed the local pressure scale height (see Section 5.1.5 and Moreno-Insertis, 1992). The domain of fully 3D MHD simulations of rising magnetic flux tubes in a model convection zone typically contains a vertical stratification of up to 3 density scale heights (see Abbett et al., 2000, 2001; Fan et al., 2003), which corresponds to the density stratification over the range from the bottom of the solar convection zone to roughly 36 Mm below the photosphere. Recent 3D spherical shell simulations of buoyantly rising flux tubes in the solar convective envelope (Fan, 2008) cover depths from the base of the solar convection zone to about 16 Mm below the surface. Because of the rapid decrease of the various scale heights in the top layer of the solar convection zone which demands increasing numerical resolution, it is not yet feasible to perform 3D MHD simulations that extend from the bottom of the convection zone all the way to the photosphere. Furthermore, there is an increased complexity in the physics of the top layer of the solar convection zone. The thermodynamics of the plasma is complicated by ionization effects and the radiative exchange is expected to play an important role in the heat transport (see review by Nordlund et al., 2009). The anelastic approximation breaks down because the plasma flow speed is no longer slow compared to the sound speed.

Rapid progress has been made in recent years in fully compressible 3D MHD simulations of magneto-convection, emerging magnetic flux, and sunspot structure within the top few Mm layer of the solar convection zone and the overlying photospheric layer, incorporating realistic physics of partial ionization of the dominant constituents and radiative transfer (see e.g. Stein and Nordlund, 2000; Cheung et al., 2007, 2008; Martínez-Sykora et al., 2008, 2009; Rempel et al., 2009; Nordlund et al., 2009). These simulations have produced results that can be directly compared with high resolution photospheric observations of the solar granulations, emerging flux regions, sunspot fine structure and the associated flows. Cheung et al. (2008) carried out 3D radiation MHD simulations of a twisted magnetic flux tube rising through the top layer of the solar convection zone (from a depth of about 5.5 Mm below the photosphere) into the photosphere. It is found that due to the strong stratification of the top layer of the convection zone, the rise of the flux tube is accompanied by a strong lateral expansion. By the time it has reached the photosphere, it appears more like a flux sheet which acts as a reservoir for small-scale flux emergence events at the granulation scale. Detailed comparisons of the simulation results of flux emergence at the photosphere layer and the new observational data from SOT of Hinode provide physical interpretations for many of the observed features in emerging flux regions (EFRs). For example, convective downflows produce serpentine-shaped emerging field lines which result in the observed mixed-polarity pattern in the interior of EFRs, where opposite-polarity flux concentrations appear to counter-stream (see Figure 35 and the associated movie).

Figure 35:
figure 35

mpg-Movie (7920 KB) Still from a movie showing Continuum intensity images and synthetic magnetograms of the simulated emerging flux region resulting from a simulation of flux emergence from the top layer of the solar convection zone into the solar photosphere. From (Cheung et al., 2008). Figure and movie reproduced with permission of the AAS (For video see appendix).

At each of the two opposite edges of the EFR, flux of one sign tends to coalesce. This may eventually lead to the formation of solar pores and sunspots, although the simulations so far have not been able to run long enough to see this actually happening.

Another interesting observed feature reproduced by the simulations is the presence of supersonic downflows at some flux “cancellation” sites. This is revealed by the simulations to correspond to the retraction of inverted U-loops due to the magnetic tension force. Simulations also found examples of surface flux concentrations undergoing convective intensification leading to the formation of Kilogauss fields and associated bright points. Such events have recently been directly observed by the high resolution Solar Optical Telescope (SOT) of the Hinode satellite (e.g. Nagata et al., 2008; Fischer et al., 2009). The intensification process is consistent with the basic theory of “convective collapse”, which results from the convective instability of plasma inside the vertical thin flux tubes in the top few hundred kilometers of the solar convection zone (see Parker, 1978; Spruit, 1979; Spruit and Zweibel, 1979), but is operating under the more realistic conditions with the effect of radiative energy transfer included (e.g. Cheung et al., 2008)

Realistic magneto-convection simulations of the evolution of emerging active region scale flux tubes in the top ∼ 20 Mm layer of the solar convection zone are yet to be carried out. It remains an open question how the top of the emerging Ω-shaped tube intensifies to form sunspots with Kilogauss field strength and β ∼ 1 at the photosphere.

8.2 Flux emergence into the solar atmosphere and the corona

Understanding how twisted magnetic fields emerge from the dense, convectively unstable solar convection zone into the stably stratified, rarefied solar atmosphere and corona is fundamentally important for understanding the formation of solar active regions and the development of precursor structures for solar eruptions such as flares and coronal mass ejections. Pioneering work in this area was carried out by Shibata and collaborators (see Shibata et al., 1989) who showed that the magnetic buoyancy instability is a mechanism through which magnetic flux reaching the photosphere can expand dynamically into the stably stratified solar atmosphere.

In recent years, a large body of 3D MHD simulations (e.g. Fan, 2001b; Magara and Longcope, 2001, 2003; Magara, 2004; Manchester IV et al., 2004; Archontis et al., 2004, 2005, 2006; Galsgaard et al., 2005, 2007; Murray et al., 2006; Magara et al., 2005; Magara, 2006, 2007, 2008; Archontis and Hood, 2008; Archontis and Török, 2008; Archontis et al., 2009; Fan, 2009b) have been carried out to model the emergence of a twisted magnetic flux tube through a multi-layered atmosphere that includes a polytropic layer representing the top of the solar convection zone, an isothermal layer representing the photosphere and chromosphere, connecting to another isothermal layer representing the million-degree corona. In these simulations, a twisted flux tube is initially placed in the convection zone at a depth of a couple of Mms and the central segment of the twisted tube rises buoyantly towards the photosphere. When the top of the tube enters the photosphere the rise velocity at the upper boundary of the tube first slows down as it encounter the stable stratification of the photosphere. As a result, magnetic flux begins to pile up at the photosphere and a steep magnetic pressure gradient is established. It is found that subsequently, the magnetic flux undergoes a run-away expansion into the atmosphere due to the non-linear growth of the magnetic buoyancy instability. Archontis et al. (2004) and Murray et al. (2006) quantitatively examined their simulation data using the criterion for the onset of the undulatory magnetic buoyancy instability (Newcomb, 1961; Thomas and Nye, 1975; Acheson, 1979):

$$ - {H_p}{\partial \over {\partial z}}(\log B) > - {\gamma \over 2}\beta \delta + k_\parallel ^2\left({1 + {{k_ \bot ^2} \over {k_z^2}}} \right),$$
(31)

where is height, H p is the local pressure scale height at the photosphere, B is the magnetic field strength, γ is the ratio of the specific heats, β is the ratio of the plasma pressure over the magnetic pressure, δ ≡ ∇ − ∇ad is the superadiabaticity of the atmosphere, which is −0.4 for the isothermal stratification of the photosphere, and k , k and k z are the three components of the local perturbation wave vector (normalized by 1/H p), with k and k being the two horizontal components parallel and perpendicular to the local magnetic field direction, and k z being the component. They found that the run-away expansion takes place at the time when the above critical condition is met. Furthermore, Murray et al. (2006) found that if either the field strength or the twist of the subsurface flux tube is too low, the magnetic pressure build-up at the photosphere may not achieve the critical condition given above and flux emergence into the atmosphere may fail to take place.

During the run-away expansion of the magnetic flux tube into the solar atmosphere, strong diverging downflows are found along the emerged field lines, forming shock fronts just above the photosphere (e.g. Fan, 2001b; Magara and Longcope, 2003; Archontis et al., 2004). Due to the magnetic tension force associated with the twisted field lines, a shear flow pattern immediately develops on the photosphere with the plasma on the two sides of the polarity inversion line moving oppositely in the east-west direction (Manchester IV, 2001; Fan, 2001b; Magara and Longcope, 2003; Manchester IV et al., 2004; Archontis et al., 2004; Manchester IV, 2007). The effect of the tension force driven shear flow is to transport the axial magnetic flux upward into the expanding portion of the flux tube in the solar atmosphere (Manchester IV et al., 2004; Manchester IV, 2007). Newly developing active regions on the photosphere often exhibit such a shear flow pattern (Zirin, 1988; Strous et al., 1996). Since the subsurface flux tubes are assumed to be significantly twisted in these simulations, the two polarities of the emerging region on the photosphere initially emerge as north-south oriented, and then they undergo shearing motion along the polarity inversion line and separate in the east-west direction (e.g. Fan, 2001b; Magara and Longcope, 2003; Manchester IV et al., 2004; Archontis et al., 2004; Magara, 2006; Fan, 2009b).

It is further found (Magara, 2006; Fan, 2009b) that, as the two polarity flux concentrations become separated, a prominent rotational or vortical flow develops within each polarity, centered on the peak of the vertical magnetic flux concentration (see Figure 36), reminiscent of sunspot rotations that have been observed in many events preceding X-ray sigmoid brightening and the onset of eruptive flares (e.g Brown et al., 2003; Zhang et al., 2008).

Figure 36:
figure 36

The left panel shows the 3D coronal magnetic field produced by the emergence of a twisted magnetic flux tube from the solar interior into the solar atmosphere, resulting from a simulation of Fan (2009b). The right panel shows the z component of the vorticity w Z on the photosphere overlaid with contours of B z with solid (dotted) contours representing positive (negative) B z. It shows counter-clockwise vortical motion (i.e. positive w Z) centered on the peaks of the vertical flux concentrations of the two polarities of the emerging region. Figures reproduced with permission of the AAS.

The vortical motions are counter-clockwise (clockwise) for a left-hand-twisted (right-hand-twisted) emerging flux tube. Fan (2009b) showed that the vortical motions in the two polarity flux concentrations twist up the inner emerged field lines, causing them to rotate in the atmosphere from an initial normal configuration (i.e. arching over the polarity inversion line from the positive to the negative polarity) into an inverse configuration (directed from the negative polarity to the positive polarity over the neutral line) (see the movie associated with Figure 37).

Figure 37:
figure 37

mpg-Movie (2983 KB) Still from a movie showing 3D evolution of a set of tracked field lines as they are being twisted up and rotate in the atmosphere due to the shear and rotational motions at their footpoints on the photosphere. In these images and the movie, a field line of a particular color corresponds to the same field line carrying the same plasma. The black field line corresponds to the original tube axis, and all the other field lines have their mid points (at x = 0, y = 0) above the mid point of the black field line (reddish field lines have mid points above bluish field lines). From (Fan, 2009b). Figure and movie reproduced with permission of the AAS (For video see appendix).

In this manner, a flux rope with sigmoid-shaped dipped core fields forms in the atmosphere, with the center of the flux rope, as represented by the O-point of the magnetic field in the central cross-section, rising into the corona (see left panel of Figure 36).

Fan (2009b) showed that the rotational or vortical motions centered on the two polarity flux concentrations are a manifestation of non-linear torsional Alfvén waves propagating along the flux tube, consistent with what has been predicted by an earlier idealized analytical model of Longcope and Welsch (2000). Due to the rapid stretching of the emerged magnetic field in the atmosphere and corona, the magnitude of α = J · B/B 2 along the coronal field lines, which is a measure of the rate of twist per unit length, drastically decreases. As a result, a gradient of α is established along the field lines of the flux tube from the interior into the atmosphere with the interior portion having a much higher magnitude of (Figure 38).

Figure 38:
figure 38

The variation of α ≡ (∇ × B) · B/B 2 as a function of z along three field lines: the black field line shown in Figure 37, which is the original tube axis, and its two neighboring blue field lines at time t = 118. The α values are plotted along these three field lines (using the same colors for the data points as those of the corresponding field Figures 37) as a function of z, from their left ends to the left apices in the atmosphere. From Fan (2009b). Figure reproduced with permission of the AAS.

This gradient of α drives torsional Alfvén waves along the flux tube, transporting twist from the interior highly twisted portion into the expanded coronal portion (Longcope and Welsch, 2000). The rotational motion will continue until the coronal equilibrates with the interior α along the field lines. The time scale for establishing the equilibrium is on the order of the Alfvén transient time along the interior flux tube, which means that the rotational motion can persists for a few days after the initial emergence. Magara and Longcope (2003) and Fan (2009b) found that the helicity flux due to flux emergence is dominating only for a brief initial period, and then horizontal shearing and rotational motions at the footpoints of the emerged fields provide the dominant and steady source of helicity flux into the atmosphere. These results suggest that the observed sunspot rotations are due to non-linear torsional Alfvén waves naturally occurring during the emergence of a twisted flux tube, and is an important means whereby twist is transported from the interior into the corona, driving the development of a coronal flux rope as a precursor structure for solar eruptions. The horizontal vortical motion and its subsurface extension corresponding to the torsional Alfvén waves may be detectable by local helioseismology techniques (see review by Gizon and Birch, 2005).

The simulations consistently show that the subsurface twisted flux tube does not rise bodily into the corona as a whole due to the heavy plasma that is trapped at the bottom concave (or U-shaped) portions of the winding field lines. While the upper parts of the helical field lines of the twisted flux tube expand into the atmosphere and the corona, the U-shaped parts of the winding field lines remain largely trapped at and below the photosphere layer, and the center of the original tube axis ceases to rise a couple of pressure scale heights above the photosphere. Nevertheless in the end, a twisted coronal flux rope with sigmoid-shaped concave upturning field lines threading under a central axis that rises into coronal heights is found to develop in many simulations (Magara and Longcope, 2001, 2003; Magara, 2004; Manchester IV et al., 2004; Magara, 2006; Manchester IV, 2007; Archontis and Hood, 2008; Archontis and Török, 2008; Archontis et al., 2009; Fan, 2009b). By Larangian tracking of the evolution of the emerged field lines as their footpoints are undergoing shearing and vortical motions (see Figure 37), Fan (2009b) found that the inner emerged field lines (including the original tube axis) rotate in the atmosphere into an inverse configuration, and that is how the sigmoid-shaped, dipped core fields of the new coronal flux rope come into being. Driven by the continued vortical motions at the footpoints of the emerged field lines, the newly formed coronal flux rope accelerates upwards, and a current sheet of an overall sigmoid morphology develops in the lower atmosphere, extending up to the base of the corona (Manchester IV et al., 2004; Magara, 2004, 2006; Archontis and Török, 2008; Archontis et al., 2009; Fan, 2009b). Figure 39 shows such an example.

Figure 39:
figure 39

Emerging magnetic field in the solar atmosphere resulting from the 3D simulation of the emergence of a left-hand-twisted magnetic flux tube by Magara (2004). The colors of the field lines represent the square value of the current density at their footpoints on a chromospheric plane located at z = 5. Top left: Top view of the magnetic field lines. Note the inverse-S shape of the brighter field lines, which is consistent with the X-ray sigmoid morphology preferentially seen in the northern hemisphere. Top right: The square of the current density (color image) and vertical magnetic flux (contours) at the chromospheric plane. Bottom left: Side view of the magnetic field lines. Bottom right: Another perspective view of the magnetic field.

It is found that the field lines going through the current sheet all show a sigmoid-shape. Thus the heating associated with the current sheet may cause these sigmoid-shaped field lines to preferentially brighten up in soft X-ray, giving rise to the observed X-ray sigmoid loops in an active region (e.g. Rust and Kumar, 1996; Canfield et al., 1999).

Archontis et al. (2009) further examined the detailed evolution of the (overall) sigmoid-shaped current layers in their simulation of an emerging twisted flux tube, and found that it agrees qualitatively with the complex features and evolution of a coronal sigmoid observed by the XRT of Hinode, which is characteristic of many evolving sigmoids that result in eruptions (see Figure 40).

Figure 40:
figure 40

Comparison between the results from a simulation of emerging flux tube (left and middle columns) and the XRT/Hinode observations (right column). The left column shows the evolution of the constant current surfaces, the middle one shows the result computed from a heating proxy, and the right column shows XRT images at three different times during the evolution of the sigmoid structure. From Archontis et al. (2009). Figure reproduced by permission of the AAS.

The simulation shows that earlier in the emergence, the current layers and the associated baldpatch separatrix field lines form two “J”-like structures. With time the electric current becomes more rich in structure, with additional thin current layers forming. The current layers and the associated reconnected field lines form an overall sigmoid-shaped structure. Then in the last phase, a central current sheet develops in the middle of the sigmoid structure, accompanied by the rapid rise of the coronal flux rope and the appearance of post-reconnection loops under the current sheet. Observations of evolving coronal X-ray sigmoids often show a similar transition from an initial morphology of two “J”-like bundles, to the development of a central brightening in the middle of the sigmoid at the onset of an eruptive flare.

The effects of the presence of a simple horizontal coronal magnetic field and the interaction and reconnection of the pre-existing coronal field with the emerging flux tube have been investigated in a series of work (Archontis et al., 2004, 2005, 2006; Galsgaard et al., 2005, 2007; Archontis and Török, 2008). It is found that the dynamic interaction and the rate of magnetic reconnection depends sensitively on the relative orientations between the upcoming emerging flux and the pre-existing coronal magnetic field. When the two flux systems coming into contact have a relative angle above 90°, the simulations show immediate and substantial magnetic reconnection, producing collimated high-speed and high-temperature jets from the reconnection site. However, the cases that have a more parallel orientation of the flux systems show very limited reconnection and none of the associated features. Most recently, Archontis and Török (2008) found that rapid reconnections with a pre-existing horizontal coronal magnetic field that is nearly anti-parallel with the top of the emerging flux system can greatly enhance the eruption of the newly formed coronal flux rope. The reconnections remove flux from both the pre-existing field and the outer emerged fields, which are arcade-like fields acting to confine the sigmoid core fields of the newly formed coronal flux rope. Thus the coronal flux rope shows an enhanced acceleration, reaching a maximum speed of about 240 kms−1. In contrast, in the case of emergence into a field-free atmosphere, the coronal flux rope that forms generally reaches a maximum rise speed of a few tens of kms−1 (e.g. Manchester IV et al., 2004; Archontis and Török, 2008; Fan, 2009b).

Isobe et al. (2005) performed high resolution 3D simulations of emerging flux and its reconnection with pre-existing magnetic field in the corona. In this calculation the emerging flux develops from an initial horizontal magnetic flux sheet (with uni-direction field) situated in the top of the convection zone. Due to the growth of the three-dimensional magnetic buoyancy instability of modes with high wave number in the direction perpendicular to the field (or magnetic Rayleigh-Taylor instability), undulating flux bundles with fine spatial scales rise into the corona and reconnect with the pre-existing coronal magnetic field in a spatially intermittent way. This results in the formation of filamentary structure resembling the observed arch filament systems in emerging flux regions (EFRs). The spatially intermittent reconnection and heating also explain the coexistence of many hot and cold loops and the jets being ejected from the loop footpoints in EFRS observed in EUV by the TRACE satellite.

Due to the need to resolve the photosphere pressure scale height (∼ 150 km), the 3D simulations of the emergence of twisted flux tubes from the interior into the solar atmosphere and the corona described above are done for a domain size of up to a few tens of Mm, significantly smaller than the size of a typical active region. The implicit assumption is that the qualitative dynamical behavior of the smaller emerging flux tubes modeled in these simulations is representative of that of the larger active region scale emerging tubes. Furthermore, these simulations, which focus on understanding how twisted active region flux tubes emerge dynamically into the atmosphere and the corona, typically ignore convection in the interior layer. They also grossly simplify the treatment of the thermodynamics by assuming an ideal gas equation of state and an adiabatic evolution of the plasma, which are not appropriate for the atmosphere layers considered. Thus these simulations of flux emergence into the atmosphere do not study the magneto-convection process of sunspot formation, which depends critically on the convective and radiative energy transport at the photospheric layer. The assumed adiabatic expansion of the tube plasma emerging into the atmosphere results in unphysically low temperatures. Radiative heat exchange, thermal conduction, and the still uncertain process of coronal heating driven by the mechanical energy of convective motions, all play an important role in the thermal energy evolution of the plasma in emerging flux regions, and in maintaining the observed temperature profile of the atmosphere layers.

Progress is being made towards building a full radiation MHD model of active region flux emergence, encompassing both the magneto-convection process of the formation of sunspots from rising flux tubes and the emergence of the active region flux into the stably stratified solar atmosphere and the corona with realistic treatment of the thermodynamics. Prototypes of such numerical models include, for example, Abbett (2007), Amari et al. (2008), and Martínez-Sykora et al. (2008, 2009). The simulations of Martínez-Sykora et al. (2008, 2009) are at present the most sophisticated in terms of the physics included. They include realistic magneto-convection in the interior layer, with a realistic treatment of radiative energy transport in the interior and the overlying atmosphere layers, and taking into account thermo-conduction along magnetic field lines. The hot corona is self-consistently maintained by the magnetic energy dissipation driven by the stresses applied to the fields due to the photospheric convective motions (Gudiksen and Nordlund, 2002). In these simulations, a quiet-sun 3D domain is first initialized which contains a convecting interior layer that realistically represent the top layer of the solar convection zone, a photosphere and chromosphere, a transition region, and a hot corona, and with a pre-existing magnetic field of varying mean strength treading through the layers. A horizontal flux tube with varying degree of twist and field strength is then transported though the lower boundary into the domain, and its subsequent rise through the convection zone and the emergence into the atmosphere layers is studied. The domain sizes considered in these simulations are similar to the idealized simulations described above. The rise of the twisted flux tube through the interior layer and its emergence at the photosphere are in agreement to what have been found in Cheung et al. (2007). Similar to previous findings, a large amount of horizontal flux is retained at the photosphere, and with greater twist and/or greater field strength, more flux crosses the photosphere into the chromosphere and the corona. Compared to previous ideal MHD simulations which ignored convection in the interior layer, the rising flux tube is more fragmented, and the dynamic expansion of magnetic flux into the upper chromosphere and corona is more patchy and slower. Plasma of chromosphere temperature is found to be lifted into great heights (far greater than that expected for hydrostatic equilibrium) due to magnetic support by the emerging fields. A low density structure of the emerged fields is found to form in the corona with an underlying high density structure of chromospheric plasma supported by the magnetic fields. These are suggestive of the formation of a filament, although the scales of the structures in the simulations are far smaller than their realistic counterparts (due to the limitation in numerical resources which prevent this type of simulations from being carried out on realistic scales of active regions and filaments). Associated with the high density structure, low density inclusions or voids are found to rise in the high density structure as shown in the synthetic Ca II images at the limb, which resembles the dynamics observed in quiescent prominences with the Hinode Spacecraft (Berger et al., 2008). New results are also found with regard to the spatial temporal evolution of Joule heating during flux emergence, and how it affects the energetics of the chromosphere and corona, and changes the emission of various transition region and coronal lines. Further advancement in computing power in the next decade is likely to enable such realistic radiation MHD simulations of flux emergence on a domain size-scale that is approaching a realistic active region.

8.3 Post-emergence evolution of subsurface fields

This is so far a largely unexplored area with mostly speculations and very few quantitative calculations. One outstanding question is whether and to what extent sunspots and active region fields at the surface remain connected to the subsurface flux tubes in the deep solar convection zone. Clearly the answer to this question is intimately related to understanding the process of sunspot formation from rising flux tubes as described in

If the surface fields remain connected to the subsurface flux tubes, then based on the shape of the Ω-loop as shown in Figure 17, it seems difficult to understand why the separation of the two polarities of a solar active region stops at the scales of ≲ 100 Mm, without continuously increasing to at least the footpoint separation of the Ω-loop at the base of the convection zone, which is ∼ 1000 Mm! Due to the angles between the legs of the Ω-loop and the surface, the magnetic tension of the legs is expected to be pulling the two poles of the bipolar region apart until the legs become vertical. However there is also an attracting force between the two poles due to the emerged magnetic field above the photosphere. If we assume that the field above the photosphere is potential, then the attracting force of the two poles due to the field above can be estimated by considering the attraction between two opposite electric charges of ±q, where q = Φ/2π and Φ is the total magnetic flux of the active region. This yields (Φ/2π)2/d 2 for the attracting force, where is the separation between the poles. On the other hand, the lateral force from the subsurface field that acts to pull the two poles apart can be estimated by integrating over the area of each pole the Maxwell stress exerted by the subsurface field from below. This leads to an estimated lateral force of \((B_{n}^{2}/4\pi){\rm ctg}\theta\Delta S\), where B n and ΔS are the photosphere normal field strength and area of each magnetic pole respectively (note that we have B nΔS = Φ), and θ is the angle between the slanted subsurface tube and the surface. Balancing the forces from above and below the photosphere (taking a representative value of θ = 45°), we find that it is possible for the two poles to reach a lateral force balance if πd 2 = ΔS. This suggests that the two poles of the bipolar region need to be very close, with the polarity separation no greater than the radius of each pole, in order for the attracting force from above to balance the force that drives the separation from below. Furthermore it appears that the force balance is an unstable equilibrium because of the 1/d 2 dependence of the attracting force which weakens in response to an increase in d.

One possible resolution of the above problem is that some other processes that can form Ω-loops of much shorter length scale compared to the undulatory buoyancy instability (Section 4.1) is responsible for the formation of active region emerging tubes, although it is not immediately clear what these processes might be. Another alternative is that the active region magnetic fields on the photosphere become dynamically disconnected from the interior flux tubes. Fan et al. (1994) speculated on the process of “dynamic disconnection” which has the same physical cause as that of the so called flux tube “explosion” described in Section 7. It can be seen from Figure 33, that even for an emerging flux loop with a field strength at the convection zone base that is as high as 105 G, the flux tube is expected to lose pressure confinement and hence “explode” at a height of about 10 Mm below the surface as a result of plasma establishing hydrostatic equilibrium (HE) along the tube. While a flux loop is rising, especially in the top few tens of Mm of the convection zone where the rise speed becomes comparable to the Alfvén speed, the variation with depth of the magnetic field strength can deviate significantly from the conditions of HE. However, after flux emergence at the surface, the plasma inside the submerged flux tube will try to establish HE. Tube plasma near the surface can cool through radiation and probably undergoes “convective collapse”, forming sunspots and pores (see Spruit and Zweibel, 1979; Stein and Nordlund, 2000) at the photosphere. However in deeper layers where plasma evolution is still nearly adiabatic, a catastrophic weakening of the magnetic field or flux tube “explosion” (Section 7) may occur, which is caused by an upflow of high entropy plasma as it tries to establish HE along the tube. This weakening of the field leads to the effective dynamic disconnection of the active region fields on the photosphere from the interior flux tubes. Fan et al. (1994) suggested that “dynamic disconnection” can explain

  1. 1.

    why the active region tilt persists once an active region has formed (see more recent work by Longcope and Choudhuri, 2002) and

  2. 2.

    why the 2D surface flux transport models based on diffusion by supergranular motion (see Wang et al., 1989) can explain well the observed dispersal of active regions and the migration of active region flux across latitudes.

A first quantitative calculation of the above process of “dynamic disconnection” has been carried out by Schüssler and Rempel (2005), using a 1-D self-similar vertical flux tube model whose top end has reached the photosphere. The model computes the quasi-static evolution of the flux tube under the effects of radiative cooling, convective energy transport, and a pressure buildup by a prescribed inflow at the bottom of the tube due to the high entropy tube plasma flowing upward to establish hydrostatic equilibrium along the field. The calculation shows that after emergence, the radiative losses near the surface drives an inward propagating cooling front accompanied by a downflow, which leads to a decrease of the gas pressure and an intensification of the magnetic field in the surface layers. In the mean time, the convergence of the radiative cooling driven downflow and the buoyancy driven upflow results in an increase of the gas pressure and a weakening of the magnetic field in the tube at a depth of a few Mms. It is found that for a reasonable range of the upflow speed, the magnetic field weakens to fall below the local equipartition value (i.e. dynamic disconnection takes place) at a depth between 2 and 6 Mm, in a time scale of up to 3 days. This is consistent with the time scale over which the evolution of a newly emerging active region changes from an “active phase” of growth with increasing polarity separation to a “passive phase” which can be well represented by transport by near surface flows and supergranular diffusion.

9 Summary and Discussion

9.1 Subsurface evolution of active region magnetic fields

If we believe that the magnetic field seen in sunspots and active regions on the solar surface originates from a large-scale, predominantly toroidal magnetic field generated at the base of the solar convection zone, then the process of how active region flux is transported through the convection zone and emerge into the solar atmosphere becomes a crucial link in the solar cycle dynamo puzzle. Theoretical, numerical, and observational studies in the past 3 decades have greatly improved our understanding of the process, but have also raised many new questions. One of the key questions remaining unanswered is the field strength B of the toroidal magnetic field at the base of the convection zone generated by the solar cycle dynamo. Whether B is significantly greater than or comparable to B eq ∼ 104 G, which is the field strength in equipartition with the kinetic energy density of the convective motions, will critically affect the formation, dynamic evolution, and properties of the emerging active region flux tubes. As was described in Section 5.7, for a buoyant flux tube with B ≳ (H p/a)1/2B eq ∼ 3B eq, the magnetic buoyancy force of the flux tube dominates the hydrodynamic force of the strongest downdrafts, and the flux tube can rise unaffected by convection. On the other hand, for B eqB ≲ 3 ×B eq the strong convective downdrafts can overcome the magnetic buoyancy and hence plays a significant role in affecting the dynamics and structure of the emerging Ω-shaped tubes.

9.1.1 Magnetic buoyancy dominated regime

Thin flux tube simulations of emerging flux loops in the solar convective envelope (Section 5.1) suggest that the toroidal magnetic fields at the base of the solar convection zone is in the range of about 3 × 104 to 105 G, significantly higher than the equipartition field strength, in order for the emerging loops to be consistent with the observed properties of solar active regions. In this field strength range, it is found that Coriolis force acting on rising Ω-shaped loops produce several asymmetries between the leading and the following sides of the emerging loops, providing explanations for the observed asymmetries between the leading and the following polarities of bipolar active regions, e.g. the active region tilts described by Joy’s law and the observed asymmetric proper motions of the two polarities of newly developing active regions.

At B ∼ 105 G, toroidal magnetic fields stored in the weakly subadiabatic overshoot region are preferably in an equilibrium state of neutral buoyancy, with the magnetic curvature force balanced by the Coriolis force due to a prograde toroidal flow (Section 3.1). This is true regardless of whether the field is in the state of an extended magnetic layer or isolated flux tubes. Detecting this prograde toroidal flow through helioseismology may be a way to probe and measure the toroidal magnetic fields stored in the tachocline region. Isolated toroidal flux tubes stored in the weakly subadiabatic overshoot region should experience a radiative heating due to the non-zero divergence of the radiative heat flux (Section 3.2). This radiative heating causes a quasi-static upward drift of the toroidal flux tubes. In order to maintain toroidal flux tubes in the overshoot region for a time scale comparable to the solar cycle period, a rather strong subadiabaticity of δ < −10−4 is needed, which is significantly more subadiabatic than the values obtained by most of the overshoot models based on the non-local mixing length theory. This strong subadiabaticity may be achieved as a result of some level of suppression of convective motions by the toroidal magnetic flux tubes themselves. Furthermore, a semi-analytical model of convective overshoot (Rempel, 2004) — based on the assumption that the convective energy transport is governed by coherent downflow plumes — shows that the overshoot region can have a subadiabaticity of δ ≲ −10−4 if the downflow filling factor at the base of the convection zone is ≲ 10−4.

Neutrally buoyant toroidal flux tubes stored in the weakly subadiabatic overshoot region are subject to the onset of the undulatory buoyancy instability depending on the field strength and the value of the subadiabaticity (Section 4.1). The toroidal flux tubes become buoyantly unstable and develop Ω-shaped emerging flux loops when the field strength becomes sufficiently large, or when the quasi-static upward drift due to radiative heating brings the tubes out to regions of sufficiently weak subadiabaticity or into the convection zone. A neutrally buoyant equilibrium magnetic layer is also subject to the same type of undulatory buoyancy instability and 3D MHD simulations show that arched buoyant flux tubes form as a result of the non-linear growth of the instability (Section 4.2).

There are several major difficulties associated with explaining active regions as rising flux tubes in the magnetic buoyancy dominated regime (with B ∼ 105 G at the base of the convection zone) as listed in the following. (1) Recent solar cycle dynamo models which take into account the dynamic effects of the Lorentz force from the large-scale mean fields have suggested that the toroidal magnetic field generated at the base of the solar convection zone is ∼ 1.5 × 104 G (Rempel, 2006b). (2) The long length scales of the Ω-tubes required for the onset of the magnetic buoyancy instability are too large compared to the observed longitudinal extent of the active regions (see Section 8.3). (3) The high twist rate required to counteract the vorticity generation by the magnetic buoyancy, which tends to break up the rising flux tube, is inconsistent with the observed twists and tilts in the majority of solar active regions on the surface (Sections 5.4 and 5.5). Possible solutions for resolving the above difficulties have been proposed. For (1), for example, the amplification of a toroidal magnetic field by conversion of potential energy associated with the superadiabatic stratification of the convection zone may be a means to reach field strength that is significantly above the equipartition value (Section 7). For (2), a possible solution may be the dynamic disconnection of surface active regions from their parent tubes in the deep interior (Section 8.3). For (3), further 3D simulations of the formation and rise of buoyant flux tubes from toroidal magnetic fields initially in mechanical equilibrium in the presence of solar rotation are necessary. Another possible solution for (3) is the much lower rate of twist of the emerged magnetic field in the solar atmosphere compared to that in the interior portion of the tube as a result of the rapid and extreme stretching of the field in the atmosphere (see Figure 38 and the discussion in Section 8.2).

9.1.2 The effects of convection

The range of field strengths of 104 G ≲ B ≲ 3 × 104 G for toroidal flux tubes at the base of the solar convection is an interesting regime where the strong convective downdrafts can overcome the magnetic buoyancy of the tube, and hence play a significant role in affecting the formation, dynamics, and structure of the emerging Ω-tubes (Section 5.7). Some simulations have found that magnetic fields of BB eq are preferentially transported downward against their magnetic buoyancy out of the turbulent convection zone into the stably stratified overshoot region by compressible penetrative convection (Section 6). This is also the range of field strength for the toroidal magnetic field suggested by the recent dynamic mean field dynamo model for the solar cycle (Rempel, 2006a,b). The rise of flux tubes with these field strengths through the solar convective envelope, under the influence of both the Coriolis force and the turbulent convective flows has not been well studied.

Previous thin flux tubes simulations (Section 5.1) which neglect the effect of convection, have shown that tubes at these field strengths tend to be significantly deflected poleward by the Coriolis force during their rise, and thus have difficulty reproducing the emergence of active regions at low latitudes. They are also found to produce tilt angles that are inconsistent with Joy’s law. Therefore, an important question is whether the the convective flows in the solar convective envelope can modify the dynamics of the rising Ω-tubes such that they emerge with properties that are consistent with the observed properties of solar active regions. At these field strengths, the convective downdrafts are capable of pinning down the flux tube while the broad helical upflows in between the downdrafts can boost the rise and thus possibly produce emerging Ω-loops which are consistent with the properties of the majority of solar active regions.

Some of the major difficulties associated with explaining active regions as rising flux tubes in the magnetic buoyancy dominated regime (with B ∼ 105 G at the base of the convection zone) discussed above may be ameliorated. For such equipartition or weakly super-equipartition field strengths, the length scale of the Ω-tubes may be largely defined by the separations of convective downdrafts extending into the deep convection zone (see e.g. Miesch et al., 2008), which are comparable to the size scale of solar active regions. Furthermore, significant twist of the tube may no longer be required to obtain reasonably cohesive emerging tubes, since the generation of strong vortex tubes by magnetic buoyancy is suppressed by the tension force resulting from the pinning-down of the flux tube at short intervals by the convective downdrafts (Section 5.7). These possibilities need to be examined by self-consistent 3D spherical-shell simulations of rising flux tubes in a rotating solar convective envelope in the presence of convective flows and the associated large-scale mean flows. Important initial steps in this area have been carried out by Jouve and Brun (2009), although the simulations presented so far have considered buoyant toroidal flux tubes whose initial field strengths and fluxes are too large compared to the values expected for active region scale flux tubes.

9.2 The twist of active region magnetic fields

Vector magnetic field observations on the photosphere show that on average solar active regions have a small but statistically significant mean twist that is left-handed sense in the northern hemisphere and right-handed sense in the southern hemisphere. The origin of this hemispheric dependent twist in solar active regions is not clear. The twist may be due to the current helicity in the dynamo generated, predominantly toroidal magnetic field at the base of the convection zone, from which buoyant flux tubes form, or it may be acquired during the rise of the flux tubes through the convection zone as a result of buffeting by the helical convective motions (called the Σ-effect) and also through the accretion of the mean poloidal magnetic field in the convection zone onto the rising flux tube (Section 5.3). Observational studies of the correlation between active region twist and tilt angles indicate that the Σ-effect cannot be the only source for the twist (Section 5.3).

It is found that in order for the tilt of the Ω-shaped emerging tube to be consistent with Joy’s law of active region mean tilt, the magnitude of the twist in the rising flux tube cannot be too high (Section 5.5). On the other hand, the rise of highly twisted, kink unstable magnetic flux tubes can produce kinked or knotted emerging tubes which may explain the origin of the unusual class of flare productive active regions called δ-sunspots which are observed to be highly twisted and often show polarity order inverted from the Hale polarity rule (Section 5.6).

9.3 Active region flux emergence into the atmosphere

The evolution of active region scale flux tubes in the top layer of the solar convection zone, the process of sunspot formation at the surface, and the post-emergence evolution of subsurface flux tubes are not well understood (Section 8). However, rapid progress has been made in recent years in realistic radiation MHD simulations of magneto-convection in the top few Mm layer of the solar convection zone and the overlying photospheric layers. These simulations are able to reproduce many observed features seen in high resolution photospheric observations of the solar granulations, emerging flux regions, and sunspot fine structure (Section 8.1).

With regard to flux emergence into the solar atmosphere and the corona, it has been shown that the magnetic buoyancy instability is a mechanism through which magnetic flux reaching the photosphere can expand dynamically into the stably stratified solar atmosphere. Recent 3D MHD simulations of magnetic flux emergence into the solar atmosphere and reconnection of emerging flux with a pre-existing coronal field are able to reproduce many observed features in newly emerging active regions (Section 8.2).

These simulations suggest that a twisted subsurface flux tube does not rise bodily into the corona as a whole due to the heavy plasma that is trapped at the bottom concave portions of the helical field lines. Shear and rotational flows on the photosphere driven by the Lorentz force of the twisted flux tube during flux emergence are the crucial means whereby twist is transported from the interior into the solar corona, leading to the formation of a coronal flux rope with sigmoid-shaped, dipped core fields. The models suggest that sunspot rotations are a manifestation of nonlinear torsional Alfvén waves propagating along the emerging flux tube, transporting twist from the tube’s interior portion, where the rate of twist is high, towards the greatly stretched coronal portion, where the rate of twist is low. Driven by the continued vortical motions at the foot points of the emerged field lines, the newly formed coronal flux rope accelerates upwards and a current sheet of an overall sigmoid morphology develops. Heating associated with current sheet may cause the sigmoid-shaped field lines going through the current sheet to brighten up as the observed X-ray sigmoid loops in an active region (Section 8.2).