Sunspot Modeling: From Simplified Models to Radiative MHD Simulations

We review our current understanding of sunspots from the scales of their fine structure to their large scale (global) structure including the processes of their formation and decay. Recently, sunspot models have undergone a dramatic change. In the past, several aspects of sunspot structure have been addressed by static MHD models with parametrized energy transport. Models of sunspot fine structure have been relying heavily on strong assumptions about flow and field geometry (e.g., flux-tubes, “gaps”, convective rolls), which were motivated in part by the observed filamentary structure of penumbrae or the necessity of explaining the substantial energy transport required to maintain the penumbral brightness. However, none of these models could self-consistently explain all aspects of penumbral structure (energy transport, filamentation, Evershed flow). In recent years, 3D radiative MHD simulations have been advanced dramatically to the point at which models of complete sunspots with sufficient resolution to capture sunspot fine structure are feasible. Here, overturning convection is the central element responsible for energy transport, filamentation leading to fine structure, and the driving of strong outflows. On the larger scale these models are also in the progress of addressing the subsurface structure of sunspots as well as sunspot formation. With this shift in modeling capabilities and the recent advances in high resolution observations, the future research will be guided by comparing observation and theory.


Introduction
Magnetic fields on the Sun exist in a large variety of phenomena and interact in various ways with plasma and radiation. In the convection zone large and small scale magnetic fields are generated. These magnetic fields are partially transported into the outer layers of the Sun, i.e., into the chromosphere and the corona. The most prominent example of a magnetic phenomenon is a sunspot as seen in the photosphere. A typical sunspot has a lifetime of a few weeks and has a size of about 30 granules. The magnetic field strength spans from 1000 to 3000 Gauss in the deep photosphere, summing up to a magnetic flux of some 10 22 Mx, typically. For an extensive review of the sunspot structure, we refer the reader to an instructive overview by Solanki (2003).
The magnetic field of a sunspot extends into the interior as well as into the outer layers of the Sun. The most detailed information of sunspots is obtained in the photosphere. The topology of the magnetic field above and beneath the photosphere is poorly understood. In particular our knowledge of the magnetic field extension into the interior presents a theoretical challenge. Direct measurements of the sub-photospheric structure are impossible, but at least for the larger scales, indirect methods are being explored in the framework of local helioseismology (cf. Gizon and Birch, 2005).
Sunspots are central to our understanding of solar magnetism in several aspects. Sunspots are the most prominent manifestation of the large scale cyclic solar magnetic field. Understanding their subsurface structure as well as the processes of formation, dynamic evolution, and decay is crucial for connecting them to the dynamo and flux emergence processes in the solar convection zone (Fan, 2009;Charbonneau, 2010). On smaller scales sunspots provide an ideal environment for studying magnetoconvection for a variety of different field configurations. While quiet Sun and plage regions have been modelled very successfully for almost 3 decades using 3D radiative MHD simulations (see the review by Nordlund et al., 2009), such models were only applied to sunspots in the past five years. The combination of detailed models with the wealth of high resolution observations has substantially advanced our understanding of sunspot structure over the past decade.
In this paper we aim to review our current understanding of sunspots. We approach the problem from different perspectives. In Section 2 we start out by characterizing the global structure of a sunspot. We describe model ideas for the sunspot structure, and describe how sunspots can be treated in a static (non-dynamic) configuration. In Section 3 we discuss the dynamic fine structure of umbra and penumbra. We summarize the key observational facts models have to explain, give an overview about several idealized models in use and summarize the recent progress in radiative MHD modeling of sunspots. In Section 4 the paradigm of sunspot formation by rising and emerging magnetic flux tubes is addressed. We also summarize observations and models of the moat region surrounding sunspots. In Section 5 we give a brief summary of our current understanding of subsurface structure and flow fields surrounding sunspots as derived from helioseismic inversions. We summarize our knowledge on sunspot modeling and conclude on our present understanding in Section 6.

Global Sunspot Structure Revisited
One of the elementary problems for sunspot models is our ignorance on the sub-photospheric morphology. We have no means to directly observe the sub-photospheric structure of sunspots. The photons that reach us originate in the photosphere, which is the transition layer between the convectively unstable interior and the convectively stable atmosphere. Knowledge about the sub-photospheric structure is expected from local helioseismology. However, current inversions techniques do not yet allow to probe the internal sunspot structure with sufficient confidence (Gizon et al., 2009;Moradi et al., 2010). Therefore, attempts to model sunspots rely on theoretical expectations and numerical simulation of magnetohydrodynamics. Any model can only be tested by comparing the (photospheric) surface signatures with observations. The fundamental observational requirement of a sunspot model is that sunspots are manifested on large spatial scales and are dynamically stable.

Sunspot time scales
Sunspots are long-lived (typically weeks) relative to a dynamical time scale of about one hour; one hour corresponds to the travel time of a magneto-acoustic wave across a sunspot in the photosphere. Even though sunspots are globally stable, high resolution observations readily demonstrate that they are not static on small scales. Down to the present-day resolution limit of some 0.1 arcsec, corresponding to 70 km on the solar surface, sunspots have dynamically evolving features that are known as umbral dots, light bridges, bright and dark penumbral filaments, penumbral grains, dark-cored bright filaments, and penumbral twists, just to mention the most prominent structures. Hence, a sunspot is a coherent phenomenon on large spatial and temporal scales. However, it seems crucial to realize that a sunspot is not static, but finds a dynamical equilibrium: A variety of small-scale features evolve on a dynamic time scale to produce a large scale coherent structure on long time scales.
In other words, the dynamic fine structure forms a globally stable sunspot and it is the goal of sunspot physics to understand how an ensemble of short-lived features with small scales is organized to form a coherent large and long-living sunspot. This "fine structure" is complex and is seen in white light images. This fine structure must be associated with magneto-convective motions on the small scales, and it is unquestioned that the small scale motions are the key to understand the energy transport within and the structure of sunspots. But it is a challenge to understand how the small-scale features form a stable sunspot. Hale (1908a) discovered the existence of (vertical) magnetic fields in sunspots. Studying such "vortices", it "occured" to him (Hale, 1908b) that sunspots result from a solar tornado that sucks plasma into higher layers, which obscures the solar surface. Such a tornado would be visible from the earth as a dark spot. The circular flow field would separate the ions and electrons by centrifugal forces such that the circular flow of electrons produces the magnetic field that is perpendicular to the surface, consistent with his observed Zeeman signatures. But only one year later, the British astromoner Evershed (1909) observing in Kodaikanal (Tamil Nadu, India) found that the flow field within the penumbra of a spot is radially outwards with respect to the center of the spot on the solar disk. His result "seems entirely out of harmony with the splendid discovery [...] by Professor Hale", "[...] and it was only after a considerable amount of evidence had accumulated that the preconceived conviction that the motion must be circular was abandoned" (Evershed, 1909).

Suppressed convection
The next theory, which still plays a role in our present understanding, came from Biermann (1941) and Alfvén (1942): In a highly ionized plasma the electric conductivity can be so large that the magnetic fields are frozen-in to the plasma. Biermann realized that the magnetic field in sunspots itself could be the reason for the spot coolness: Outside a sunspot energy is transported to the surface by overturning convection. In a sunspot the convective flow is inhibited by magnetic tension. Hence, a sunspot is dark because it is cooler, and it is cooler because the magnetic field suppresses the heat transport by convection. Hence, the darkness of a spot is due to a decreased surface brightness. Compared to the sunspot surroundings, which has a heat flux of 6.31 Ö 10 7 W/m 2 , the average penumbra heat flux is reduced by 25% and the umbra heat flux by 77%, respectively (e.g., Jahn and Schmidt, 1994). These values for the heat flux correspond to effective temperatures of 5777 K 1 , 5275 K, and 4000 K of quiet Sun, penumbra, and umbra, respectively.

Heat flux dilution by funneling
Another concept, discussed by Jahn (1992) and relying on an idea from Hoyle (1949), also plays an important role in our understanding: plasma = 8 2 of a monolithic bundle of magnetic field lines increases with depth beneath the photosphere. Hence, to balance the hydrostatic stratifications in and outside the spot and to conserve the magnetic flux, must increase with depth and the spot diameter must decrease. With this geometrical funnel and the assumption that the umbra is thermally isolated from the penumbra and the surroundings, Hoyle constructed the following explanation for the umbral darkness: At a certain depth a given heat flux enters the umbra from below. This entering heat flux is the same as in the sunspot surroundings. Then, even if all that energy is transported (by magneto-convection) to the surface it will dilute because the umbral area at the surface is larger. And as we will discuss further down, this effect is present in the tripartite models of Jahn and Schmidt (1994). Deinzer (1965) 2 pointed out that convection cannot be suppressed completely since sunspots umbrae are still as hot as 4000 K, and radiative heating or heat conduction cannot supply the necessary heat flux. From this theoretical argument it must be concluded that energy transport by convection must exist in sunspots and, in particular, in sunspots umbrae. Meyer et al. (1974) studied the possible modes of magneto-convection and found that for the first 2000 km beneath the surface, the magnetic diffusivity, , is smaller than thermal diffusivity, . In this surface region, convection sets in as overstable oscillations (corresponding to standing Alfvén waves). The latter was also proposed earlier by Savage (1969) and has been found in idealized magneto-convection models, as we describe in Section 3.5.

Modified convection and funneling
In deeper layers, down to 20 000 km overturning convection takes place. Hence, heat transport by (magneto-) convection is essential in sunspots, meaning that magnetic fields cannot inhibit convection. It is true, however, that the magnetic field modifies the mode of convection. In Sections 3.5 and 3.6.2, we will describe simulations of magneto-convection in strong magnetic fields, but it should be noted that until today there is no parametrized theory of heat transport in magnetoconvection as there is the mixing length theory (or more sophisticated moment approaches) for convection without magnetic field.
Realizing that global sunspot models must rely on energy transport by convection, a mixing length theory with a depth dependent mixing length parameter was developed (Deinzer, 1965;Jahn, 1989). A model that could be compared quantitatively with observed sunspots was presented by Jahn and Schmidt (1994). Most notably it turned out that these models rely on Hoyle's concept: The 'tripartite' model has three stratifications: the thermally isolated umbra, the penumbra, and the quiet Sun. In a typical model, the umbra diameter doubles between 15 Mm depth and the surface, i.e., the area increases by a factor of 4 and the heat flux is reduced by the same amount. In some of these models, the regular solar heat flux enters the umbra at a depth of 15 Mm leading to the same umbral surface brightness as it is observed, i.e., Hoyle's concept is built-in in the tripartite models. Note that the depth dependence of the mixing length parameter decides about the temperature gradient of the umbral depth stratification.

Why there is no bright ring around a funnel-shaped spot
For a long time another expectation of a funnel-shaped spot was discussed: It was thought that the funnel blocks heat, which accumulates beneath the spot such that a bright ring around the spot should be produced. Several observational investigations have found some evidence for the existence of such bright rings: Fowler et al. (1983) found bright rings in the 0.1 -0.3% range, the more recent study by Rast et al. (2001) brightness enhancements of 0.5 -1%. A major challenge in these observational studies is the proper separation from the dominant effect of facular brightening, which is a near surface effect unrelated to the presence of a sunspot. Whether there is observational evidence for a bright ring independent from facular brigthening remains an open issue. It is, however, clear from the present observational constraints that the amplitude of such a bright ring would be insufficient to account for the heat flux blocked by the sunspot. Spruit (1977) resolved this puzzle: The thermal conductivity of the convection zone is so large that thermal disturbances are smoothed out so effectively that the remaining temperature fluctuations are too small to be observed. In a later review paper, Spruit (1992) made this effect "more understandable with a kitchen analogy": If you put a small piece of a thermally insulating material on top of a electrically heated copper plate, there will not be an enhanced temperature in the vicinity of the insulator, since the blocked heat will be distributed very effectively across the copper plate. The insulator will be fairly cool at the top, but the copper plate, which has a high heat conductivity, will essentially have the same high temperature everywhere. This idea has been applied subsequently in more detailed sunspot models by Spruit (1982a), Spruit (1982b, Chiang and Foukal (1985) and was recently reviewed by Spruit (2000). Also MHD models with radiative transfer do not show evidence for bright rings (Rempel, 2011c).

Subsurface morphology
In the preceding Section 2.2, we implicitly assumed that a sunspot is a funnel-shaped object consisting of magnetized plasma. This class of models is termed monolithic. One can also envisage that bundles of field lines separate beneath the photosphere such that gaps/columns of field-free plasma exist just below the surface. The latter type has been referred to as the jelly fish or cluster or spaghetti models (Parker, 1979b,c;Spruit, 1981b;Choudhuri, 1992). The subsurface morphology relates also to the question of what determines the apparent stability of sunspots over long time scales despite more rapid evolution on smaller scales. We will come back to this aspect in Section 4.4.
Spaghetti model and field-free gaps: The field-free columns close at the photosphere, but broaden with depth, since plasma = 8 / 2 is about unity in the photosphere and increases with depth: At the surface the magnetic pressure is strong enough to squeeze the field-free gaps such that Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 only small umbral dots remain at the surface as the imprint of the field-free convective columns. An umbral dot would correspond to the peak of a field-free gap. More recently, Spruit and Scharmer (2006) suggested that such field-free gaps in the inclined magnetic field of the penumbra may result in elongated bright filaments, instead of in point-like dots, thereby proposing an explanation for the brightness of the penumbra. The surplus brightness of the penumbra relative to the umbra would then be due to the fact that the convective cell can become larger in the more inclined and weaker magnetic field as in the less inclined (more vertical) and stronger field of the umbra.
Are spots monolithic or spaghetti-like? At present, we do not have the means to give a rigorous answer. As discussed in Jahn (1992), fine structure phenomena like umbral dots can be explained in both frameworks. The same is true for the umbral heat flux: it is either done by fieldfree convection in gaps or by a magneto-convective process. The spaghetti model is compelling, since it explains naturally the existence of umbral dots and bright filaments, but it is unknown (a) what forces would keep the magnetic field lines together in the photosphere so that they can form stable long-lived sunspots, and (b) how it happens that the spot has a magnetic field strength that decreases systematically from spot center. In the spaghetti model, magnetic and non-magnetic convecting plasma are distinct from each other, while a monolithic model forms an entity in which magneto-convection takes place. Theoretical considerations like the one mentioned in Section 2.2.3 (Meyer et al., 1974) predict the existence of overturning magneto-convective processes for depths larger than 2000 km in which magnetized and non-magnetized plasma gets mixed, and hence, favor the monolithic approach. These monolithic models have the advantage that they can be investigated in a quantitative manner, at least if one assumes that the magnetic field can be described by a mean magnetic field, which finds itself in a global magnetostatic equilibrium.
Obviously, the sunspot brightness and its subsurface structure are intimately connected. In summary, it can be said that the coolness of sunspots relative to the surrounding quiet Sun is explained by the tension of the magnetic field, which tends to suppress convective motions and/or the channeling of heat in an expanding fan, thereby diluting the heat flux. It is more difficult to understand why sunspots are as hot as they are. Since neither radiative transport nor heat conduction can account for the surface brightness of sunspots, the energy must be transported by (magneto-)convective flows. Indeed, the fine structure manifests the inhomogeneities of the magnetic and velocity field and testifies that the energy transport in sunspots happens on small spatial scales by the motion of plasma. The essential question is how magnetized these convective flows are. Observations point toward reduced magnetic field strength in umbral dots (see Section 3.1), which agrees also with recent results from 3D radiative MHD simulations (see Section 3.6.2). Overall, magneto-convective models point toward the presence of overturning convection in both umbra and penumbra that takes place in regions with reduced, albeit not zero field strength. Most importantly, it is found that this mode of convective energy transport can originate in initially monolithic magnetic field, i.e., the presence of sunspot fine structure cannot be taken as support for the "spaghetti-like" subsurface structure.

Stability of monolithic models
Sunspot models have recently been reviewed in detail by Moradi et al. (2010), but see also Spruit (1981a) and Jahn (1997). Here we consider the issue of sunspot stability in monolithic models: Sunspots are stable relative to the dynamical time, i.e., Alfvén waves are estimated to travel across a spot in about 1 hour, while the lifetime is in order of weeks. How can it be that all this dynamic fine structure constitutes a spot that is globally stable?
The question of stability can be addressed if one assumes a monolithic vertical magnetohydrostatic magnetic flux tube that fans out with height, being bounded by current sheets (Pizzo, 1986;Jahn, 1989;Pizzo, 1990;Jahn and Schmidt, 1994). In such models it is implicitly assumed Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 that the heat transport is attributed to magneto-convection. Yet, the heat transport is not described dynamically, but parametrized in the context of mixing length theory (e.g., Hansen and Kawaler, 1994;Stix, 2004) with a reduced mixing length parameter (Deinzer, 1965;Jahn, 1989). Thereby the dynamic fine structure is ignored and only their averaged effect on the stratification for umbra and penumbra is accounted for. As a consequence, the magnetic field in the models does not really consist out of simple and straight field lines, but must be represented by a mean field, which result from averaging the small scales. In monolithic models it is then implicitly assumed that this mean field follows the rules of magnetostatics.  Jahn and Schmidt (1994) consists of three stratifications as sketched in the left panel: An umbra, a thick penumbra, and the quiet Sun. The umbra is thermally isolated and confined by a current sheet called peripatopause, while the penumbra is confined by a current sheet called magnetopause. For depths smaller than bp the magnetopause transmits 70% of the quiet Sun heat flux, thereby accounting for the surplus brightness of the penumbra when compared to the umbra. The right panel shows a surface plot of the (logarithmic) magnetic pressure. Credit: Jahn and Schmidt (1994), reproduced with permission of the ESO.
Exemplary for the class of monolithic models, we consider the tripartite sunspot model presented by Jahn and Schmidt (1994) (but see also Pizzo, 1990, for a nicely constructed model). It is configured to be in magneto-static equilibrium with a total pressure balance horizontally and a hydrostatic equilibrium vertically. The three stratification are separated by two current sheets between the umbra, and penumbra, and between the penumbra and the quiet Sun horizontally, as sketched in the left panel of Figure 1. The gas pressure jumps are balanced by magnetic pressure, which is shown logarithmically as a surface plot in the right panel of Figure 1. The tripartite model describes a sunspot down to 15 Mm beneath the surface. This configuration can be stable against the interchange instability (Meyer et al., 1977;Schüssler, 1984;Buente et al., 1993) in the first 5 Mm or so beneath the photosphere (Jahn, 1997). In these upper layers of the convection zone the inclination of the magnetopause, i.e., the interface between spot and surrounding, is so large that buoyancy forces make the spot to float on the granulation. In deeper layers, beyond 5 Mm, the inclination of the outermost magnetic field line, i.e., the magnetopause, is small relative to the vertical. There, interchange (fluting) instability is no longer suppressed by buoyancy effects, and the magnetic configuration of a monolithic sunspot is unstable. In these depths one would expect that strands of field lines separate to form a spaghetti configuration. However, in this depth range is larger than , and magneto-convection is active, such that those spaghetti are far from being Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 al dente: The magnetic and non-magnetic plasma is expected to mix. Indeed, it has been proposed that the magnetic field strength progressively weakens in these deep layers shortly after the formation of a sunspot. The decreasing field strength, the convective motions, and the interchange instability dynamically disrupt the sunspot magnetic field from the deeper roots (Schüssler and Rempel, 2005). Hence, the magnetic field in the deeper layers may be dispersed, but the floating part of the sunspot is stable.
As mentioned before, it is essential to realize that monolithic models cannot be static. To transport sufficient energy, magneto-convection must also be present in the upper layers where it is stable against fluting. This seems feasible and we see no necessity to discard monolithic models and to prefer jelly fish models, instead. Jelly fish models may be compelling to explain umbral dots and bright filaments, but we are faced with the following problems: (1) How is the separation of plasma into magnetic and non-magnetic components maintained? How could that be achieved in a turbulent and convectively unstable stratification? At least in layers deeper than 2 Mm, magneto-convection is expected to mix the plasma. (2) What makes the sunspot to be stable and to behave like a coherent structure in which the field strength and inclination monotonically decreases outwards from spot center?

Sunspot Fine Structure
Reviews on sunspot fine structure have been nicely presented recently by Tritschler (2009), Borrero (2009.

Observations of umbral dots
Sunspots umbrae show about 20% 3 of the quiet Sun luminosity, a value which cannot entirely be explained through radiative energy transport. Indeed, observations with a spatial resolution better than 0.5 arcsec reveal that the umbra harbors dynamic inhomogeneities, which have been called umbral dots by Danielson (1964). They are observed as bright dot-like spots with typical sizes of half an arcsec or less, embedded in a more uniform and darker background. These umbral dots seem to be present in all sunspots, although their intensity varies significantly. In some spots they can be almost as bright as bright penumbral filaments, in other spots their intensity is much smaller. In the latter case, the dot-like intensity variations occurs in a background that shows a lower intensity.
Umbral dots are an obvious signature of convection, yet it is not so obvious to understand the type of convection that leads to umbral dots. In the field-free gap idea of Parker, the convection is confined by the strong surrounding magnetic field, such that the column of convection narrows upwards and only a small brightening is seen at the surface. It is established observationally that the magnetic field in umbral dots is weaker than in the surroundings and that an upflow of at least a few hundred m s -1 is associated with them (Socas-Navarro et al., 2004;Rimmele, 2004Rimmele, , 2008Bharti et al., 2007). The latter two observations also establish the presence of dark lanes across umbral dots. From observations it is still unclear how small the field strength in umbral dots really is, since it is difficult to resolve the field strength gradients along the LOS and because there are stray light issues. Numerical models presented in Section 3.6 point toward a dramatic reduction to almost zero field strength in sub-photospheric layers, however, the reduction in line forming regions is less pronounced.

Morphological description
The penumbra is a manifestation of small-scale structure. The variety of the intensity fine structure has recently been reviewed by Solanki (2003), Bellot Rubio (  Intensity pattern: In essence, there are bright and dark filaments, as well as penumbral grains. These features are barely visible at a resolution of about one arcsec. At smaller scales, these features exhibit structure on spatial scales of 0.1 arcsec in the inner penumbra (see, e.g., Scharmer et al., 2002) and of 0.35 arcsec in the mid and outer penumbra (Sütterlin, 2001): Bright filaments have dark cores (Scharmer et al., 2002;Sütterlin et al., 2004) and display inclined dark stripes 4 along filaments (Ichimoto et al., 2007b). Spruit et al. (2010) termed these stripes striations and interpreted them as corrugations, which form as a consequence of the fluting instability at the interface between the ambient magnetic field and the filament. In their explanation the overturning convection inside a field-free gap drags the inclined field lines downward such that the stripe migrates outwards as it is observed. Of course another scenario based on flow channels is also possible: The fluting instability could produce such migrating stripes at the interface between a magnetized flow channel and the ambient magnetic field, which wraps around the flow channel. In this scenario, the outward (Evershed) flow could drag the striations to migrate along the flow channel, and since the ambient magnetic field is inclined, the stripes appear inclined relative to the flow direction. The challenge to understand the intensity fine structure consists in measuring their spectroscopic and spectropolarimetric signatures with the goal to derive their thermodynamic properties as well as their velocity and magnetic field on scales as small as possible. Only recently, with the technological advance of adaptive optics and with observations from space, it has become possible to acquire data with a spatial resolution of 0.3 arcsec for exposure times as long as 5 sec or more. This is a necessity to collect enough photons to have high spatial, spectral, and polarimetric resolution.
At a spatial resolution of better than half an arcsec, it can be demonstrated that not only the intensity and velocity, but also the magnetic field consists of a filamentary structure (Title et al., 1993;Langhans et al., 2005;Tritschler et al., 2007;Ichimoto et al., 2007aIchimoto et al., , 2008a. Actually, at a spatial resolution of better than half an arcsec, all physical quantities in the penumbra show small-scale variations and predominantly filamentary (radially elongated) features. The latter is demonstrated in Figure 2 in which maps of intensity, velocity, and circular polarization are displayed at the SP/Hinode resolution of about 0.3 arcsec. Therefore, it is obvious that there is an intimate interaction between the convective flows and the magnetic fields.
However, at a spatial resolution worse than 1 arcsec, the penumbra looks fairly uniform and is on average brighter than the umbra, but less bright than in the surrounding granulation. But even if the penumbra is less bright on average, the small scale peak-to-peak intensity variation in the penumbra is larger and the spatial scales of the variations are smaller than in the granulation. The same is true for velocities in the penumbra. Line-of-sight (LOS) velocities in the penumbra of more than 5 km s -1 have been derived from Doppler shifts of photospheric lines (e.g., Wiehr, 1995;Franz and Schlichenmaier, 2009) and radial flow channels with widths of less than half an arcsec are observed (e.g., Tritschler et al., 2004;Rimmele and Marino, 2006).

Evershed flow, uncombed magnetic field, and net circular polarization (NCP)
To understand the nature of the penumbral fine structure, it is essential to know the topology of the velocity field and the magnetic field. The first attempt to measure the flow field was undertaken by Evershed in 1908(see Evershed, 1909 in order to test Hale's tornado theory of sunspots. Yet, instead of a circular flow, Evershed found a radial outflow of plasma. Until today we lack a fully consistent theory for sunspots, although substantial progress in modeling the characteristic features of the penumbra has been made in recent years. We summarize first important observational aspects, before discussing models for sunspot fine structure. Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Overview: Before we review this topic, we think it is instructive to describe Figure 3, adhering to Franz (2011): The figure shows the velocity maps (upper row, determines from line wing bisector shift in Fe i 630.15 nm) and 'Ichimotograms' (bottom row) of a penumbra close to disk center (left column, spot 1) and a penumbra that is seen at a heliocentric angle of 47°(right column, spot 2): Spot 1: At disk center, horizontal velocity components do not contribute to the LOS velocity, such that the vertical velocity component is prominent. The blue elongated patches, which dominate the inner and mid penumbra, are associated with upflows of up to 1.5 km s -1 . The yellow and red downflow patches are dominant in the very outer penumbra, but also exist in the mid and a few even in the inner penumbra. In the outer penumbra the downflow velocities exceed 5 km s -1 , but are clipped at 1.5 km s -1 in the figure so that smaller velocities are better visible. The Ichimotogram is an image of Stokes-V in the red line wing at 0 + 211 nm (such images were discussed by Ichimoto et al., 2007a) 5 . In white patches the V-signal has opposite sign as in the overall spot and, hence, trace locations where the polarity is opposite to the spot polarity. These opposite polarity patches are co-spatial with downflow patches. Spot 2: At large heliocentric angles the horizontal velocity component becomes visible. Since the horizontal flow speeds are larger than the vertical velocities (see, e.g., Schlichenmaier and Schmidt, 2000), the velocity map is dominated by the horizontal (Evershed) flow, producing a blueshift on the center-side penumbra, and a redshift on the limb-side penumbra. The flow velocity are clipped at ± 3.5 km s -1 . Velocities exceeding 4 km s -1 are common on both penumbral sides.
The flow field: With high spatial resolution, it is now established that the flow has a filamentary structure Rimmele and Marino, 2006). On azimuthal average, the flow is essentially horizontal with a small upward component in the inner and a small downward component in the outer penumbra Schmidt and Schlichenmaier, 2000;Tritschler et al., 2004;Langhans et al., 2005). Recent observations have revealed that radially aligned up-and downflows exist on small scales next to each other (Sainz Dalda and Bellot Rubio, 2008;Franz and Schlichenmaier, 2010). Regarding the height dependence of the flow, St John (1913) (see also Ichimoto, 1987Ichimoto, , 1988 found that the flow velocity decreases with the formation height of the absorption line. In chromospheric lines the flow reverses its sign, which is being referred to as the inverse Evershed flow. Studying the line asymmetries of photospheric lines, one finds convincing evidence that the flow is predominantly present in the very deep photosphere, i.e., beneath = 0.1 (Maltby, 1964;Schlichenmaier et al., 2004;. The peaks of the flow velocities measured in the penumbra are substantially larger than what is measured in the granulation. Individual penumbral profiles exhibit line satellites that are Doppler shifted by up to 8 km s -1 (e.g., Wiehr, 1995). From inversions, velocities well above 10 km s -1 have been found by del Toro Iniesta et al. The magnetic field: Attempts to describe the magnetic field as being uniform along the lineof-sight are clearly inconsistent with the measured Stokes ( ), ( ), and ( ) profiles (e.g., Westendorp Plaza et al., 2001a,b). In particular, the penumbral V-profiles with three or more lobes and non-vanishing NCP-values (see below) cannot be explained by one component Schlichenmaier and Collados (2002): Gradients and/or discontinuities along the line-of-sight must be present. Therefore, it was proposed that the magnetic field is interlocked or in other words uncombed (Solanki and Montavon, 1993). In order to keep things as simple as possible, the magnetic field is assumed to have two components with different inclinations. Indeed, if the observed Stokes profiles with a spatial resolution of about 1 arcsec are interpreted with two components by means of inversions techniques, the fits to the observations are substantially better and reproduce essential features of the line asymmetries, which is not possible with only one component Borrero et al., 2004Borrero et al., , 2005Beck, 2008). Such inversions yield one less inclined magnetic component that is only slightly Doppler shifted, and a second component with somewhat weaker and more inclined, i.e., approximately horizontal field. This second component carries the Evershed flow, with spatially averaged flow speeds of about 6.5 km s -1 .
These inversions also show that the magnetic field of the second component is aligned with the associated flow, pointing slightly upwards in the inner and slightly downward in the outer penumbra. The inclination angle relative to the local vertical of the first magnetic field component increases from some 30°at the umbral-penumbral boundary to some 60°at the outer penumbral boundary. Inversions, which are optimized to locate the height of the flow layer, find that the flow is present in the very deep atmosphere, in the continuum forming layers Borrero et al., 2006;Jurčák and Sobotka, 2007;Jurčák and Bellot Rubio, 2008).
At 0.3 arcsec spatial resolution, spectropolarimetric measurements reveal that, at least in the inner penumbra, the more inclined magnetic component, which carries the flow, is associated with the dark cored bright filaments. Individual dark cores have a smaller degree of circular polarization than their lateral brightenings (Langhans et al., 2007). A thorough analysis shows that the latter statement is also true for the total polarization and that the dark core magnetic field is weaker and more inclined than in the lateral brightenings . Additionally, these studies confirm that the dark cores harbor strong Evershed flows .
The magnetic canopy: Outside the white-light boundary of the penumbra, the inclined magnetic field continues into the chromosphere. It forms a magnetic canopy in the surroundings of the sunspot, rising with increasing distance up to a height of approximately 800 km (Solanki et al., 1992). In the canopy a radial outflow is present, which is interpreted as the continuation of the Evershed flow (Solanki et al., 1992;Rezaei et al., 2006). However, it is estimated that only a few tenth of the mass flows into the canopy. The rest of the penumbral Evershed flow must disappear within the penumbral downflow regions.
The net circular polarization (NCP): The NCP, ∫︀ ( ) d , is a quantity that intimately links the flow and the magnetic field: NCP can only be non-zero, if and only if velocity gradients along the line-of-sight are present (e.g., Sánchez Almeida and Lites, 1992). The magnitude and the size of the NCP depends on the gradient of the line-of-sight velocity, but also on the gradients in the magnetic field strength, inclination, and azimuth (Landolfi and Landi Degl'Innocenti, 1996;Müller et al., 2002Müller et al., , 2006. A predominantly horizontal flow channel embedded in a less inclined background magnetic field successfully explains symmetry properties of NCP maps of sunspots  as well as the center to limb variation of NCP (Martínez Pillet, 2000;Borrero et al., 2007).
Yet, some recent interpretations of NCP maps require that the flow component should be associated with stronger magnetic field (Tritschler et al., 2007;Ichimoto et al., 2008b), rather than being associated with the same or weaker magnetic field in the flow channels, as we would expect from the models. Since there are also other indications for these stronger magnetic fields (e.g., Bellot Cabrera Solana et al., 2008;, the concept of embedded flow channels will need to be reviewed taking into account these new measurements. The magnetoconvective models described in Section 3.6.4 lead to flow channels with enhanced horizontal magnetic field strength consistent with these recent observational findings. Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Magnetized or non-magnetized flow: In terms of modeling the Evershed flow, it is crucial to know whether or not the flow is magnetized. While NCP can be generated by a field-free flow in a magnetized environment (e.g., Steiner, 2000), the observed profiles in certain locations in the penumbra show more than two lobes (e.g., Schlichenmaier and Collados, 2002;Ichimoto et al., 2007a;Beck, 2008Beck, , 2011Franz and Schlichenmaier, 2010).
Two example profiles of Stokes-V, which have more than two lobes, are displayed in Figure 4. The profiles come from a spot close to disk center ( = 3 ∘ ). The left profile shows a bump in the blue wing (of the two iron lines at 630.15 nm and 630.25 nm). This bump is attributed to an magnetic upflow that is superposed with a contribution from the ambient magnetic field, which is roughly at rest. Both components have the same magnetic polarity. It can be shown that a non-magnetic upflow does not produce such a bump in Stokes-V, even though it makes a bump in Stokes-I. The right panel profile is crucial to understand the penumbra: It shows an extra bump of opposite sign in the red wing of Stokes-V. This bump must be generated by a redshifted component of opposite polarity! Hence, this profile can only be reproduced if (at least) two magnetic magnetic components of opposite polarity contribute to the profile (e.g., Franz, 2011) 6 . Note that a 'typical' magnetogram, which measures Stokes-V only at two wavelength positions, would show the spot polarity, and would miss the opposite polarity part of the profile. For a spot at disk center, Hinode/SP profiles show clear signatures of opposite polarity for about 40% of all downflow patches, i.e., for about 17% of all penumbral pixels. This is a lower limit, since these signatures are blended by noise.
Another clear evidence for the magnetic nature of the Evershed flow comes from the inversion results based on two components, which we described above. They also demonstrate that the Doppler shifted second component is magnetized. This result is inferred from spots off disk center, typically at ≈ 30 ∘ .

Penumbral modeling
The previous section stresses the point that the penumbra is a phenomenon of complex interaction of magneto-convective forces and radiation in a regime of inclined magnetic field of intermediate strength.
One simplified view on this problem is to consider a separation between convective plumes and a magnetic configuration as it is done in the field-free gap model. Another simplified view is by dealing with the problem in ideal MHD, in which the thin flux tube approximation is applicable. The latter perspective is taken in the siphon flow model and the dynamic extension of it, the moving tube model. Yet, for a full understanding it seems necessary to take into account dissipative magneto-convection driven by radiation. But we want to stress that simplified models often help to isolate the dominating physical processes, and to understand the essentials.

Convective models: field-free gap and convective rolls
Originally proposed by Parker (1979b) to explain the umbral dots, Spruit and Scharmer (2006) and Scharmer and Spruit (2006) extended the concept of the field-free gaps to explain the bright penumbral filaments and they realized that such a configuration may also produce the dark cores within bright filaments, caused by a subtle radiative effect at the top of the field-free gap. The idea of field-free gaps in the penumbra is that the inclined penumbral magnetic field produces bright elongations instead of dots. The gaps are supposed to be void of magnetic field and to be connected to the surrounding quiet Sun. Within the gaps, overturning convection transports ample amounts of heat, which would account for the brightness of the penumbra. The convective flow field is directed upwards along the central lane of the filament and downward at the edges of the long sides of the filaments. Within the field-free gap there may exist a radial outflow that corresponds to the Evershed flow. The problem with this description is that the Evershed flow, which is observed to be magnetized, needs to be non-magnetized in the field-free gap. This is, however, not an intrinsic limitation of the model since the field-free gap model could be generalized to accommodate a horizontal field.
The field-free gap model shares some similarity with the model of the convective rolls proposed by Danielson (1961) (see also Grosser, 1991, for a numerical investigation on this model). Here convective rolls lie radially aligned next to each other, two such rolls would form one filament as they rotate in opposite direction, producing an upflow in the central lane and a downflow at the lateral lanes. While the resulting flow near the = 1 level is very similar to the flow assumed in the gap model and essentially indistinguishable in observations, there is nevertheless a profound difference. The upflow in the gap model is very deep reaching and connects, therefore, to layers with large heat capacity, while shallow convection rolls alone cannot provide the brightness of the penumbra over long time scales.
Danielson assumed that a horizontal magnetic field component would be associated with the rolls. This model has been discarded for two reasons: (1) There was no evidence for the corresponding convective flow field, and (2) a major fraction of the magnetic flux in the penumbra is directed upwards, and not horizontal. However, reason (1) depends on spatial resolution and the issue is not settled yet, as we cannot rule out the existence small amplitude vertical motions of a few hundred m s -1 . Reason (2) could be overcome by assuming that the rolls are separated by less inclined (more vertical) magnetic field lines, which constitute a more or less static background magnetic field. And, magnetized rolls interlaced by a static background field that is less inclined relative to the vertical would also meet the observational requirements of two magnetic components in the penumbra. In this respect, at least in principle, it is possible that the horizontal magnetic component carries an Evershed flow.
The problem here is that -up to now -there is only little support for downflows along the Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 edges of bright filaments. Rimmele (2008) finds a convective-roll like flow field in a filament that extends into the umbra for a sunspot close to disk center. However, he uses filtergrams with only two spectral positions at modest spectral resolution. Bellot Rubio et al. (2010), who acquired spectroscopic data at high spatial and spectral resolution, did not find indications for up-and downflows associated with a dark-cored bright filament at disk center. Still, the crucial question of vertical flows in the penumbra is not settled. In order to minimize the effects of the horizontal radial outflow and of possible flows in azimuthal direction, sunspot observations at disk center are needed to learn about the presumably small vertical flow component. See also the discussion in Section 3.7 for further detail.

Flux tube models: stationary and dynamic
Another class of simplified models is based on the assumption that penumbral filaments can be identified with primarily horizontal flux tubes near the = 1 surface. The properties of flux tubes have been studied in stationary setups as well as dynamic configurations.
In the self-consistent magneto-static tripartite sunspot model of Jahn and Schmidt (1994) the surplus brightness of the penumbra relative to the umbra is produced by a heat transfer through the magnetopause, i.e., through the interface between the quiet Sun and the penumbra. This additional heat is thought to be distributed horizontally by interchange convection of magnetic flux tubes. The idea of dynamic magnetic flux tubes is compatible with the observationally finding of multiple magnetic components in the penumbra.
This motivated the study of the dynamics of a single thin magnetic flux tube as it evolves in a 2D static model background Schlichenmaier et al. (1998a,b). However, these studies did not confirm the concept of interchanging magnetic flux tubes that distribute heat horizontally. Instead, these studies created a new picture: The simulated tube lies along the magnetopause of the tripartite sunspot model and is taken to be a bundle of magnetic field lines with penumbral properties. Initially the tube is in magneto-static equilibrium. However, at a magnetopause that is sufficiently inclined, radiative heat exchange between the tube and the hotter quiet Sun triggers an instability: A thin magnetic flux tube that initially lies along the magnetopause, (a) feels the hotter quiet Sun, (b) heats up by radiation most effectively just beneath the photosphere, (c) expands, (d) rises through the sub-photospheric convectively unstable stratification, and (e) develops an upflow along the tube, which brings hot sub-photospheric plasma into the photosphere. (f) This hot upflow cools radiatively in the photosphere and streams radially outwards with supercritical velocity. The radiative cooling sustains the gas pressure gradient that drives the flow. (g) The outflow intrudes the convectively stable photosphere up to a height of some 50 to 100 km. The equilibrium height is determined by the balance of the diamagnetic force, which pulls the conducting tube upwards toward decreasing magnetic field strength and the downward acting buoyancy, which increases as the tube is being pulled up in a convectively stable stratification.
Weak magnetic field at footpoint: The gas pressure gradient that drives the flow is caused by a surplus gas pressure building up inside the part of the tube that rises through the subphotospheric stratification. At the footpoint, i.e., the intersection of the tube with the transition layer from convectively unstable to stable, the gas pressure is high, and in order to balance the total pressure with the surroundings, the magnetic field strength is strongly decreased relative to the surroundings. In this sense the upflow footpoints can be considered as regions of weak magnetic field strength. In other words, the moving tube model is a magneto-convective mode, which consists of a region of weak-field plasma that harbors hot upflows and that travels inwards.
In principle, the effect leading to the up-and outflow works like an inverse convective collapse: In the classical convective collapse the plasma in the tube is cooled and a downflow occurs. Here, the heating of the plasma results in an upflow, and consequently the magnetic field strength in the Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 tube decreases as the flow continues. In the photosphere, the gas pressure gradient is sustained by radiative cooling.
The moving tube scenario successfully explains a number of observational findings: (i) Penumbral grains are the photospheric footpoints of the tube, where the hot and bright plasma enters the photosphere. (ii) The upflow turns horizontally outwards in the photosphere and cools radiatively until it reaches temperature equilibrium. This determines the length of the penumbral grains (Schlichenmaier et al., 1999). (iii) The footpoints migrate inward, as many observed penumbral grains do (e.g., Sobotka and Sütterlin, 2001). (iv) The horizontal outflow corresponds to the Evershed flow. (v) The tube constitutes a flow channel being embedded in a background magnetic field. This is in agreement with the uncombed penumbra, and produces realistic maps of NCP.
Magneto-convective overshoot: An interesting effect that can be studied with the idealized moving tube model, is related to overshooting (Schlichenmaier, 2002. The upflow shoots into the convectively stable photosphere, and is turned horizontally by the magnetic curvature forces along the tube. The dominating forces here are the centrifugal force of the flow, 2 and the magnetic curvature force, 2 /(4 ), with being the curvature. In equilibrium equals A , with A being the Alfvén velocity. During the evolution of the tube the velocity is roughly constant, but the magnetic field strength and, hence, A decreases, leading to an overshoot into the photosphere that is convectively stable. As the plasma overshoots the density increases and buoyancy forces pull the flow back down. This results in an oscillation of the outflow around its equilibrium position such that the tube adopts a wave-like shape, i.e., the plasma first shoots up and then down, again passing the equilibrium position. Such a wave can be considered quasi-stationary, and the crest of such a wave can be compared with the properties of a siphon flow (see below). Hence, the flow yields a serpentine shape, looking like a sea serpent, and evidence for such radially aligned up and downflows has been presented by Sainz Dalda and Bellot Rubio (2008). The amplitude of this wave increases as the magnetic field strength decreases, and eventually the downflow part dives in the sub-photosphere. There the stratification is convectively unstable and the magnetic flux tube experiences a dynamic evolution, that produces outward propagating waves. This scenario produces downflows and makes the tube to disappear within the penumbra. Thereby it would solve a problem of the moving tube model: the out-flow would not extend into the surrounding canopy, but would disappear within the penumbra, as it is observed.
Serpentine flow: Such a two-dimensional serpentine solution was criticized to be unstable in three dimensions (Thomas, 2005), arguing that buoyancy forces make the wavy tube to fall over sideways. But this argument is not valid, since the influence of the upflow at the footpoint of the tube is not taken into account. At the footpoint the plasma is ejected upwards into the photosphere and due to conservation of momentum, the plasma overshoots and follows an up-and down wavy behavior. The fact that the density at the upper crest is larger than in the surroundings does not make the tube to fall over. As an analogy, one may think of a jet of water directed upwards with a garden hose. As long as the jet is pointing upwards with the hose (at the footpoint), the jet of water will not fall over. The jet of water will not fall over, even though the density of water is larger than that of the surrounding air. Since the footpoint of the upflowing flux tube and its inclination is constrained, the boundary condition circumvents the wavy flow to fall over. Therefore, the argument of Thomas (2005) is only true for a serpentine flow without a footpoint, and is not applicable here.
Siphon flows: Siphon flow arches are stationary magnetic flux tube models, which were proposed to explain the Evershed flow (e.g., Meyer and Schmidt, 1968;Thomas, 1988;Degenhardt, 1991;Thomas and Montesinos, 1991). This class of models makes the ad hoc assumption of different magnetic field strengths at the two footpoints of a magnetic arch, which is responsible for a gas Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 pressure gradient along the tube driving the flow. In the dynamic sea-serpent solutions (see above) a quasi-stationary solution exists . This solution corresponds to one (out of four) particular siphon flow solution: a flow with a supercritical flow speed along the arch. Since the pressure difference at the footpoints has to be assumed a-priori, this is not a fully self-consistent explanation of the Evershed flow, but these models allow to describe properties of the flow in the visible layers of the sun.
Heat transport: Temporal measurements of the intensity evolution rule out the existence of interchange convection (Solanki and Rüedi, 2003) and, also, the numerical work of the moving tube model did not confirm the concept of interchange convection of magnetic flux tubes as the heating mechanism for the surplus brightness of the penumbra: A crucial result of the numerical investigation is that a tube rises and develops an upflow, but the upflow does not stop nor does the tube sink back down to the magnetopause. Hence, instead of interchange convection the moving tube simulations suggests that the heating occurs in form of upflow channels along magnetic field lines. Ruiz Cobo and Bellot Rubio (2008) demonstrate that such an upflow is capable to account for the brightness of the penumbra and that such upflows can produce dark-cored bright filaments with a length of up to 3 Mm. Yet, even if such upflows can transport enough heat to account for the brightness of the penumbra, Schlichenmaier and Solanki (2003) have shown that downflows within the penumbra are obligatory: If each hot upflow along a tube turns horizontal and stays horizontal, there is not enough space for all the horizontal tubes. Hence, the horizontal magnetic tubes must turn downwards after a few Mms to make room for more tubes with more upflows. The downward directed part of the tube is then associated with a downflow. In this respect, the overshoot scenario (serpentine flow) may help: the hot upflow cools and the cool downflow heats up in the hot sub-photosphere, and re-enters the photosphere as a hot upflow. Hence, the moving tube scenario encounters problems in accounting for sufficient heat transport, but there are ways to solve the heat transport problem with channeled flows. And these channeled flows are driven by radiative cooling.
In present-day MHD simulations the energy transport in umbral dots, light bridges, and filaments in the inner penumbra is accomplished by a magneto-convective mode, which may be characterized as convective elongated cells, very similar to the energy transport in granulation. The following section will summarize models of sunspot fine structure based on MHD simulations with radiative transfer.

Idealized magneto-convection simulations
Magneto-convection has been investigated in idealized setups for several decades and have been reviewed by Proctor and Weiss (1982), , Schüssler (2001) and Weiss (2002). Apart from the Rayleigh number characterizing the degree of convective stability, these studies typically focus on how magneto-convection patterns change with imposed field strength (typically expressed through the Chandrasekhar number), the ratio of magnetic to thermal diffusivity = / as well as the field inclination angle. Early studies summarized by Proctor and Weiss (1982) were based on the Boussinesq approximation, followed by 2D compressible studies (Hurlburt and Toomre, 1988;Weiss et al., 1990). Here, it was found that for < 1 convection in strong field regions is oscillatory, while steady overturning motions are present for > 1. It has been conjectured by Weiss et al. (1990) that umbral dots can be explained through oscillatory magneto-convection in the sub-photospheric layers ( < 1 is realized in the upper most 2 Mm of a sunspot umbra). Matthews et al. (1995) and Weiss et al. (1996) expanded this work to 3D where convection takes place in a lattice of pulsating dots. The regime of moderately strong field was studied by Tao et al. (1998), here magnetic field separates from convective motions (flux-separation), which is realized in granulation and plage regions. Weiss et al. (2002) presented 3D studies of flux separation in Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 photospheric convection, varying the field strength while using a fixed profile that varies (bottom to top) between 2.2 to 0.2. Based on this work it was also suggested that in certain regions of the umbra an intermediate regime with flux separation on small scales is realized. Hurlburt et al. (1996) investigated 2D magneto convection in inclined field. Here oscillatory convection transitions to traveling waves that can lead to both pattern motion and average horizontal flows near the top boundary.  presented the corresponding 3D traveling wave pattern for different inclination angles. They found convection cells with a pattern motion toward the umbra, while fluid is rapidly moving outward in the wake of the traveling convection cells. It has been speculated by the authors that several aspects of penumbral structure and flows are represented by traveling wave magneto-convection.
While idealized simulations point toward oscillatory and traveling wave like convection under the condition < 1, which is realized about 2 Mm beneath the photosphere, MHD simulations with radiative transfer and a realistic equation of state (described in Section 3.6.2) show the immediate transition to overturning convection in umbra as well as penumbra. To our knowledge it has not been thoroughly studied which additional ingredient (radiative transfer, partial ionization, location of photospheric boundary away from domain boundary allowing for convective overshoot) is responsible for the change of behavior compared to the idealized models summarized above.
Recent magneto-convection studies by Thomas et al. (2002a,b), Weiss et al. (2004), and Brummell et al. (2008) focused on the role of turbulent magnetic pumping for the formation and maintenance of a sunspot penumbra. Overall, pumping was found to be very efficient in the idealized setups to hold down magnetic field lines near the outer edge of the penumbra and it was conjectured that this process together with a convective fluting instability is responsible for the formation of penumbrae as well as the Evershed flow in terms of a siphon flow in the overarching flux loops resulting from this process.
3.6 Radiative magneto-convection 3.6.1 Numerical challenges Significant progress in our ability to simulate sunspots using realistic MHD simulations (i.e., MHD simulations that include the solar equation of state and multidimensional radiative transfer) was only possible over the past couple years. This is primarily due to the fact that pursuing radiative MHD simulations on the scale of sunspots with sufficient resolution to still capture the essential scales of magneto-convective energy transport requires fairly large computational domains and accordingly computing power. Additional to the computational domain size also the physical parameters encountered in and above the umbral region of a sunspot pose significant numerical challenges. The combination of several kG magnetic field with the rather small density scale height leads to a steep increase of the Alfvén velocity above the sunspot umbra reaching values in excess of a few 1000 km s -1 . Such high velocities lead to severe time step constraints for explicit codes that make such a simulation almost impractical. A very low value of also leads to problems in codes that use the conservative formulation of the energy equation, since the determination of the internal energy requires to compute the small difference between the nearly equal values of the total and the magnetic energies. Hence, progress in this field calls for either restricting simulations to small domain sizes with focus on certain aspects of the problem, or by relaxing some of the most severe numerical constraints, primarily with respect to very large Alfvén velocities (see, e.g, appendix of Rempel et al., 2009b).

Magneto-convective modeling of umbral dots
The first realistic radiative MHD simulation of a sunspot umbra was presented by Schüssler and Vögler (2006  vertical extent located entirely within the umbra region. In order to circumvent the severe Alfvénic time step constraint, the upper boundary was placed about 400 km above the Rosseland optical depth of unity ( = 1) in the umbra (and, therefore, still beneath the level of = 1 in the quiet Sun). The average vertical field strength in the box was chosen to be 2500 G and the initial stratification as well as entropy at the bottom boundary condition were adjusted to yield a thermally relaxed state with about 17 -18% of the quiet Sun radiative energy flux.
The simulations show a magneto-convective mode that consists of instationary almost field-free upflow plumes that form naturally even in a monolithic strong field region. The plumes start like an oscillatory convection mode as it has been proposed before. But owing to the sharp drop in density near the solar photosphere, upflowing plasma has to expand and weakens the vertical magnetic field to a degree that overturning convection can set in. The result are narrow almost field-free gaps filled with overturning convection that extent about 500 km downward in depth. The typical lifetime of about 30 min is closely related to the size of the energy reservoir these plumes have access to. The energy transported by the overturning plume and radiated away near = 1 leads to a decrease of the superadiabatic gradient of the stratification, which was the original driver of the instability. Once the stratification is sufficiently stable the overturning convection stops and the field-free gap closes again. The typical lifespan and photospheric appearance of an upflow plume is shown in Figure 5.
The photospheric manifestation of the upflow plumes is a umbral dot, visible for about 25 min, showing a substructure with a central dark lane, in some cases also a more complicated y-shape dark lane. The dark lane originates from a local density and pressure enhancement above the stagnation point at which most of the mass flux turn around. As a consequence the = 1 levels are elevated, leading to a shift of the line formation height into a lower temperature region.
While the overall photospheric appearance of this convection mode might look very similar to Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 the intrusion of unmagnetized plasma in a shallow cluster model as proposed by Parker (1979a), the physical origin is fundamentally different: Almost field-free upflow plumes can form within a strong monolithic magnetic field as a consequence of convective instability (energy reservoir), strong stratification (expansion of upflows leading to weakening of field), and flux expulsion by convection within the upflow plume (cf. Weiss, 1964). In that sense the presence of umbral dots can not be used as a support for the cluster model. On the contrary, we may conceive that the umbra is an overall monolithic fully magnetic structure, in which the fine structure is a local disturbance. The dots are produced locally by magneto-convection processes, the latter being a necessity for the large amount of energy that needs to be transported.

Magneto-convective modeling of umbra/penumbra transition and sections of penumbrae
In order to simulate the transition from umbra to penumbra, it is required to laterally extend the computational domain and also include regions of granulation. The inclusion of granulation automatically requires that the top boundary of the domain is placed a few 100 km above the quiet Sun photosphere and therefore almost 1000 km above the = 1 level in the umbra. Under these conditions a severe Alfvénic time step constraint becomes unavoidable and all of the magnetoconvection simulations discussed here relax this constraint by artificially limiting the Lorentz-force in low regions. As long as the reduced Lorentz-force remains dominant over pressure forces, the overall field topology is not affected, but the computational expense is reduced by about 2 orders of magnitude.
To keep the computational expense at a minimum it is also beneficial to consider sunspots in 'slab' geometry, i.e., modeling a slender rectangular section of a sunspot. The first realistic radiative MHD simulation based on this concept was performed by Heinemann et al. (2007). They considered a rectangular section of a (slab-like) small sunspot of about 4 Mm diameter. The main result of this simulation is the formation of filamentary structures in the outer part of the spot, various properties of which (such as dark cores, inward propagation during formation phase, outflows, and strongly inclined magnetic field) are consistent with observational results However, the filaments found by Heinemann et al. (2007) are much shorter than the typical lengths of real penumbral filaments, and the overall extension of the simulated penumbra is very small. While the magneto-convection mode they identified shares many similarities with the mode found by Schüssler and Vögler (2006), umbral dots were not present in their simulation. Rempel et al. (2009b) presented a slab simulation based on the same concept, but allowing for a substantially larger sunspot of about 20 Mm diameter. Due to larger overall extent, the filamentary structure of the inner penumbra is more pronounced and individual filaments reach a length of up to 3 Mm and exhibit a clearly visible central dark lane (see Figure 6).
The umbral region also shows the formation of umbral dots similar to those found by Schüssler and Vögler (2006), which allows to clearly identify the common magneto-convective origin of both structures. Figure 7 summarizes the properties of magnetic field, temperature and velocity on a plane perpendicular to the filament highlighted in Figure 6. The vertical extent of the displayed region is 1.3 Mm and almost the same as in Figure 5, and the field and flow properties can be compared accordingly. In both cases energy transport takes place in form of hot rising plumes that lead to a strong reduction of the vertical magnetic field component, however, in the case of penumbral filaments the vertical extent of the region with reduced field strength and vertical flow is larger. The inclined magnetic field near the periphery of the spot causes a symmetry breaking, which leads to elongated filaments with strong outflows along flow tubes of nearly horizontal field near optical depth unity. In addition to the flow along the filament, the upflow also turns over into a motion perpendicular to the filament axis. Dark lanes appear above the strongest upflows owing to the upward bulging of the surface of optical depth unity and the piling up of plasma in Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Most of the energy radiated away in the photosphere is provided by the deep reaching central upflow pattern, which connects to layers with substantially larger heat capacity. Shallow roll-type convection as indicated in Figure 7 is only of secondary importance for the energy transport and is observationally indistinguishable from the deep reaching component of the flow pattern.
Based on the simulations of Heinemann et al. (2007) it has been suggested by Scharmer et al. (2008) that the Evershed flow can be identified with the horizontal flow seen in the filaments and that it is identical with the horizontal flow component of the magneto-convection in the filament channels. An average outflow of about 1 -2 km s -1 in the penumbral region is also present in the simulations of Rempel et al. (2009b). However, the overall appearance of the penumbral region seen in Heinemann et al. (2007) and Rempel et al. (2009b) seems more representative of an inner sunspot penumbra with highly intermittent filaments penetrating partially into the umbra region. Rempel et al. (2009b) also found weak inflows at the edges of filaments where the vertical flow is downward directed (see Figure 7). The indications of an outer penumbra with extended regions of horizontal field and strong radial outflows with large filling factor are very weak in slab simulations, however, they are present in simulations of full sunspots described in the following section.
Both, the simulations of Heinemann et al. (2007) and Rempel et al. (2009b) show some evidence of a weak moat flow diverging from the sunspot and transporting magnetic flux away from the sunspot. Due to the limited horizontal extent of the domain this flux accumulates near the domain boundaries and forms small pores as visible in Figure 6. Kitiashvili et al. (2009) presented recently a numerical setup that can be considered as a model for magneto-convection within a sunspot penumbra. They focused on a small section of magnetoconvection in strongly inclined magnetic field, where the inclination and field strength is imposed through the initial state and maintained through the boundary condition. They found a strong dependence of the average horizontal outflow on field strength and inclination angle. Average outflow speeds in the 1 -2 km s -1 range required about 1.5 kG field strength combined with an imposed mean inclination of 85°. In addition, they reported on temporal variations of the simulated Evershed flow in the range from 15 -40 min, that can be associated with Evershed clouds (Shine Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Living Reviews in Solar Physics http://www.livingreviews.org/lrsp- -3 et al., 1994Rimmele, 1994;Cabrera Solana et al., 2007). A further analysis of these simulations with focus on magnetic flux returning beneath the photosphere and forming so called "sea-serpent magnetic structures" was presented by Kitiashvili et al. (2010a).

Magneto-convective modeling of full sunspots
While the slab geometry can capture many aspects of the transition from umbra to inner penumbra with a minimal computational expense, simulations of circular sunspots are a better setup for realistic simulations of the outer penumbra as well as the surrounding moat region. Since the slab geometry only allows expansion of the horizontal field in one dimension there is a tendency of too strong horizontal magnetic field in the penumbral regions when applied to large spots with diameters around or in excess of 20 Mm.
A simulation of a complete sunspots is much more demanding since the combination of large domain size and the required high resolution requires large computational grid sizes. This task was recently accomplished by Rempel et al. (2009a), who performed a simulation of a opposite polarity sunspot pair in a domain of 98.304 Ö 49.152 Ö 6.144 Mm at a resolution of 32 km in the horizontal directions and 16 km in the vertical requiring a total of 1.8 billion grid points. A pair of opposite polarity sunspots was chosen as setup in order to cover a variety of combinations of field strength and inclination angles in a single simulation run.  Rempel et al. (2009a). Displayed is the sunspot pair after about 5.75 hours of temporal evolution (3.75 hours in high resolution). The sunspot on the left has more umbral dots due to the overall weaker field strength (3 kG), most extended outer penumbra can be found in the center region with the shortest separation between the opposite polarity spots. (To watch the movie, please go to the online version of this review article at http://www.livingreviews.org/lrsp-2011-3.) Figure 8 shows an intensity image from that simulation after about 5.75 hours of temporal evolution. While both spots have an identical flux (1.6 Ö 10 22 Mx), the spot on the left (right) has a central field strength of about 3 (4) kG; owing to the periodic boundary conditions the magnetic field is more inclined along the x-direction. This underlying field geometry is also manifested in the fine structure: The weaker spot shows more umbral dots, most extended penumbrae are found in the x-direction, preferentially in the center of the simulation domain. The latter is due to the fact Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3  Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 that the spot separation in the box is 43 Mm, while the separation across the periodic boundaries is 55 Mm. Figure 9 shows a magnetogram as well as subsurface field strength for the same snapshot.
In addition to umbral dots and dark cored filaments, which have been studied before in slab geometry, this simulation presents for the first time an extended outer penumbra with a strong radial outflow that has a filling factor close to unity and average velocities of up to 5 km s -1 (peak flow speeds can reach 14 km s -1 ). The location of regions with radial outflows is strongly related to the average inclination angle of the magnetic field. Figure 10 presents the relation between average flow velocity in the x-direction and field inclination in the center of the domain averaged over 8 Mm in the y-direction and 1 hour in time. Coherent outflows from the spots start where the inclination angle exceeds 45°. Note that this relation holds for the 4 penumbral regions present in this plot, despite their different appearance in intensity ( Figure 8) and is, thus, of more universal nature. The Evershed flows are also clearly visible in the animation that is provided with Figure 9.
Recently, Rempel (2011a) presented a detailed analysis of this simulation with regard to the photospheric appearance of sunspot fine structure and the physical origin of large scale outflows. Figure 11 summarizes the fine structure of the penumbra at = 1. Radially aligned filaments in intensity correspond to regions with moderately enhanced radial field and strongly reduced vertical field strength resulting in a strong variation of the inclination angle. Note that the enhancement of horizontal field is restricted to a narrow boundary layer forming around = 1. Fast horizontal outflows are present along horizontal stretches of the field. Overturning convection is found everywhere in the penumbra, upflows are preferentially located in the center of filaments. The correlations between Evershed flow, intensity, and field strength are presented in Figure 12 and are discussed in more detail in Rempel (2011b) and Rempel (2011a). While in the inner penumbra strong outflows are found preferentially in bright filaments, the correlation disappears in the outer penumbra. This is consistent with findings from recent high resolution observations Ichimoto et al., 2007a). The correlation with field strength is negative in the inner and positive in the outer penumbra, the latter was also proposed by Tritschler et al. (2007) and Ichimoto et al. (2008b) as explanation for changes in the net circular polarization (NCP) at different viewing angles.
It is found that maintaining the penumbra brightness requires overturning motions at the = 1 level with a RMS value of 1 km s -1 , about half of what is found in granulation (simulation results point toward an approximate relation of the form ∝ √︀ RMS ( = 1) between azimuthally averaged intensity and vertical RMS velocity). Fast horizontal outflows are primarily driven in a narrow boundary layer found beneath = 1 in the central upflow regions of filaments. In this boundary layer the Lorentz force facilitates the energy exchange between pressure driving in upflows and the acceleration of fluid that takes place primarily in the horizontal direction. This process is not limited to filament heads and takes place along the full length of filaments. The resulting outflow reaches a maximum velocity close to the local Alfvén velocity in the deep photosphere and declines rapidly with height (see Rempel, 2011a, for details).
Most of these features are found to be robust with respect to numerical resolution within the currently accessible range. A recent simulation with 16 km horizontal and 12 km in the vertical resolution is presented in Figure 13 and further described in Rempel (2011b).
Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Figure 11: Still from a movie showing Penumbral fine structure at = 1 level. The intensity image (panel a) shows radially aligned filaments that correspond to regions with moderately enhanced horizontal field (b) and strongly reduced vertical field (c). The result is a strong variation in field inclination (d). Fast radial outflows (e) with velocities of about 8 km s -1 within flow channels are present along horizontal stretches of field. The vertical velocity (f) shows upflows mostly along the center of individual filaments and downflows near their edges in lateral direction. Downflows are also found near the outer end, where the inclination exceeds 90°, indicating magnetic flux returning beneath the surface. Contour lines in (e) and (f) highlight outflows faster than 10 km s -1 and downflows faster than 5 km s -1 . This figure is reproduced from Rempel (2011a) by permission of the AAS. (To watch the movie, please go to the online version of this review article at http://www.livingreviews.org/lrsp-2011-3.) Figure 12: Correlations between horizontal (Evershed) flow and intensity (a) and magnetic field strength (b). In the inner penumbra outflows are preferentially found in bright features, the correlation essentially disappears toward outer penumbra and turns weakly negative. The magnetic field strength is reduced in flow channels in the inner penumbra, but enhanced in the outer penumbra.
Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 3.6.5 Unified picture of magneto-convection Altogether, the results of Schüssler and Vögler (2006), Heinemann et al. (2007), Scharmer et al. (2008), Rempel et al. (2009a,b), and Rempel (2011a,b) indicate a new level of realism in the theoretical modeling of sunspot structure. The basic properties of the simulated umbral dots and penumbral filaments are consistent with a variety of observational results and provide a basis for a physical understanding of umbral and penumbral structure in terms of a common magnetoconvective process that is modulated by the varying inclination angle of the magnetic background field.
In the almost vertical magnetic field of the umbra upflow plumes strongly expand as consequence of stratification and reduce the field strength to a degree that overturning convection sets in. The observable manifestations are bright umbral dots with short dark lanes above the central upflow. With increasing inclination angle the symmetry between the horizontal directions becomes broken resulting in elongated peripheral umbral dots and dark cored filaments of the inner penumbra. The vertical magnetic field component is diminished to almost zero field strength similar to umbral dots as a consequence of expansion and overturning motions. In contrast to that the horizontal field component is moderately enhanced in a shallow boundary layer leading to the formation of an uncombed structure of filaments with almost horizontal field and more vertical background field in between. The central upflows in the filaments are deflected outward by the inclined field in a narrow boundary layer leading to outflows with almost Alfvénic velocity in the deep photosphere.
The common element in all these regimes is overturning convection on scales much shorter than the radial extent of the penumbra, which is responsible for most of the energy and mass transport. Deep reaching upflow plumes connect the photosphere to layers with substantial heat capacity and maintain the substantial brightness of the penumbra. The presence of the inclined magnetic field imposes a large degree of anisotropy as well as preferred direction leading to the appearance of large scale organized horizontal outflows.

Critical assessment of MHD simulations to date
While radiative MHD simulations capture the essential building blocks of penumbral fine structure, they are not fully realistic yet. Major shortcomings are related to the currently affordable resolution, which is sufficient to resolve individual filaments, but insufficient to properly resolve turbulent flows within filaments, i.e., these flows are currently forced to be more or less laminar. Also the initialization of these models is quite arbitrary since a self-consistent simulation of sunspot formation is still out of reach. Since there is strong evidence from observations that many aspects of penumbral structure and large scale flows around sunspots show hysteresis, it is likely that certain aspects of penumbral structure found in MHD simulations are tied to the initial state, which is typically an axisymmetric self-similar field configuration. Current limitations with regard to resolution and initial state will be relaxed in the future with more computing power becoming available.
Overall, MHD simulations are currently very successful in explaining the penumbral brightness, the filamentation resulting in a magnetic structure comparable to the uncombed penumbra (Solanki and Montavon, 1993), and the fast Evershed outflow along the almost horizontal component of the field. Since the outflow is essentially a convective flow (Scharmer et al., 2008;Rempel et al., 2009a,b;Rempel, 2011a,b) simulations show an Evershed flow in the deep photosphere, reaching its peak velocity near = 1 and falling off rapidly with height. This is in disagreement with the investigation of Rimmele (1995) and Stanchfield II et al. (1997), which found evidence for elevated flow channels. On the other hand, the bisector of spectral lines points toward outflows in the deep photosphere that decline with height . The currently most controversial aspect is the evidence (or lack of evidence) for overturning convection in observations of the penumbra. Most of the evidence in favor of overturning convection is based on twisting Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 motions of filaments that are interpreted as overturning convection (Ichimoto et al., 2007b;Bharti et al., 2010). Márquez et al. (2006)  Recently, Bharti et al. (2011) analyzed the visibility of overturning motions in a MHD penumbra model through forward modeling of spectral lines and their degradation to the spatial and spectral resolution of observations. They concluded that the visibility of overturning motions depends strongly on the used spectral line and can be easily masked by projected Evershed flows, even if a sunspot is only a few degrees away from disk center. From their synthetic analysis, they propose that the C i 538.0 nm line is most promising to detect overturning motions. However, they did not consider the effect of line blends, present in the cooler parts of the penumbra, nor the fact that the line depression can be smaller than 5% in dark penumbral filaments, which makes a detection very difficult. Scharmer et al. (2011) andJoshi et al. (2011) claim direct evidence for overturning motions inside the penumbra using C i 538.0. Scharmer et al. (2011) finds a similar level of correlation between intensity and vertical velocity as found in the quiet Sun photosphere. The vertical RMS velocity is about 1.2 km s -1 and very similar to predictions from the above mentioned MHD simulations. While these recent results point clearly toward a convergence between observations and theoretical predictions, in both paper the results are only obtained after a stray-light correction was applied. Therefore rigorous evidence for the existence of overturning convection in penumbra is still missing.

Sunspot Formation and Evolution
The previous section focused on the detailed structure of sunspots down to the smallest currently observable scales. While many aspects of sunspot fine structure can be understood to a large degree separated from the overall evolution of sunspots, a comprehensive understanding of sunspots requires to model the processes of their formation, dynamic evolution and decay. Models that focus on the latter have to address much longer time and length scales. As a consequence, such models require major simplifications or are computationally expensive. It is currently not possible to describe the entire process of sunspot formation and decay consistently within a single model. Typically, processes in the lower and upper convection zone are treated independently due to the vast separation of time and length scales. The evolution of sunspots after the emergence process has been modeled primarily through simplified models.
For more details on the flux emergence process we refer here also to the reviews of Moreno-Insertis (1997), Fisher et al. (2000), Fan (2009), and further references therein.

Flux emergence in lower convection zone
The solar convection zone encompasses a density contrast of about 10 6 , leading to vast range of length and time scales as well as convection regimes. While the bottom of the convection zone hosts strongly subsonic flows (Ma ∼ 10 −4 ), convective motions in the solar photosphere turn supersonic. The pressure scale height ranges from about 50 Mm at the base of the convection zone to about 100 km in the photosphere, convective overturning time scales range from weeks to minutes. Modeling the deeper layers of the convection requires filtering out sound waves (to avoid overly severe time step constraints) while fully accounting for stratification, which is achieved through the anelastic approximation (see Glatzmaier, 1984, for full 3D approaches) or the thin flux tube approximation in simplified models. The assumptions underlying the anelastic as well as thin flux tube approximation break down in the upper most 10 -20 Mm of the convection zone, which requires taking compressibility fully into account.
Early models of the flux emergence process (Choudhuri and Gilman, 1987;Fan et al., 1993Fan et al., , 1994Moreno-Insertis et al., 1994;Schüssler et al., 1994;Caligari et al., 1995) were based on the thin flux tube approximation. Studying the time evolution of a closed 1-dimensional flux loop in a background stratification taken from a solar convection zone model, these studies were able to explain large scale properties of active regions, such as the low latitude of emergence, latitudinal trend in tilt angles as well as asymmetries between leading and following spots. Necessary condition for this agreement was an initial field strength at the base of the convection zone in the 100 -150 kG range. Very similar values for the field strength were also found in independent studies of flux storage in the subadiabatic overshoot region Schüssler, 1993, 1995).
Based on two-dimensional MHD simulations it was early realized by Schüssler (1979) that untwisted magnetic flux tubes cannot rise coherently and fragment. It was shown later by Moreno-Insertis and Emonet (1996) and Emonet and Moreno-Insertis (1998) that this fragmentation can be alleviated provided that flux tubes have enough initial twist.
More recently also 3D MHD simulations of rising flux tubes based on the anelastic approximation have become possible (Fan, 2008) and give support for results from earlier simulations based on the thin flux tube approximation. It was, however, found by Fan (2008) that there is a very delicate balance between the amount of twist required for a coherent rise and the amount of twist allowed to be in agreement with observations of sunspot tilt angles (twist with the observed sign produces a tilt opposite to the effect of Coriolis forces on rising tubes).
The simulations presented above consider the flux emergence process decoupled from convection. First attempts to address flux emergence in global simulations of the convection zone were made recently by Brun (2007, 2009). Understanding the interaction of emerging flux Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 with the ambient convective motions in the convective envelope is a crucial step toward more realism; however, currently the focus on the global scale limits the resolution required to resolve this interaction in detail. A complementary approach that circumvents the resolution problem was recently taken by Weber et al. (2011). They simulated the rise of a thin flux tube through a convection zone taken from a 3D global convection model, where the coupling between flux tube and surrounding convective flow field is accomplished through the drag force. Overall, the interaction between rising flux and convective flows reduces the sensitivity of results with respect to the initial field strength, although the best agreement with observed active region properties is found for about 50 kG initial field strength. The latter is a factor of about 2 -3 lower than values previously inferred from simulations not considering the interaction with convective motions.

Flux emergence in upper convection zone
Models for the flux emergence process in the upper most layers of the convection zone require fully compressible MHD since flows are approaching and exceeding the speed of sound in the photosphere. In addition, fully realistic simulations have to take into account partial ionization in the equation of state and radiative energy transport in the photospheric layers. Over the past few decades, realistic MHD simulations on granular scale evolved to a degree that allows a detailed comparison with high resolution observations (see review by Nordlund et al., 2009, and further references therein). Realistic MHD simulations of sunspot formation are much more demanding due to the substantially longer times and length scales inherent to this problem. While there has been in recent years a large body of idealized 3D MHD simulations addressing the emergence of flux from the photosphere into the corona, only very recently MHD simulations with radiative transfer and a realistic equation of state are utilized to address the last stages of flux emergence in the subphotospheric layers. Cheung et al. (2007) considered the emergence of small flux concentrations (10 18 -10 19 Mx) in the upper most 2 Mm of the convection zone, and Cheung et al. (2008) studied the flux emergence of (10 20 Mx) tubes in the uppermost 5 Mm of the convection zone. While all these cases show signatures in the photosphere that are in agreement with observations, such as highly distorted granulation, formation of supersonic downflows, separation of mixed polarities, they were not yet able to produce "pore-like" flux concentrations. In part this was due to the amount of flux considered, in part it was due to time-step and stability constraints imposed once the flux reached the photosphere and layers above, which did not allow running these simulations sufficiently long. The formation of a strong flux concentration was finally achieved by Cheung et al. (2010) building on code improvements introduced by Rempel et al. (2009a,b) and utilizing a new bottom boundary condition that accounts for a half-torus shaped flux tube entering the domain from beneath following Fan and Gibson (2003). In this case the rise of a twisted flux tube with 7 Ö 10 21 Mx flux in the upper most 7.5 Mm of the convection zone lead to the formation of 2 pores, showing fine structure in form of umbral dots as well as light bridges (see Figure 14). All these simulations clearly show a strong influence of the near surface convection on the flux emergence process, leading to the formation of small scale loop like structures that emerge on a granular scale.
A complementary approach was recently introduced by Stein et al. (2011a). Instead of prescribing initially a concentrated flux tube they emerged horizontal magnetic field across the bottom boundary in inflow regions in a 20 Mm deep domain. Experiments with different initial field strength revealed that only runs with less than 20 kG field strength at 20 Mm depth lead to flux emergence signatures in agreement with observational constraints. In this setup the overall amount of flux reaching the photosphere is too small to form an active region. Recently, Stein et al. (2011b) reported on an experiment in which the magnetic field was in addition scaled up everywhere in the domain proportional to B on a time scale of 30 min to produce larger flux concentrations leading to a complex active region.
Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Expanding the simulation domains to also include chromospheric layers and the lower corona is computationally very expensive. As a consequence many emergence simulations do not include convection and focus entirely on the buoyant rise. Recently, Abbett (2007) developed a compressible MHD code in which the radiative energy transport is approximated through a combination of radiative diffusion beneath and cooling functions above the photosphere. This simplification allows to cover the full range from the upper convection zone into the corona and allows to simulate flux emergence and interaction with convective motions. A flux emergence simulation with realistic photospheric convection but simplified thermodynamics in the chromosphere was presented by Tortosa-Andreu and Moreno-Insertis (2009). In contrast to this Martínez-Sykora et al. (2008) focus on detailed physics of time dependent ionization and non-LTE radiative transfer. While the former approaches are more efficient to address large active region scale problems, the latter is required for in-depth comparison with observations.

Formation of a penumbra
The formation of a penumbra is an integral part of the flux emergence process. To date most numerical models do not show clear evidence for penumbra formation, in part since the grid resolution in most flux emergence simulations is insufficient, the total emerging flux falls short of that of solar active regions, and as described below boundary conditions matter. The recent experiment by Stein et al. (2011b) shows some evidence of penumbra formation, although here the magnetic field strength was scaled up artificially during the simulation to reach flux concentrations of sufficient size.
Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 Based on more idealized models it has been suggested by Weiss et al. (2004) and Brummell et al. (2008) that turbulent pumping near the outer edge of a penumbra plays a key role in maintaining the overall field geometry required for the existence of a penumbra. While such processes are certainly present in most radiative MHD simulations of sunspots, their overall role has not been quantified. Most recent MHD models point to a strong role of magnetoconvection within the penumbra (not just near the edge) in maintaining the finestructure of sunspots locally. In addition, in most MHD models the extent of penumbrae is subject to boundary conditions (see, e.g., Figure 8). Rempel et al. (2009a) found extended penumbrae only inbetween (relatively close) opposite polarity spots, similarly the setup of Stein et al. (2011b) also involves opposite polarity flux in close proximity. Rempel (2011b) artificially increased the field inclination near the top boundary to obtain an extended penumbra for an individual sunspot in a periodic domain. The difficulty of obtaining penumbrae with realistic extent is most likely an artifact of the use of periodic boundary conditions in horizontal directions that put strong constraints on the allowed global field structure. A more realistic setup including parts of the corona and relaxing horizontal periodicity might be required to address this aspect self-consistently.
On the observational side high resolution observations of flux emergence including the transition from a proto-spot into a penumbra confirm the picture that emerging bi-poles separate, and that the flux patches of the proper polarity merge. Figure 15 displays snapshots that trace the growth of a sunspot as it develops a penumbra. Towards the lower right, elongated granules mark emerging bi-poles. As they separate, the flux patches with the spot polarity migrate toward the spot. As the spot increases in size it forms a penumbra (see Schlichenmaier et al., 2010a,b). Figure 15: Snapshots in the G band of a spot on July 4, 2009, as it develops a penumbra. The spot is located at = 28°and was observed from 08:32 UT until 13:03 UT at the German VTT in Tenerife (Schlichenmaier et al., 2010b). Flux emergence takes place in the lower right corner of the images. Flux patches of the spot polarity are observed to migrate towards and merge with the spot. As the area and the magnetic flux of the spot increases the penumbra forms on the side opposite to the flux emergence.

Sunspot evolution past emergence process
After the appearance of an active region in the photosphere, sunspots show a transition from an active initial phase toward a more passive later phase. In the active phase both polarities continue to separate and move away from the original emergence site, in the later passive phase sunspots do not show strong motions with respect to ambient plasma (Švanda et al., 2009). It has been suggested that this change of behavior is related to the subsurface connectivity of sunspots. The early stage with the strong separation of both polarities is the natural consequence of an Ω-shaped loop emerging into the photosphere. The related asymmetries between the leading and following spots have been well studied in thin flux tube simulations (Fan et al., 1993;Caligari et al., 1995) and are also reproduced by 3D simulations (Fan, 2008). Extrapolating these results past the appearance time of active regions in the photosphere would predict rather strong separation of polarities if the magnetic field remains rooted in the tachocline. A process that would stop the separation process, requires a dramatic change in field connectivity. Whatever the underlying process is, it has to be very reliable since we do not observe active regions with peculiar behavior in that respect. The two possible scenarios discussed in literature are reconnection of sunspots leading to shallow U-loops (Schrijver and Title, 1999) and dynamical disconnection (Fan et al., 1994;Schüssler and Rempel, 2005). The latter leads to rather shallow sunspots in which the magnetic field strength drops to sub-equipartition values at about 5 -10 Mm depth. The magnetic field would become passive to turbulent motions there, which effectively disconnects the sunspot from the magnetic root at the base of the convection zone. While this process was simulated within a time evolved onedimensional self-similar sunspot model, a self-consistent ab initio 3D simulation is still outstanding. On the other hand, it is also not clear whether a rather shallow sunspot would be stable enough to explain live times of several weeks, which is certainly longer than the overturning time scale of convection in the surface layers. Simulations of sunspots in 6 -8 Mm deep domains (Rempel et al., 2009b,a) require fixing the magnetic field at the bottom boundary to prevent a rapid decay within a few hours. Recent simulations by Rempel (2011c) in up to 16 Mm deep domains show strong evidence that this constraint can be significantly relaxed since the intrinsic convective time scales increase dramatically with depth. Rempel (2011c) found that sunspot lifetimes of about 1 -2 days can be achieved in a 16 Mm deep domain regardless of the bottom boundary condition; extrapolating this result should yield lifetimes of about 10 days for sunspots anchored in 50 Mm depth. The latter is substantially deeper than the disconnection depth that was suggested by Fan et al. (1994), Schrijver and Title (1999) and Schüssler and Rempel (2005). Addressing this problem fully consistently requires to model all stages of the flux emergence process from the base of the convection zone into the photosphere and beyond, which will likely become feasible within the next decade.
To date the decay of sunspots has been primarily addressed by simplified turbulent decay models, such as Petrovay and van Driel-Gesztelyi (1997), Rüdiger and Kitchatinov (2000) and Kitchatinov and Olemskoi (2006). With proper assumptions about the quenching of the turbulent magnetic diffusivity these models can describe the overall flux loss rates and shape of the decay curve, a more detailed modeling requires taking fully into account the interaction of the sunspot with the surrounding convection zone as well large scale flows (moat). The 3D interaction of sunspots with convective flows was recently modeled by Botha et al. (2011) in an idealized setup and by Rempel (2011c) using radiative MHD models.

Moat flows
Apart from the Evershed flow within the spot, there is also a large scale outflow pattern surrounding spots. Figure 16 shows a 15 minute time average of a sunspot flow map taken with HMI/SDO (http://sdo.gsfc.nasa.gov/) 7 . The spot is off disk center and viewed at an angle of 35°relative the local vertical. Within the spot (bounded by black contour) the Evershed flow reaches values exceeding 1.6 km s -1 . Outside the spot, the granulation shows an overall redshift (blueshift) on the limb (center) side of about 0.5 km s -1 , extending to about twice the spot radius. The area occupied by this flow patterns is called "moat region" and has typically extents out to about 2 sunspot radii. This flow patterns was first found by Sheeley Jr (1969) by tracking photospheric bright points and further studied through Doppler measurements by Sheeley Jr (1972). Harvey and Harvey (1973) studied further the behavior of the outward moving bright points and introduced the term "moving magnetic features (MMF)" to describe the time dependent magnetic structure of the moat region. They found that MMF activity is mostly limited to decaying sunspots and plays a role in transport of flux away from the spot. Most MMFs are dipolar but have on average a net flux with the polarity of the sunspot. No significant correlation with the presence of a penumbra was found. A very extensive study of the energy and mass flux of the moat region was conducted by Brickhouse and LaBonte (1988), who found that the radial extent of the moat region is about twice the spot diameter. The average moat flow velocity from their study is about 500 m s -1 . Recently, the connection between the Evershed flow and moat flow has been addressed by several authors (see, e.g., Sainz Dalda and Martínez Pillet, 2005;Cabrera Solana et al., 2006;Vargas Domínguez et al., 2008;Zuccarello et al., 2009;Vargas Domínguez et al., 2010), however, so far the observational evidence is not conclusive to either proof or disproof a connection.
Helioseismic measurements of moat flows were performed by Gizon et al. (2000). They used inversions of the f -mode, which senses mainly the upper most 2 Mm of the convection zone and found good agreement with photospheric measurements. A more sophisticated inversion based on all ridges from f to p4 was presented by Gizon et al. (2009Gizon et al. ( , 2010b. The near surface flow agrees well with MMF tracking, and outflows with comparable amplitude were found down to a depth of 4.5 Mm. Recently, Featherstone et al. (2011) presented a helioseismic inversion combining ringdiagrams with three different resolution levels. They found that the moat flow consists of a two components, a superficial photospheric flow and a deeper reaching flow peaking at about 5 Mm depth.
On the theoretical side, Meyer et al. (1974) suggested that the moat flow is essentially a supergranular flow, which arises as sunspots are preferentially formed in a supergranular downflow vertex. Heat flux blockage by the extended penumbra leads to a reversal of the flow direction in the proximity of the spot turning the initially converging flows into outflows.
2D axisymmetric simulations by , Botha et al. (2006), and Botha et al. (2008) typically produce a 2 cell flow around sunspots, a converging flow in their proximity (which has been referred to as "collar flow") and a diverging "moat" flow further out. Recently, this work has been expanded to 3D simulations in cylindrical geometry by Botha et al. (2011). While the basic results from the 2D axisymmetric simulations are confirmed, the collar flow breaks up into several cells in the azimuthal direction and allows for sunspot decay due to turbulent erosion. Since both, the 2D axisymmetric and 3D cylindrical simulations do not capture the effect of a penumbra, it is not clear how this collar flow would be further modified in a more realistic sunspot model. Zhao et al. (2010) suggested that a superficial Evershed flow is just added on top of it.
Large scale flows are also present in recent 3D MHD simulations of sunspots (Heinemann et al., 2007;Rempel et al., 2009a,b;Rempel, 2011a,c). Here outflows dominate the picture at all depth levels, i.e., there is no compelling evidence for the presence of a collar flow in simulations of sunspots with penumbrae. An outflow is also the flow response one would expect from the penumbral blockage of heat flux as suggested by Meyer et al. (1974). Converging collar flows are, however, found in simulations of pores (Cameron et al., 2007b;Kitiashvili et al., 2010b) and play a role in their formation and stabilization against decay. A converging flow around pores is also seen in observations (Wang and Zirin, 1992;Sobotka et al., 1999;Vargas Domínguez et al., 2010).
Recently, Cheung et al. (2010) simulated the formation of a pair of small (3 Ö 10 21 Mx) sunspots through a flux emergence process in the upper most 7.5 Mm beneath the solar photosphere. Due to the large amount of high entropy plasma that emerges, a buoyantly driven large scale outflow surrounds the sunspots from early on in the simulation. Interestingly a coherent pair of spots was forming out of initially dispersed field despite diverging (azimuthally averaged) horizontal mean velocities. The analysis of the simulation revealed that mostly the transport of field due to the electromotive force resulting from small scale fluctuating motions dominated the spot formation, while large scale flows had a destructive effect throughout the process. The small-scale velocity field correlations allow in principle for a transition between the growing and decaying phase of sunspots without requiring changes in the large scale flow pattern as originally proposed by Meyer et al. (1974). However, this is still an untested hypothesis and numerical simulations are currently just progressing to the point at which processes on the scale of the moat region or entire active regions can be properly addressed. Substantial progress is likely over the next decade.

Sunspot Models and Helioseismic Constraints
Since the advent of local helioseismology, the subsurface structure of sunspots (thermal, magnetic, and flows) is among the prime targets. Independent investigations of subsurface flow and field structure are of great interest, since we cannot infer from surface observations alone whether a sunspot is monolithic or better described by a cluster model. In this section we only highlight very briefly recent developments and refer to the reviews by Gizon and Birch (2005) and Gizon et al. (2010a) for details.
A major difficulty in applying helioseismic inversions to sunspots comes with the fact that sunspots are not small perturbations with respect to the surrounding quiet Sun -at least not for the upper most few Mm. This intrinsically limits the applicability of linear inversion methods most commonly used. In addition, the effect of magnetic field cannot be easily captured by a simple isotropic wave speed perturbation due to its directional character. Direct application of time distance helioseismology to sunspots revealed a two layer structure of sunspots. Interpreting travel time shifts as wave speed perturbation, Kosovichev et al. (2000) found a reduction of wave speed in the upper most 4 Mm and an enhancement of wave speed further down, extending in depth beyond 10 Mm. The relative amplitude of wave speed perturbation was found to be up to 5% translating to either a substantial temperature change of more than 2000 K or a magnetic field almost reaching 20 kG just a few Mm beneath the photosphere. A similar 2 layer pattern was also found by Zhao et al. (2001) with respect to subsurface flows: an inflow from about 1 -4 Mm depth and an outflow beneath. These results were also summarized in a recent review by Kosovichev (2006). Complementary, also ring diagram analysis has been applied to active regions (Basu et al., 2004;Bogart et al., 2008), providing results that are qualitatively comparable to the previously mentioned time-distance inversions. However, the transition from reduced to enhanced wave speed is found at a depth of more than 7 Mm.
These results have been challenged recently by a long list of comparative studies summarized in Gizon et al. (2009Gizon et al. ( , 2010b and Moradi et al. (2010). Most of the observed travel time shifts can be attributed to the upper most 2 Mm, where relative perturbations are large with respect to the quiet Sun reference stratification. Independent support for a rather superficial thermal anomaly in sunspots comes also from the quasi-1D approach of Lindsey et al. (2010). Their method focuses on the frequency dependence of travel time shifts for waves that propagate almost vertical in sunspots (and are reflected in the umbra) to separate thermal and magnetic effects. A very promising new approach that has been taken by several groups in parallel is the forward modeling of wave propagation through stationary sunspot models that automatically addresses the complicated nature of the wave propagation in inclined magnetic field (Khomenko and Collados, 2006;Cameron et al., 2007a;Parchevsky and Kosovichev, 2007;Hanasoge, 2008;Shelyag et al., 2009). By modeling the propagation of f, p1, and p2 -modes, Cameron et al. (2011) demonstrated that the observed helioseismic signatures of sunspots can be well captured by a rather shallow model. Evidence for a dominant influence from near surface effects was reported previously by Braun and Birch (2006) based on observations of the frequency dependence of p-mode travel times in sunspots. Also recent 3D simulations of sunspots (Rempel et al., 2009b,a) can be explored to infer helioseismic signatures. In these simulations the full spectrum (to the degree it is allowed by the finite size of the domain) of modes is naturally excited by convection and can be analyzed similar to observed modes by extracting artificial Dopplergrams on constant height or constant optical depth surfaces. As presented in Moradi et al. (2010) and Braun et al. (2011), fully compressible 3D radiative MHD models of sunspots predict a rather shallow wave speed perturbation and are in terms of travel times consistent with observations. Overall, there is now very strong evidence, both observationally and theoretically, that a detailed understanding of near surface effects will be essential for making progress in sunspot seismology.

Summary and Conclusions
Over the past few decades there has been a tremendous improvement of sunspot models, in part driven by improvements of instrumentation leading to unprecedented observational constraints on models, in part driven by the increase of computing power finally allowing radiative MHD simulations of entire sunspots. While for a long time penumbral fine structure was only accessible through simplified models (assuming certain field geometries or convective flow structures), in recent years numerical simulations evolved to a degree where they start capturing the essential elements of sunspot fine structure.
To explain different aspects, two "simple" model classes have been proposed for the penumbra: (a) Flux tube models, which assume that the magnetic field structure associated with penumbral filaments can be described as flux tubes located close to the = 1 surface. In these models flows channeled by magnetic fields account for the filamentation, the Evershed flow, and the line asymmetries. While these models are most successful in reproducing spectropolarimetric measurements they have difficulties to reproduce the overall downflow in the outer penumbra and to account for all of the energy transport needed to explain the penumbral brightness. (b) The gappy penumbra model, which makes the strong assumption that regions of strong magnetic field and convective energy transport are spatially separated. This model assumes field-free gaps in which energy transport by sufficiently deep reaching elongated convection cells. This model does not give an explanation for the Evershed flow (i.e., the model does not require the existence of strong horizontal flows), but such could be included. However, these flows would be essentially nonmagnetic and it has not yet been demonstrated whether the gappy model can explain the observed spectropolarimetric signature. With respect to the energy transport this class of models differs fundamentally from flux tube models. While in the former only the footpoints provide hot plasma to the photosphere, in the latter upflows along the full length of the filaments are responsible for the energy supply.
Radiative MHD simulations reveal a full spectrum of magneto-convection regimes depending on the strength and inclination angle of the background field. Some of the convection patterns found in simulations have some resemblance of features found in simpler models. However, none of the simpler models can give a good characterization of magneto-convection in the penumbra as a whole. While umbral dots are essentially field-free upflow plumes (with much similarity to the gappy models), filaments in the inner penumbra remain magnetized with field strength around 1 kG. The outer penumbra with strong horizontal outflows is best characterized by anisotropic magnetoconvection with no clear separation between convection and magnetic field. The manifestation of field inclination, field strength and horizontal flows close to photospheric layers show much similarity with flux tube models. However, the latter do not include overturning convective motions throughout the penumbra (except for the footpoints), which is the key element in radiative MHD simulations (as well as the gappy model).
Yet, at this point, it has not been fully investigated if radiative MHD simulations are consistent with all observational facts. In particular, it remains to be seen whether the flow pattern of convective rolls can be measured, and whether the observed penumbral line asymmetries in the Stokes parameters including the NCP can be reproduced by such models. Spectropolarimetric measurements need to have a spatial resolution of better than 0.1 arcsec to be comparable to the models. On the other hand, models need to be able to explain the substantial brightness of the penumbra, which has not been achieved by models not relying on overturning motions throughout the penumbra. Making an unambiguous connection between model results and observations requires the forward modeling of spectral lines and detailed comparison with observations, work in that direction is currently in progress.
Detailed modeling of formation, evolution, and decay of sunspots remains at this point an open problem. However, MHD simulations are now in the process of approaching the relevant scales Living Reviews in Solar Physics http://www.livingreviews.org/lrsp-2011-3 and it can be expected that within the next decade substantial process will be made in this regard. Photospheric MHD simulations will reach domain sizes that can host entire active regions and domain depths, which allow coupling to flux emergence simulations in the lower convection zone.