Coupling the solar dynamo and the corona:
wind properties, mass and momentum losses during an activity cycle
We study the connections between the sun’s convection zone and the evolution of the solar wind and corona. We let the magnetic fields generated by a 2.5D axisymmetric kinematic dynamo code (STELEM) evolve in a 2.5D axisymmetric coronal isothermal MHD code (DIP). The computations cover an 11 year activity cycle. The solar wind’s asymptotic velocity varies in latitude and in time in good agreement with the available observations. The magnetic polarity reversal happens at different paces at different coronal heights. Overall sun’s mass loss rate, momentum flux and magnetic braking torque vary considerably throughout the cycle. This cyclic modulation is determined by the latitudinal distribution of the sources of open flux and solar wind and the geometry of the Alfvén surface. Wind sources and braking torque application zones also vary accordingly.
Subject headings:solar dynamo — solar wind — Sun: corona — Sun: magnetic fields — Sun: convective zone
The sun’s magnetic field varies in time following a -year cycle. The most visible manifestation of this periodic behaviour is the well-known years sunspot cycle, the sunspots themselves corresponding to high local concentrations (up to about wide) of vertical magnetic field at the photosphere. The sun’s photosphere is in reality permeated by magnetic flux at several scales, and at all latitudes and longitudes. This fact brings forward the idea that the magnetic field plays a role in coupling the dynamics of the inner layers of the sun and its atmosphere. Studies about the convection zone (hereafter CZ) and atmospheric magnetic fields have mostly been done separately, nonetheless. The first reason for this is that the photosphere essentially separates two regions of plasma with different regimes ( and ) and disparate characteristic time and length scales, making any numerical investigations comprising them both a real challenge. The second reason is that the actual coupling and/or cross-transport phenomena are still poorly understood (except perhaps at very small-scale, or in simple scenarios). Finally, the third reason is linked to the use of different observational techniques for gathering data relating to the dense layers below the optically thick photosphere and to the essentially optically thin and rarefied atmosphere (see review by Zurbuchen, 2007).
The most remarkable exceptions are the numerical studies made in small cartesian domains which include the photospheric layers in the domain, the typical numerical domain size being of order of a few tens of (Rempel et al. 2009; Leenaarts et al. 2009; Martínez-Sykora et al. 2008; Vögler et al. 2005, among others; see also review by Wedemeyer-Böhm et al. 2009). These studies are unable, though, to capture the slowly varying magnetic field and solar wind’s properties at a global scale.
The structure of the magnetic field inside the sun is unreachable by direct observation and/or measurement. Helioseismology techniques have been used, though, to deduce some of the properties of the magnetic field in the sun (Kosovichev, 2006; Kosovichev & Duvall, 2006; Antia et al., 2000, 2003). On the other hand, dynamo models are able to reproduce photospheric observables in greater detail.
Currently our understanding of the inner solar magnetism, its structure, evolution and origin rely on multi-D models of solar dynamo. This “fluid” dynamo is at the origin of the intense magnetic activity of the sun. In particular, the solar interface dynamo paradigm has received much attention (Parker, 1993). This model assumes that the locations of generation of the toroidal and poloidal global field are separated, with the tachocline playing a central role in organising toroidal field.
The classical explanation for the cyclic activity of the large scale magnetic field is that a dynamo process acts in the solar interior to regenerate the three components of the magnetic field and sustain them against ohmic dissipation. The inductive action of the complex fluid motions would thus be responsible for the vigorous regeneration of magnetic fields and for its nonlinear evolution in the solar interior (see Charbonneau 2010 and Miesch 2005 for recent reviews on the subject). Understanding how these complex physical processes operating in the solar turbulent plasma non-linearly interact is very challenging. One successful and powerful approach is to rely on multi-dimensional magnetohydrodynamic (MHD) simulations. In this context, two types of numerical experiments have been performed since the 70’s: kinematic mean-field axisymmetric dynamo models which solve only the mean induction equation (Steenbeck & Krause, 1969; Roberts, 1972; Stix, 1976; Krause & Radler, 1980) and full 3D global models which explicitly solve the full set of MHD equations (Gilman, 1983; Glatzmaier, 1985; Cattaneo, 1999; Brun et al., 2004).
Clearly, both approaches are complementary and are needed to better understand the magnetic solar activity. Recent progress have been made with 3D numerical models of magnetic stars. Large-scale magnetic cycles are starting to be found in simulations of the Sun, such as those performed by Ghizaru et al. (2010), or of solar-like stars (Brown et al., 2011). However, butterfly diagrams similar to the solar observations — with the poloidal field field reversal happening when the toroidal field is at its maximum — are still difficult to reproduce. Kinematic mean-field models and their associated simplifying assumptions have thus been used extensively to reproduce several features of the large-scale solar cycle. In particular the use of a differential rotation profile inferred from helioseismology associated with an alpha-effect (due to the helical turbulence of the stellar convective envelope) antisymmetric with respect to the equator enabled Charbonneau & MacGregor (1997) to produce a solar-like butterfly diagram with the ingredients of the modern interface dynamo. A similar model will be considered in this work to catch the large-scale behaviour of the inner solar magnetic field.
Above the photosphere, a very complex structure comprising open magnetic flux-tubes (coronal holes) and magnetic loops with different length-scales arises (magnetic carpet, canopy, coronal loops) (Stix, 2002; Aschwanden, 2005, and references therein). The wealth of ground and space based observational data gathered in the latest years has been bringing up a great insight on the magnetic structure in the lower atmospheric layers (chromosphere and lower corona). The properties of the magnetic field in the (optically thin) higher corona are harder to determine observationally, though. Space-borne in-situ measurements may introduce further constraints about how the solar wind properties connect to the coronal structures. Remarkable examples of these were given by the Ulysses consecutive polar orbits (McComas et al., 2003; Issautier et al., 2004, 2008). The SoHO spacecraft further provided complementary in-situ and solar surface (and low coronal) observations. Future missions such as the Solar Orbiter should provide more refined data to complete the scenario. These should combining high resolution imagery and spectroscopic data, magnetogram and in-situ measurements of the wind’s properties in an orbit with very close helio-synchronous passes and varying orbital inclinations.
The solar wind’s outflow properties depend on the particular geometry of the coronal field at each moment of the solar cycle. The amplitude and latitudinal distribution of the solar wind velocity and mass flux depend on parameters such as the positions of the wind sources at the surface, the local magnetic field strength and the expansion factor of each particular coronal flux-tube the wind flows along. The details about the coronal heating mechanisms, most notably the amplitude and location of energy dissipation, are also of importance. Leer & Holzer (1980) pointed out that heating below or above the critical sonic point may produce different effects on the wind velocity, and that high speed winds require energy deposition in the supersonic region. Hansteen & Leer (1995) showed that the mass flux does not depend strongly on the mode and location of energy deposition, but rather on the amplitude of the energy flux. Note also that the heat flux profile could be related to the magnetic field’s amplitude and geometry (Cranmer et al., 2007; Pinto et al., 2009). The resulting angular momentum losses and wind’s braking torque (applied on the sun) depend on wind’s velocity and mass flux together with the (time varying) geometry of the Alfvén surface (the surface at which the wind’s velocity equals the Alfvén velocity; see definitions in §5). Some authors have studied the influence of the topology of coronal field in the properties of stellar winds and resulting braking torques, but mostly by using simplified configurations (e.g, dipolar versus quadrupolar field, as in Matt & Pudritz 2008). Others have furthermore considered how the presence of strong magnetic spots (localised magnetic flux enhancements) specifically affects the angular momentum loss rate (e.g, Aibéo et al., 2007; Cohen et al., 2009).
More detailed studies of the global solar magnetic structure usually adopt surface magnetogram data for the radial component of the field as a lower boundary. The atmospheric magnetic field’s geometry is deduced using potential field extrapolation techniques. The coronal and heliospheric hydromagnetic conditions are then deduced either by using a set of semi-empirical relations (Arge & Pizzo, 2000; Luhmann et al., 2002; Wang & Sheeley, 1992; Schrijver & DeRosa, 2003), or by finding stable MHD solutions for the mapped fields (e.g, Usmanov, 1993; Hu et al., 2008). The usage of magnetogram data as lower boundary conditions encounters a few difficulties. Namely, only the line-of-sight component of the field is available, making the estimation of the radial components difficult at the polar regions. The standard PFSS technique (standing for “potential field source surface”) further assumes the surface field is strictly radial. One still needs to assume an outer boundary condition which imposes the magnetic field to be purely radial as well at a source-surface (typically, a spherical surface placed at about ). This condition emulates the field-line opening usually caused by the solar wind flow, while keeping the global field at its lowest (current-free) energy state. These studies essentially provide quasi-stationary predictions of the coronal fields at a given moment. Despite all the simplifications and assumptions this family of models has been successful in predicting qualitatively the coronal magnetic field’s topology (observed, for example, in white light during solar eclipses) and became a standard to which other types of model can be compared. Their MHD counterparts are limited in different ways. Beyond the much increased computational power they require, they rely on phenomenological assumptions for the heat transport in the corona (as the heating processes are still under debate). These models show, nevertheless, some quantitative differences when compared to the PFSS ones. The relative sizes of the streamers can be different (usually more elongated) and the heliospheric magnetic field tends to be more uniform (as it is known to be). See, for example, Riley et al. (2006) for a comparison between these two families of models. A third type of models — nonlinear force-free field (NLFFF) models — places itself somehow in between the two previously described, in the sense that it tries to add some complexity to the PFSS scenario while avoiding the issues found by global MHD computations. These require, though, good quality vector magnetogram data at the surface of the sun and still find some difficulty in accurately predicting real solar features (DeRosa et al., 2009).
In this paper, we investigate the influence of the cyclic evolution of the large-scale magnetic field produced by the solar dynamo on the solar wind properties. We chose not to perform this study based on surface magnetogram magnetic field data (as in Wang et al., 2003; Hu et al., 2008; Luhmann et al., 2009, for example), but rather pursue a different path: our source fields are the result of well-tested and validated kinematic dynamo models (Jouve et al., 2008). An isothermal MHD model of the solar wind and corona is used to produce a temporal sequence of steady-states spanning a complete activity cycle. The focus here is on estimating the wind’s velocity, mass and angular momentum flux spatial profiles as they vary during an activity cycle in response to the variations of the magnetic field’s topology, rather than reproducing a particular solar cycle. The variability of global properties such as the sun’s mass loss rate and wind’s breaking torque are studied in regards to the activity cycle. Also, we expect to gain a deeper insight on the connections between the sub-surface and coronal physical processes by proceeding this way.
The remainder of this manuscript is organised as follows: the methods and numerical codes used are presented in §2. The results are described thoroughly in §3 (coronal magnetic field evolution), §4 (solar wind speed) and §5 (mass and angular momentum flux, breaking torque). A discussion follows in §6 and a brief summary in §7.
2. Coupling solar dynamo and wind models
We used two 2.5D axisymmetric MHD codes. The first one — STELEM (Jouve & Brun, 2007, and §2.1) — computes mean-field kinematic MHD dynamo solutions given the meridional circulation and differential rotation profiles in the solar convection zone. A Babcock-Leighton or an source terms generate the poloidal magnetic field. The second code — DIP (Grappin et al., 2000, and §2.2) — computes the temporal evolution of an MHD solar corona with a self-consistent wind. In the latter, the magnetic field separates into two components: an imposed stationary external potential field (whose sources lie within the sun), and an induced component produced by the flows within the numerical domain. The dynamo magnetic fields produced by STELEM match a potential field at the surface of the sun. This potential source field can then be directly transmitted to DIP as the external component of the coronal magnetic field. Further details about both numerical codes are given hereafter, in §2.1 and §2.2. Details about the coupling of the two codes are given in §2.3.
2.1. STELEM: a representative dynamo solution
To investigate the global solar cycle features produced by dynamo models, we start from the hydromagnetic induction equation, governing the evolution of the magnetic field in response to advection by a flow field and resistive dissipation.
As we are working in the framework of mean-field theory, we express both magnetic and velocity fields as a sum of large-scale (that will correspond to the mean field) and small-scale (associated with fluid turbulence) contributions. Averaging over some suitably chosen intermediate scale makes it possible to write two distinct induction equations for the mean and the fluctuating parts of the magnetic field.
A closure relation is then used to express the mean electromotive force in terms of the mean magnetic field, leading to a simplified mean-field equation (Moffatt, 1978). The prescribed mean velocity field will here only consist in its longitudinal component (the solar differential rotation) and the magnetic diffusivity is assumed to be constant. The source term for the poloidal field is linked to the turbulent helical motions within the convection zone. We thus obtain a simple dynamo model with constant diffusivity.
Working in spherical coordinates and under the assumption of axisymmetry, we write the total mean magnetic field B as
Reintroducing this poloidal/toroidal decomposition of the field in the mean-field induction equation, we get two coupled partial differential equations, one involving the poloidal potential and the other concerning the toroidal field
where , is the differential rotation and is the -effect.
All these quantities are now dimensionless, thanks to the presence of the two Reynolds numbers , and , being the turbulent magnetic diffusivity, a measure of the rotation rate at the equator and a measure of the intensity of the alpha-effect. The product of these two numbers gives us the dynamo number which measures the efficiency of the source terms to make the magnetic energy grow in time against Ohmic diffusion.
Equations (1) and (2) are solved in an annular meridional cut with the colatitude and the radius i.e from slightly below the tachocline (e.g. ) up to the solar surface (see Figure 1, left panel). At and boundaries, both and are set to 0. At , we compute a perfect conductor condition. At the upper boundary, we smoothly match our solution to an external potential field, i.e. we have vacuum for .
The principle of the finite element method is to look for solutions of the weak formulation of the equations, here equations (1) and (2). Those solutions are taken to be linear combinations of well-chosen trial functions. In our case, they are Lagrange polynomials of degree 1 and the elements are rectangles in the plane. Applying this spatial method to our equations results in a system of first order ODEs in time governing the evolution of the coefficients of the linear combinations. The temporal scheme that we use is adapted from Spalart et al. (1991). It is similar to a Runge-Kutta 3 method and is thus explicit and of order 3.
The STELEM code is here used to produce a cyclic dynamo field within the model convection zone, whose values at the top of our domain will be reintroduced in the DIP code described below. We chose a simple dynamo model to produce our representative solution, with physical source terms consistent with observations. The rotation profile is deduced from helioseismic inversions and the -effect is antisymmetric with respect to the equator and positive in the Northern hemisphere (in agreement with the preferred handedness of the turbulent helical convective motions). This is similar to case A of (Jouve et al., 2008). The following values were used for the various parameters: (critical value for which dynamo action just starts to compete against Ohmic diffusion), (corresponding to a value of , in agreement with observations). The spatial resolution is .
Figure 2 shows the temporal evolution of both the latitudinal and the radial field at the surface, at all latitudes. We have here the evidence of a cyclic magnetic field: the surface poloidal field changes sign regularly and at different times depending on the latitude of interest. This simple model does not aim at reproducing the solar dynamo in its finest details. The goal here is rather to assess the influence of a cyclic poloidal background field (having the same general properties as that in the sun) on the wind’s properties. The dashed white lines labelled and in the figure indicate the beginning and the end of the particular time interval used to compute the solar wind evolution. This time interval spans an entire cycle period and is roughly located between two activity minima.
2.2. DIP: wind model
The DIP code is a 2.5 D axisymmetric model of the solar corona obeying the compressible MHD equations for a one-fluid, isothermal and fully ionised plasma. The continuity and momentum equations are
We set and a uniform temperature in all the corona. The gas pressure is deduced from the equation of state
valid for a fully ionised hydrogen plasma. The isothermal approximation discards a complete treatment of the energy fluxes in the corona. But we should note firstly that the corona is nearly isothermal (at least in the first ). Secondly, assuming is a proxy to the actual thermal state of the coronal plasma, summing-up the combined effects of the thermal conduction and the (still debated) heating source. The choice of a particular value for the coronal temperature is somewhat arbitrary, though. The temperature chosen here is justified both by representing an average coronal value in the first (see, e.g, results by Hansteen & Leer, 1995; Endeve et al., 2003; Guhathakurta et al., 2006; Pinto et al., 2009, among others) and empirically as it produces correct wind solutions. As discussed later in the text, the wind’s mass flux depends both on the coronal temperature and on wind geometrical expansion factors. We focus here on the geometrical effects linked to the activity cycle and keep our reference temperature fixed throughout the simulations.
The magnetic field decomposes into a potential external component and on a component induced by the flows. The former only has poloidal components (), while the latter comprises both poloidal and toroidal components. The total magnetic field is then defined as
Furthermore, . We integrate the evolution equation for the magnetic potential
The poloidal components of are then computed as
The azimuthal component of the induced field comes directly from the induction equation (in the direction only), that is
This method guarantees that the solenoidal condition is satisfied at all times. This papers focus only on the poloidal components of the coronal field, though (as explained hereafter, in §2.3).
The diffusive terms are adapted so that grid scale () fluctuations are correctly damped. The kinematic viscosity is defined as , typically with and . The simulations presented in this paper were performed with a vanishing , so as to approach as much as possible the limit of ideal non-resistive MHD. Additional diffusive terms are used. These are implicit numerical filters (applied over , and ) which dissipate mostly at the grid scale and minimise the dissipation of large scales fluctuations (Lele, 1992). This filtering scheme allows the diffusive parameters and to be lowered while avoiding spurious Gibbs fluctuations, and allow for lower mid-scale damping (Grappin & Léorat, 2001). Note that actual kinetic dissipation should happen at scales much smaller than the grid size, anyway.
We used a grid which is uniform in latitude and non-uniform in radius. The grid’s cells radial extent is at the lower boundary and at the upper boundary.
Both the upper and the lower numerical boundaries (respectively at and ) are transparent, i.e, mass can flow through them and waves are not spuriously reflected there (see Figure 1, right panel). This is achieved by writing the MHD equations in their characteristic form. In simple terms, this consists in projecting the system of equations in terms of the primitive quantities into an equivalent system defined in terms of characteristic variables (Thompson, 1987; Vanajakshi et al., 1989; Poinsot & Lele, 1992; Roe & Balsara, 1996). The evolution equation for these variables explicitly define the time-evolution of the MHD system in terms of wave-modes propagating upwards and downwards, or equivalently, incoming and outgoing modes. We must then constraint all the incoming modes, and let the others free. Plasma is free to flow through and its properties to vary in time, the actual boundary values at each moment depending on the current state of the system. The mass fluxes and velocities (or their gradients) are, therefore, not arbitrarily set at the numerical boundaries. This is a critical feature: the actual solar wind velocities and mass flux at each point of the surface will automatically match each unique transonic and transalfvénic wind solution. There is no need to set in advance arbitrary values for the density, velocity or mass flux at the numerical boundaries. This means no spurious boundary layers will form there and that the resulting mass fluxes will not be artificially constrained. We note that the resulting surface density and mass flux do not vary by large factors (not more than in all latitudinal domain, and less than within open field regions at all moments of the cycle). Furthermore, this type of boundary conditions better describes dynamic finite-frequency phenomena than any type of non-transparent (i.e, “rigid”) conditions (Grappin et al., 2008).
The solar wind develops into a stable transonic and transalfénic solution in the open field regions after a pressure (sonic) perturbation is applied at the outer boundary. This perturbation propagates inwards until both the sonic and alfvénic surfaces appear (and are enclosed in) within the domain. From then on, the solar wind profiles evolve towards a unique stable state (Velli, 1994; Grappin et al., 1997). Note that the frontiers between open (coronal holes) and closed flux regions (streamers) are not set in advance, but they are rather a result of the competition between the magnetic tension and the dynamical pressure in the wind flow. In other words, the solar wind settles where the magnetic field is incapable of containing the outward flow and forces the magnetic field to align with the radial direction.
Other aspects of the numerical model are more thoroughly discussed in Grappin et al. (2000).
2.3. Coupling method
We used a time-series derived from a STELEM run (Jouve & Brun, 2007) which correctly produces a cyclic behaviour of the poloidal magnetic field. We used a simple dynamo model with source terms for both the toroidal field and the potential .
The time-series was scaled in order to span an -year period, sampled with a time step of months. The vertical dashed white lines labelled and in Figure 2 show the beginning and end of the time-series ( and ). These are the moments when the poloidal field is at its simplest configuration (i.e, lowest multipolar order). For simplicity, we will call these moments “activity minima”, while the magnetic field’s polarity inversion phase (between lines and in the same figure) will be called “activity maximum” hereafter. The time-series comprises equally spaced samples, but we will refer to a subset of samples ( to ) for illustrative purposes throughout the text. We let the solar wind fully develop in the domain and reach a steady state (see §2.2 for further details on the formation of the wind.) The amplitude of the external magnetic field was scaled so that the total field would match coronal amplitudes (the scaling factor being kept constant during the whole cycle). Then, for each consecutive sample, we substituted directly . The system was again let free to relax and attain a new steady state for each iteration . The stability of the relaxation method was tested so that we could be sure that each relaxed steady-state solution did not depend on the history. That is, we found no hysteresis when cycling back and forth through different states (given that we let the system relax at each stage). Different permutations of the original time-series also lead to the exact same steady-states as long as each solution’s set of parameters was kept the same. The wind solution for each instant of time depends only on the corresponding .
Some of the runs were performed at different grid resolutions, and numerical convergence was verified.
Both numerical codes are time-dependent in nature, but the result of the whole coupling procedure is a sequence of steady-state solutions. Our procedure generates a map of the poloidal coronal magnetic field and wind flow during an activity cycle rather than the dynamical evolution of a particular event.
Toroidal fields and flows will not be considered in this paper. Concerning rotation, the sun is a slow magnetic rotator, and therefore the magneto-rotational effects on the poloidal wind flow are negligible (Belcher & MacGregor, 1976).
3. Coronal magnetic field
The temporal evolution of the corona and solar wind in response to the dynamo field variations is shown in Figure 3. Only the first of the northern hemisphere are displayed at six different instants of the activity cycle, corresponding to the lines in Figure 2. That is: years. The colour-scale represents the solar wind’s velocity projected on the unit magnetic field vector in units of Mach number, that is . Orange/yellow and blue/green shades therefore trace different -field polarities in the open field regions (respectively, positive and negative). The sharp transitions between positive and negative polarities in this figure outline current sheets (note that the wind flow does not change sign across these transitions, but does).
Some elements and characteristics are observed consistently throughout the whole cycle. Higher concentrations of magnetic flux at the surface (or equivalently, of current below the surface) translate into coronal loop arcade systems. Strong flux concentrations appearing in coronal holes shape up as helmet streamers, and end in a current sheet which extends outwards. Smaller flux concentrations embedded in unipolar flux regions form nearly symmetric bipolar structures with no current sheet (as the “giant polar plumes” in Pinto et al. 2010 and the “pseudo-streamers” in Wang et al. 2007). Loop structures placed inside stagnant zones (i.e with no wind flow, as inside large streamers) mostly maintain their potential field configuration.
We start the computation at the moment when the magnetic structure of the corona is at its simplest (Figure 3, first panel, ; also first white line in Figure 2). One and only large equatorial streamer extends from the surface up to nearly , where the heliospheric current sheet starts. The frontiers between the streamer and the coronal holes cross the surface at latitudes and . The letter A indicates the latitudinal extent of the streamer. The open magnetic field has positive polarity in the northern coronal hole and negative polarity in the southern coronal hole. The streamer itself is divided into four magnetic connectivity regions around one X-type null point. The null point itself is located over the equator at . The four connectivity regions mentioned above are then the group of small equatorial loop arcades place below the null point, the group of larger equatorial arcades above the null and filling up most of the streamer, and the two groups of arcades to the north and south of the null (note that Figure 3 only shows one hemisphere, and that the system is symmetric with respect to the equator). As the solar cycle starts moving away from the minimum, new flux concentrations emerge at mid-latitudes (showing up as new groups of coronal loops inside the equatorial streamer). These new structures slowly migrate polewards, attaining the streamer boundaries at about . The equatorial streamer is disrupted at this point. A fraction of the magnetic flux will remain in the equatorial region, forming a smaller streamer (Figure 3, second panel, ). The rest of the magnetic flux reconnects and forms new plume/pseudo-streamer structures at mid-latitudes (as in Pinto et al. 2010 and Wang et al. 2007). The letter B indicates the pseudo-streamer position. New coronal holes now appear at low latitudes. The polewards progression continues at a steady pace opening up its way by reconnecting with the open magnetic flux. At about one of the magnetic arcades of each of the newly formed pseudo-streamer breaks and quickly opens up (Figure 3, third panel, , letter C). As a result, a new coronal hole with inverse magnetic polarity rapidly forms. This coronal hole will grow wider and fill up all the polar regions as the polarity inversions proceeds. The previous coronal hole’s field then closes down (Figure 3, fourth panel, ), and ends up disappearing below the surface (Figure 3, fifth panel, ). The corresponding coronal arcades (the ones closing down near the poles) are indicated by the letter D in the figure. At all traces of the previous coronal hole magnetic field have disappeared (Figure 3, last panel, ). The system is back to a state very close to its original state, but with the polarity of the magnetic field reversed (the dynamo model used here produces very regular and symmetric cycles).
Note that the polarity reversal happens quickly in the corona, even if the underlying -field evolves slowly. The wind flow is responsible for the quick “opening-up” of field lines as the magnetic flux concentrations evolve slowly at the surface, and ultimately for the change of connectivity between contiguous regions. Furthermore, at the polar axis, the magnetic field’s inversion occurs at different times at different heights. To better describe this property, the top panel in Figure 4 shows radial cuts of the polar magnetic field at the same moments as those shown in Figure 3. The bottom panel shows the temporal evolution of the polar field at different heights. The sign-switch happens first at higher altitudes and proceeds downwards. The downwards progression of the reversal of is quick ( of order of a few days) between and , but from there on it slows down considerably. The delay is between and , and about between and . The overall delay (between top and bottom) is about . This delay in the lower coronal layers is due to the slow disappearance of the coronal arcades near the poles, as shown in the last three panels of Figure 3.
Figure 5 shows as a function of at different altitudes and at different moments of the cycle. Faraway from the surface, and is nearly independent of the latitude, except in the vicinity of a current sheet. Close to the surface (beneath ), the magnetic field organises in much more complex ways. The potential component of the magnetic field dominates in the lower layers of the corona whereas the induced component largely dominates above (as the wind flow becomes stronger and approaches its asymptotic velocity). From there on, the solar wind flow and the total magnetic field essentially align with the radial direction.
The multiple current sheets shown in the intermediate panels in Figure 3 may be interpreted as a highly warped current sheet in the real non-axisymmetric corona.
At the activity minimum the open flux is restricted to large polar coronal holes (about in latitude in each hemisphere). This open flux will eventually fill all the available space at greater heights (spanning in latitude, from pole to pole), as the streamer thins out; above all the magnetic field is open. Closer to the maximum, the open flux sources are more spread in latitude but cover a smaller latitudinal extent altogether (about in latitude in each hemisphere). These multiple thin coronal holes will nevertheless grow with height, and will end up filling up all available space. Above all the magnetic field is open. In other words, the average flux-tube expansion factors will be much higher at the maximum than at the minimum.
4. Solar wind speed
We focus now on the variation of the solar wind velocity during the solar cycle, and on how it is distributed in latitude. The evolution of the coronal magnetic topology (§3) has a direct influence on the size and distribution of wind sources at the surface of sun, via changes in position and width. The local magnetic pressure dominates over the wind’s dynamical pressure at the lowest layers, and the local -field’s amplitude and inclination mostly determine whether a given field line is open or closed. Conversely, the latitudinal distribution of the solar wind faraway from the sun cannot be trivially predicted from the magnetic field’s configuration at the surface. The solar wind flows accelerates along open field lines. A fluid element of cross-section at the surface will accelerate along a magnetic flux-tube with cross-section . The final velocity profile depends on the flux-tube’s expansion factor . The expansion factor itself results from the competition between the wind’s dynamical pressure and the magnetic pressure in the corona.
The terminal wind speed and flux tube expansion factor are inversely correlated at all latitudes and at all times, agreeing with the well established Wang–Sheeley–Arge semi-empirical relation (Wang & Sheeley, 1990; Arge & Pizzo, 2000). In the initial acceleration phase (below ) the expansion factors vary considerably, both in latitude and radius. At larger radii, though, the multiple and initially independent wind flows merge into a bulk spherical outflow (, independently of the latitude). Latitudinal inhomogeneities still subsist in this flow, being a result of the evolution (acceleration) of each parcel of wind flow along its path starting from the sun’s surface.
Figure 6 shows the latitudinal distribution of the solar wind speed at (the domain’s outer boundary) during the solar cycle. The wind velocity vs. latitude diagram shows good qualitative agreement with the predictions made by Wang & Sheeley (2006) matching ULYSSES in situ measurements and the more recent multi-station IPS (interplanetary scintillation) coronal observations by Tokumaru et al. (2010). The separation between fast and slow wind is well visible in the figure. At the activity minima, fast solar wind originates essentially from high latitude regions, while the slow wind flows mostly closer to the streamer frontiers, at lower latitudes. As the activity cycle progresses from the minimum to the maximum, the slow wind expands over towards the poles and takes over most of the latitudinal domain. On the declining phase of the activity cycle, the fast wind recovers the polar regions and progressively extends towards lower latitudes, restraining the slow wind flow to the equatorial region. Some irregularities appear above this otherwise too simple scenario. Most remarkably, two short-lived wind channels appear at low latitudes between and , that is, during the polarity inversion. This can be seen in Figure 3 (panels ); note how the blue shaded wind channel appears and evolves. These flows originate in the newly formed mid-latitude coronal holes. They follow the corresponding magnetic flux-tubes, bending over from at the surface down to at the outer boundary. Although they seem to fade away quickly in Figure 6, they actually last till the declining phase of the activity cycle, but the wind flowing within these coronal holes slows down. This slowing down is due to the increasingly higher expansion factor for the coronal hole, as can be seen in Figure 3. Note how the blue shaded coronal hole expands to fill the whole hemisphere faraway from the sun while its latitudinal extent at the surface remains approximately constant (from to ). This continues while the polar closed-field regions progressively disappears; all open flux will merge into a wide polar coronal hole afterwards (). Also, the closing down of magnetic flux near the poles (during the cycle’s decay phase, at about ) is related to the appearance of fast wind flows close to the polar axis. These correspond to newly formed thin polar flux-tubes expanding almost radially, which will also merge afterwards into the wide polar coronal hole.
The wind’s velocity values presented here are expected to be lower than the values measured in situ near the Earth’s orbit. The reader should note that, on one hand, our numerical domain extends only up to , and that the wind flow has not yet reached its asymptotic velocity at this height. Nevertheless, the relative variations of in latitude should not change considerably. The wind undergoes a purely spherical expansion between and . On the other hand, this does not fully accounts for the low velocities values found, though. We could expect velocities in the order of from this model if the domain extended up to . Faster wind velocities require a complete and consistent energetic treatment, which is beyond the scope of the current paper. Besides, Pinto et al. (2009) show that the anti-correlation between terminal wind velocity and flux-tube expansion factor is still verified in self-consistent non-isothermal cases. The only strong constraint over the domain’s radial extent regarding the physical correctness of the numerical model is that it has to completely contain all critical surfaces (sonic and alfvénic; cf. §2.2), which it always does.
5. Mass flux, momentum flux and magnetic braking torque
Figure 7 shows the radial mass flux associated with the solar wind in the corona at two different instants of the solar cycle. In both cases, the net mass outflow in the polar coronal holes is higher near the streamer boundaries than closer to the poles. At the maximum of activity, thin coronal holes appear also at low latitudes (as described previously). Despite their small latitudinal extent at the surface, the associated mass flux is important when compared to that in the polar coronal holes. For an outflow with a given latitudinal extent , the actual surface area it crosses is equal to . Therefore, low latitude coronal hole are more prone to produce higher mass outflow rates. Arguably, this is a consequence of the axi-symmetrical nature of our model. Nevertheless, the real three-dimensional sun shows some degree of axi-symmetry in the distribution of the streamers and coronal holes. These low-latitude wind streams correspond to the two channels visible in Figure 6, at about .
Figure 8 shows the total mass loss rate
evaluated at the outer boundary of the numerical domain. The mass loss rate evolves in par with the activity cycle. That is, is maximal at about (during the activity maximum) and minimal during the activity minima. The amplitude of the mass loss rate varies by a factor of about in time, from at the activity minimum and at the activity maximum. This trend supports the idea that the lower latitude outflows — which appear mainly close to the activity maximum — contribute with higher net mass outflow rates, as discussed in the previous paragraph. Note that the increase in cannot be due to the variations in the velocity of the wind, as it is actually lower at almost all latitudes during the activity maxima (see Figure 6) and contributes to lowering the net mass flux. On the other hand, the variations of at the lower numerical boundary are small throughout all the cycle. The amplitude of the density fluctuations at the surface within open-field regions is always below of its average value (amounting to at most if both open and closed-field regions are considered altogether).
The way the open flux maps from the solar surface up the outer domain must therefore be the main cause for the variations found in the mass loss rate.
The solar wind outflow carries angular momentum away from the sun. The specific angular momentum flux, magnetic braking torque and spin-down time-scale are deduced from the Alfvén surface’s geometry and solar rotation rate .
The Alfvén surface is the geometric locus where the wind velocity equals the Alfvén speed . The Alfvén radius is the cylindrical distance from the rotation axis to this surface. Classical wind theory (Weber & Davis, 1967, and many others) states that the angular momentum balance problem can be simplified as follows. The plasma inside the Alfvén surface is kept in solid rotation while flowing outwards along magnetic field lines. The plasma flow then becomes super-alfvénic and its angular momentum is conserved thereafter. The physical process responsible for maintaining the solid rotation while is the magnetic tension, which works against the tendency for the field-lines to spiral backwards in the azimuthal direction (keeping them as straight as possible). Beyond the critical surface, the magnetic tension loses its efficiency.
Of course, in the real sun the transition from rigid to non-rigid rotation is a continuous and smooth one. The azimuthal magnetic field does not suddenly change from straight to spiralled. But quantitatively, it all works out as if there was a sharp transition. Following this picture, one can define an average Alfvén radius — an “effective lever arm length”. The specific angular momentum flux rate is then
The resulting torque applied on the sun is (Matt & Pudritz, 2008)
The angular momentum per unit volume of a parcel of solar wind plasma rotating with azimuthal velocity is
The angular momentum per unit volume crossing a surface element is then . Integrating over a spherical surface of radius and assuming axi-symmetry translates into
We then define the magnetic spin-down time-scale as
where is the sun’s angular momentum. We estimated
The main difficulty now lies in the definition of . Figure 9 shows the Alfvén surface at two different instants of the cycle (). This critical surface shows a regular shape for most of the activity cycle, being close to spherical at most latitudes (especially, higher latitudes). At low latitudes, though, the Alfvén surface approaches a more cylindrical shape. Some irregularities appear as inward incursions as vanishes in current sheets (there’s a small but finite outflow, so there). The most evident example of such an excursion corresponds to the equatorial streamer (e.g Figure 9, left panel). We defined here the average Alfvén radius as the average cylindrical radii of the Alfvén surface (that is, , the distance to the axis) weighted by the local mass flux crossing the surface
The sections of the Alfvén surface crossing heliospheric current sheets are discarded from the computation (that is, the inwards incursions described above). The average Alfvén radius was found to be about at the minimum and at the maximum (the time average being about ). Figure 10 (top panel) shows the average Alfvén radius as a function of time during the whole solar cycle. The Alfvén radius correlates well with the global magnetic field’s amplitude for the most part of the cycle, but shows a negative correlation with the polar surface magnetic field (compare with Figure 4)
The values found for were then used to compute the solar spin-down time-scale (eq. 11, with in eq. 10) and the specific angular momentum flux (eq. 7, using rather than ). The results are displayed in Figure 10 (middle and bottom panels) as a function of time. The spin-down time-scale due to the breaking torque exerted by the solar wind was found to vary between and (the time average being ). The quantities and show a strong negative correlation, meaning the “effective lever arm length” mostly determines the breaking efficiency. Note that the mass flux rate evolution throughout the cycle works against this trend (cf. Figure 8): the highest mass outflow rate occurs at about the activity maximum, but its effect is overwhelmed by the large variations of during the activity cycle.
We have shown that the coronal magnetic field’s topological variations due to the solar dynamo have a strong influence over the solar wind’s velocity and associated mass and momentum fluxes. The ratio of open to closed magnetic flux at the sun’s surface also varies. These elements altogether define the variations of the Alfvén radii (both in latitude and in time) and how effectively the solar wind carries angular momentum away from the sun.
The average Alfvén radius was found to be about at the activity minima and at the maximum. The global poloidal magnetic field’s amplitudes follow the same trend. Unexpectedly, though, the polar magnetic field amplitudes are anti-correlated with the global field’s amplitudes. The polar magnetic field (at the surface) happens to be stronger at the moments when the closed to open flux ratio is higher. As a consequence, the polar magnetic flux-tubes show higher expansion rates at these moments of the cycle (higher expansion factors imply stronger -field decays). This behaviour could be related to the nature of the particular dynamo model used in this work, which is a simplified one (as exposed in §2.1). We should nevertheless note that, in spite of its simplicity, the model reproduces well many solar coronal features. Namely, the topological characteristics of the poloidal magnetic field we obtained are consistent with those found by Wang et al. (2003) (using a different type of model for the surface field). Also, the latitudinal distribution of the wind velocities during the cycle we computed agree with those proposed by Wang & Sheeley (2006) (well matched by ULYSSES data; see their Figure 3) and those by Tokumaru et al. (2010) (built from IPS data collected during solar cycles).
The mass loss rate (as carried away by the wind) we computed averages out in time to , a value comparable to in-situ measures near the Earth (giving a standard value of about if the wind flow was spherically symmetric). Conversely, the solar spin-down time-scale due solely to the torque exerted by the magnetised wind outflow averages out to value of the order of , above the standard .
One reason for this disparity may lie in the fact that we used a one-fluid isothermal MHD wind model. Considering a proper (i.e, non-isothermal) treatment of the energetics of the wind could lead to different quantitative estimations on the wind velocities. For example, Pinto et al. (2009) found a broad range of wind velocities in their configurations including a chromosphere – corona transition region and a coronal heating flux dependent on the local magnetic field’s strength. The actual asymptotic wind velocities depended strongly on the radial distribution of the mechanical heat flux dissipation. Under the right conditions, a non-isothermal two-fluid approach can also produce higher bulk velocities (see Endeve & Leer, 2001; Grappin et al., 2011). An additional word of caution should be given in respect to the choice of a constant and uniform coronal temperature. It is known that the mass flux in an isothermal wind flowing within a given magnetic flux-tube depends on the coronal temperature as
(cf. Parker, 1964; Leer & Holzer, 1979; Hansteen & Leer, 1995, among others) where is the density at the base of the flux-tube, a constant and a geometrical factor describing the radial expansion of the flux-tube. The absolute values found for the wind’s mass flux are therefore very sensitive to the choice of the parameter , while the amplitude of its variations during the cycle (keeping constant) are only due to the geometrical factor and the (small) variations in . Additionally, it has been suggested that in non-isothermal scenarios the mass flux scales linearly with the mechanical heat flux injected at the surface (Leer & Holzer, 1979; Hammer, 1982; Withbroe, 1988; Hansteen & Leer, 1995). We emphasize that, nevertheless, the results described in this manuscript depend essentially on other parameters of the problem, namely the cyclic variations in the global coronal topology and the consequent open flux-tube expansion factors (cf. the function above). The isothermal approximation should in this case provide a good assessment of the variability of the global wind properties during the activity cycle. Future work will address these issues by taking into account a more complete (i.e, non-isothermal) treatment of the energetics of the wind.
Also, we neglected any effects due to the rotation on the wind’s velocity (which should be small, anyway, as the sun is a slow rotator).
The most remarkable feature found is the large temporal variability of these quantities (Alfvén radius and spin-down time-scale) over an activity cycle. The sun loses mass and momentum non steadily, which translates into variable braking torques applied to the surface. We should add that the latitude and the width of the open-field and wind sources at the surface also vary in time, meaning the torque is applied at different places at different times of the cycle. During most of the activity cycle, the open-field crosses the solar surface at high latitudes only. The exception is the interval about the polarity inversion, when the polar coronal holes are greatly reduced and low latitude open-field flux-tubes appear. We will consider whether this effect can be translated into meaningful upper boundary conditions (forcing, mass and momentum sources/sinks) for the solar convection and dynamo processes in future research.
Now, as the solar cycle progresses the effects of an ever-changing magnetic field topology add up to the variations of the global amplitude of the field. This makes it hard to distinguish between their respective contributions. Let us first clarify the status of our coronal model in regard to this issue. We ran our model through additional test cases in which we varied only the external magnetic field’s amplitude, and not its topology (i.e, the potential magnetic field only varied by a multiplicative factor between consecutive runs). In the first case, we took the first element of our time-series ( in the figures in the preceding sections) and varied the amplitude of the external field by a constant multiplicative factor. The second case corresponds to a purely dipolar external field whose amplitude we also varied, letting us perform a comparison with previous studies. We then computed the mean Alfvén radius and the mass loss rate for each new run and tested against the correlation
suggested by Matt & Pudritz (2008), where is the global field’s amplitude. This scalar quantity is well defined for any given analytical external magnetic field but not for the ones obtained numerically from the dynamo model. We chose here to associate it to the unsigned magnetic flux integrated over the solar surface (after comparing both definitions for the dipolar test cases). Figure 11 (top and middle panels) shows that the expression above fits our test cases for a power law index , which is in good agreement with what Matt & Pudritz (2008) found for their dipolar case. This result validates the (rather large) variations in Alfvén radius we found. In addition, we plotted the same relation for the computations running through the whole cycle (Figure 11, bottom panel). The same power law is still observed for most of the cycle, but it gets slightly steeper just before the polarity inversion and then deviates away from it (flipping sign) immediately after (between and ), drawing a closed cycle in the diagram. This result suggests that the Alfvén radius correlates positively with when the dipolar component of the open magnetic field is strong (which happens during most of the cycle), but may correlate negatively for some other higher order geometries.
As expressed earlier in the text, we computed a succession of steady-state wind solutions rather than the detailed dynamical evolution of a particular event. One justification for this choice is of numerical nature: the disparity of time-scales for the magneto-convective and coronal phenomena is very large. The second and perhaps more important justification is that varying the potential magnetic field in time introduces a non causal perturbation to the system. That would correspond to an instantaneous propagation of a perturbation to the field sources. This violation of causality remains even for a continuous and arbitrarily slow evolution of the background potential field. Note that localised perturbations to a MHD system are propagated away with phase velocities which correspond to the MHD wave modes. These can be arbitrarily small (even null in some places) and very anisotropic. To work around this issue one would need to include the field sources in the domain (and then self-consistently compute the system’s response to their perturbation). Alternatively, one could emulate the proper physical behaviour of such a system by injecting/propagating the hydro-magnetic perturbations through the numerical boundaries. This issue will be the subject of future work.
7. Summary and conclusions
We have performed MHD numerical simulations coupling a solar dynamo model with a corona and solar wind model throughout a complete activity cycle. In short, we found that:
The polarity reversal happens rather abruptly in the corona, in contrast with the progressive evolution of the solar wind’s velocities and of the surface magnetic field (Figure 4).
The zones of application of the braking torque due to the wind vary in time. Overall, the wind’s breaking torque should contribute to slow down the surface layers at high latitudes, but regions of application of torque appear occasionally at lower latitudes.
Future work will focus on testing other types of solar dynamo models, using Babcock-Leighton flux transport based on large scale meridional circulation (Dikpati et al., 2004; Jouve & Brun, 2007) or turbulent magnetic pumping (Guerrero & de Gouveia Dal Pino, 2008), as well as non-isothermal winds. We will furthermore consider the formation of a non-rigid coronal rotation and its effects on the solar wind properties.
http://www.stars2.eu) and the CNRS PNST Interfaces group. Computations were carried out using CNRS IDRIS and CEA’s CCRT facilities. We thank S. Matt for useful discussions.
- Aibéo et al. (2007) Aibéo, A., Ferreira, J. M., & Lima, J. J. G. 2007, Astronomy and Astrophysics, 473, 501
- Antia et al. (2000) Antia, H. M., Chitre, S. M., & Thompson, M. J. 2000, Astronomy and Astrophysics, 360, 335
- Antia et al. (2003) —. 2003, Astronomy and Astrophysics, 399, 329
- Arge & Pizzo (2000) Arge, C. N., & Pizzo, V. J. 2000, Journal of Geophysical Research, 105, 10465
- Aschwanden (2005) Aschwanden, M. J. 2005, Physics of the Solar Corona. An Introduction with Problems and Solutions (2nd edition)
- Belcher & MacGregor (1976) Belcher, J. W., & MacGregor, K. B. 1976, The Astrophysical Journal, 210, 498
- Brown et al. (2011) Brown, Miesch, Browning, Brun, & Toomre. 2011, submitted to ApJ
- Brun et al. (2002) Brun, A. S., Antia, H. M., Chitre, S. M., & Zahn, J. 2002, Astronomy and Astrophysics, 391, 725
- Brun et al. (2004) Brun, A. S., Miesch, M. S., & Toomre, J. 2004, Astrophysical Journal, 614, 1073
- Burnett (1987) Burnett, D. S. 1987, Finite Element Analysis: From Concepts to Applications. (Reading, Massachussetts: Addison-Wesley)
- Cattaneo (1999) Cattaneo, F. 1999, Astrophysical Journal, 515, L39
- Charbonneau (2010) Charbonneau, P. 2010, Living Reviews in Solar Physics, 7, 3
- Charbonneau & MacGregor (1997) Charbonneau, P., & MacGregor, K. B. 1997, The Astrophysical Journal, 486, 502
- Cohen et al. (2009) Cohen, O., Drake, J. J., Kashyap, V. L., & Gombosi, T. I. 2009, Astrophysical Journal, 699, 1501
- Cranmer et al. (2007) Cranmer, S. R., van Ballegooijen, A. A., & Edgar, R. J. 2007, Astrophysical Journal Supplement Series, 171, 520
- DeRosa et al. (2009) DeRosa, M. L., et al. 2009, The Astrophysical Journal, 696, 1780
- Dikpati et al. (2004) Dikpati, M., de Toma, G., Gilman, P. A., Arge, C. N., & White, O. R. 2004, The Astrophysical Journal, 601, 1136
- Endeve & Leer (2001) Endeve, E., & Leer, E. 2001, Solar Physics, 200, 235
- Endeve et al. (2003) Endeve, E., Leer, E., & Holzer, T. E. 2003, Astrophysical Journal, 589, 1040
- Ghizaru et al. (2010) Ghizaru, M., Charbonneau, P., & Smolarkiewicz, P. K. 2010, The Astrophysical Journal, 715, L133
- Gilman (1983) Gilman, P. A. 1983, The Astrophysical Journal Supplement Series, 53, 243
- Gilman et al. (1989) Gilman, P. A., Morrow, C. A., & Deluca, E. E. 1989, The Astrophysical Journal, 338, 528
- Glatzmaier (1985) Glatzmaier, G. A. 1985, Geophysical and Astrophysical Fluid Dynamics, 31, 137
- Grappin et al. (2008) Grappin, R., Aulanier, G., & Pinto, R. 2008, Astronomy and Astrophysics, 490, 353
- Grappin & Léorat (2001) Grappin, R., & Léorat, J. 2001, Astronomy and Astrophysics, 365, 228
- Grappin et al. (2000) Grappin, R., Léorat, J., & Buttighoffer, A. 2000, Astronomy and Astrophysics, 362, 342
- Grappin et al. (1997) Grappin, R., Leorat, J., Cavillier, E., & Prigent, G. 1997, Astronomy and Astrophysics, 317, L31
- Grappin et al. (2011) Grappin, R., Wang, Y., & Pantellini, F. 2011, The Astrophysical Journal, 727, 30
- Guerrero & de Gouveia Dal Pino (2008) Guerrero, G., & de Gouveia Dal Pino, E. M. 2008, Astronomy and Astrophysics, 485, 267
- Guhathakurta et al. (2006) Guhathakurta, M., Sittler, E. C., & Ofman, L. 2006, Journal of Geophysical Research (Space Physics), 111, 11215
- Hammer (1982) Hammer, R. 1982, Astrophysical Journal, 259, 767
- Hansteen & Leer (1995) Hansteen, V. H., & Leer, E. 1995, Journal of Geophysical Research, 100, 21577
- Hu et al. (2008) Hu, Y. Q., Feng, X. S., Wu, S. T., & Song, W. B. 2008, Journal of Geophysical Research (Space Physics), 113, 03106
- Issautier et al. (2008) Issautier, K., Chat, G. L., Meyer-Vernet, N., Moncuquet, M., Hoang, S., MacDowall, R. J., & McComas, D. J. 2008, Geophysical Research Letters, 35, 19101
- Issautier et al. (2004) Issautier, K., Moncuquet, M., & Hoang, S. 2004, Solar Physics, 221, 351
- Jouve & Brun (2007) Jouve, L., & Brun, A. S. 2007, Astronomy and Astrophysics, 474, 239
- Jouve et al. (2008) Jouve, L., et al. 2008, Astronomy & Astrophysics, 483, 949
- Kosovichev (2006) Kosovichev, A. G. 2006, Solar Physics, 238, 1
- Kosovichev & Duvall (2006) Kosovichev, A. G., & Duvall, T. L. 2006, Space Science Reviews, 124, 1
- Krause & Radler (1980) Krause, F., & Radler, K. 1980, Mean-field Magnetohydrodynamics and Dynamo Theory, 1st edn. (Elsevier)
- Leenaarts et al. (2009) Leenaarts, J., Carlsson, M., Hansteen, V., & van der Voort, L. R. 2009, Astrophysical Journal, 694, L128
- Leer & Holzer (1979) Leer, E., & Holzer, T. E. 1979, Solar Physics, 63, 143
- Leer & Holzer (1980) —. 1980, Journal of Geophysical Research, 85, 4681
- Lele (1992) Lele, S. K. 1992, Journal of Computational Physics, 103, 16
- Luhmann et al. (2002) Luhmann, J. G., Li, Y., Arge, C. N., Gazis, P. R., & Ulrich, R. 2002, Journal of Geophysical Research (Space Physics), 107, 1154
- Luhmann et al. (2009) Luhmann, J. G., et al. 2009, Solar Physics, 256, 285
- Martínez-Sykora et al. (2008) Martínez-Sykora, J., Hansteen, V., & Carlsson, M. 2008, Astrophysical Journal, 679, 871
- Matt & Pudritz (2008) Matt, S., & Pudritz, R. E. 2008, The Astrophysical Journal, 678, 1109
- McComas et al. (2003) McComas, D. J., Elliott, H. A., Schwadron, N. A., Gosling, J. T., Skoug, R. M., & Goldstein, B. E. 2003, Geophysical Research Letters, 30, 24
- Miesch (2005) Miesch, M. S. 2005, Living Reviews in Solar Physics, 2, 1
- Moffatt (1978) Moffatt, H. K. 1978, Magnetic Field Generation in Electrically Conducting Fluids, 1st edn. (Cambridge University Press)
- Morel (1997) Morel, P. 1997, Astronomy and Astrophysics Supplement Series, 124, 597
- Parker (1964) Parker, E. N. 1964, The Astrophysical Journal, 139, 72
- Parker (1993) —. 1993, The Astrophysical Journal, 408, 707
- Pinto et al. (2010) Pinto, R., Grappin, R., & Leorat, J. 2010, in Twelfth International Solar Wind Conference, Vol. 1216 (Saint-Malo, (France): AIP), 80–83
- Pinto et al. (2009) Pinto, R., Grappin, R., Wang, Y., & Léorat, J. 2009, Astronomy and Astrophysics, 497, 537
- Poinsot & Lele (1992) Poinsot, T. J., & Lele, S. K. 1992, Journal of Computational Physics, 101, 104
- Rempel et al. (2009) Rempel, M., Schüssler, M., & Knölker, M. 2009, The Astrophysical Journal, 691, 640
- Riley et al. (2006) Riley, P., Linker, J. A., Mikić, Z., Lionello, R., Ledvina, S. A., & Luhmann, J. G. 2006, The Astrophysical Journal, 653, 1510
- Roberts (1972) Roberts, P. H. 1972, Royal Society of London Philosophical Transactions Series A, 272, 663
- Roe & Balsara (1996) Roe, P. L., & Balsara, D. S. 1996, SIAM Journal on Applied Mathematics, 56, 57
- Schrijver & DeRosa (2003) Schrijver, C. J., & DeRosa, M. L. 2003, Solar Physics, 212, 165
- Spalart et al. (1991) Spalart, P. R., Moser, R. D., & Rogers, M. M. 1991, Journal of Computational Physics, 96, 297
- Steenbeck & Krause (1969) Steenbeck, M., & Krause, F. 1969, Astronomische Nachrichten, 291, 49
- Stix (1976) Stix, M. 1976, Astronomy and Astrophysics, 47, 243
- Stix (2002) —. 2002, The Sun: An introduction, 2nd edn. (Springer)
- Thompson (1987) Thompson, K. W. 1987, Journal of Computational Physics, 68, 1
- Tokumaru et al. (2010) Tokumaru, M., Kojima, M., & Fujiki, K. 2010, Journal of Geophysical Research (Space Physics), 115, 04102
- Usmanov (1993) Usmanov, A. V. 1993, Solar Physics, 146, 377
- Vanajakshi et al. (1989) Vanajakshi, T. C., Thompson, K. W., & Black, D. C. 1989, Journal of Computational Physics, 84, 343
- Velli (1994) Velli, M. 1994, Astrophysical Journal, 432, L55
- Vögler et al. (2005) Vögler, A., Shelyag, S., Schüssler, M., Cattaneo, F., Emonet, T., & Linde, T. 2005, Astronomy and Astrophysics, 429, 335
- Wang et al. (2007) Wang, Y., Biersteker, J. B., Sheeley, N. R., Koutchmy, S., Mouette, J., & Druckmüller, M. 2007, Astrophysical Journal, 660, 882
- Wang & Sheeley (1990) Wang, Y., & Sheeley, N. R. 1990, Astrophysical Journal, 355, 726
- Wang & Sheeley (1992) —. 1992, The Astrophysical Journal, 392, 310
- Wang & Sheeley (2006) —. 2006, Astrophysical Journal, 653, 708
- Wang et al. (2003) Wang, Y., Sheeley, N. R., & Jr. 2003, The Astrophysical Journal, 599, 1404
- Weber & Davis (1967) Weber, E. J., & Davis, L. 1967, Astrophysical Journal, 148, 217
- Wedemeyer-Böhm et al. (2009) Wedemeyer-Böhm, S., Lagg, A., & Å. Nordlund. 2009, Space Science Reviews, 144, 317
- Withbroe (1988) Withbroe, G. L. 1988, Astrophysical Journal, 325, 442
- Zurbuchen (2007) Zurbuchen, T. H. 2007, Annual Review of Astronomy and Astrophysics, 45, 297