.. . ..x. ..é...rwu~..._, , aucpuhuun. . . ... .v. . . $M§§Wj£1 . r1319. H. J.- rhm‘ ,. . .V F“ .1.“ ‘1“. . . .11 1.4.3.... .2: ...: ...... . “my-4J4. .411..- .. 12.1. n: .41: .. .... \t. 1:. ... . litlg. ‘~‘$“A“ 1.15.1311? .... 2.1!. .1111 um... 51‘ ...-“flan. {.11. I 1:81. . . 21....h. W..." $1.11.!uicleslx I 12 $1136. 111...! «01111:... I 1.1521... 11.3.1155. I An ...... 91. zialkfiin ....Mrifi .1“ V. :1... 5.. 3.: 22:11 hm...) . ..Ifl {1.12.1 ... 1:4! . . i.” 1' .. 1230.1: ...... 319‘}:- 1§IIA§l «. :th i . .151... 1 -..... ( .... ... -3133. 1:. .tlv East). {..J.v§1..nls )4! . 6.3.5535... .1“. .2 "ea-v.4: gaging ., . . ... b. .1 . .. Lg.» 1. fl 1W 1 . v . . .. (.... firhfifig , ,. . a... s. 13H... . 2 . , .. . . u ....uflvatmna 2.... 7 .../v This is to certify that the dissertation entitled A Numerical Investigation of the Interaction Between Convection and Magnetic in a Solar Surface Layer presented by DAVID JOHN BERCIK has been accepted towards fulfillment of the requirements for Ph. D . degree in Phfiics % 121% Mam/professor R. Stein Date February 15 , 2002 MS U it an Affirmative Action/Equal Opportunity Institution 0-12771 LIBRARY Michigan State University PLACE IN RETURN Box to remove this checkout from your record. TO AVOID FINES return on or before date due. MAY BE RECALLED with earlier due date if requested. DATE DUE DATE DUE DATE DUE Aim-10’: 91526335 6/01 c:/CIRClDateDue.p65-p.15 A NUMERICAL INVESTIGATION OF THE INTERACTION BETWEEN CONVECTION AND MAGNETIC FIELD IN A SOLAR SURFACE LAYER By David John Bercik A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Physics and Astronomy 2002 ABSTRACT A NUMERICAL INVESTIGATION OF THE INTERACTION BETWEEN CONVECTION AND MAGNETIC FIELD IN A SOLAR SURFACE LAYER By David John Bercik Three-dimensional numerical simulations are used to study the dynamic inter- action between magnetic fields and convective motions near the solar surface. The magnetic field is found to be transported by convective motions from granules to the intergranular lanes, where it collects and is compressed. A convective instabil- ity causes the upper levels of magnetic regions to be evacuated, compressing the field beyond equipartition values, and forming “flux tubes” or “flux sheets”. The degree to which the field is compressed controls how much convective transport is suppressed within the flux structure, and ultimately determines whether the mag- netic feature appears brighter or darker than its surroundings. For this reason, the continuum intensity is not a good tracer of the lifetimes of magnetic features, since their bright/dark signature is transient in nature. Larger magnetic structures form at sites where a granule submerges and the surrounding field is pushed into the resulting dark hole. These micropores are devoid of flow in their interior and cool by radiating radially. The convective downflows that collar the micropore heat its edges by lateral radiation, but fail to penetrate far enough into the interior to prevent an overall cooling, and therefore darkening, of the micropore. Magnetic features undergo numerous mergers or splittings during their lifetimes as a result of being pushed and squeezed by the expansion of adjacent granules. Larger structures survive for several convective turnover times, but smaller structures are too weak to resist convective motions, and are destroyed on a convective time scale. ACKNOWLEDGMENTS This dissertation would not have been possible without the help of a number of people. I would like to thank my advisor, Bob Stein, for his continued support, encouragement, assistance and patience throughout my graduate student career. I am grateful to Bob and Ake Nordlund for allowing me to use and modify their convection code, without which this work would not have been possible. I would also like to thank Shantanu Basu and Dali Georgobiani for their con- tributions to getting the MHD version of the code up and running. Many thanks to Axel Brandenburg for sharing his insights about magnetic fields and his hospi— tality during my trip to Newcastle England. I appreciate the conversations about theoretical astrophysics that I have had with Sydney D’Silva, Bill Abbett, Regner Trampedach and Jeff Kreissler. I am also grateful to them for the non-astrophysical conversations that helped me to retain my sanity. I am appreciative for the extended use of computer resources at Michigan State University and the National Center for Supercomputing Applications. I also would like to express my thanks to NSF and NASA for the grant support provided to work on this project. A special note of thanks is extended to my family for their unwavering support and encouragement, without which I would never have accomplished this goal. iv TABLE OF CONTENTS LIST OF TABLES vii LIST OF FIGURES viii 1 INTRODUCTION 1 2 NUMERICAL METHODOLOGY 9 2.1 Equations ................................. 10 2.2 Time Advance ............................... 13 2.3 Spatial Derivatives ............................ 15 2.4 Diffusion .................................. 17 2.5 Boundary Conditions ........................... 21 2.6 Input Physics ............................... 26 2.6.1 Equation of State ......................... 26 2.6.2 Radiative Transfer ........................ 26 2.7 Initial Conditions ............................. 29 3 GENERAL PROPERTIES OF THE MAGNETIC FIELD 32 3.1 Location of the Magnetic Field ..................... 32 3.2 Magnetic Field Distribution ....................... 38 3.3 Magnetic Field Correlations ....................... 43 4 CONVECTION VS. MAGNETO-CONVECTION 53 5 MAGNETIC FLUX TUBES 62 5.1 Structure ................................. 62 5.2 Energy Budget .............................. 78 6 MICROPORES 86 6.1 Formation ................................. 87 6.2 Evolution ................................. 91 7 CONCLUSIONS 101 7.1 Summary ................................. 101 7.2 Future Work ................................ 103 A Derivative Schemes 105 A.1 First Derivatives ............................. 105 A11 Horizontal Direction: Uniform, Periodic Grid .......... 105 A.1.2 Vertical Direction: Non-Uniform, Non—Periodic Grid ...... 109 A.2 Second Derivatives ............................ 117 A.2.1 Horizontal Direction: Uniform, Periodic Grid .......... 117 A22 Vertical Direction: Non-Uniform, Non-Periodic Grid ...... 120 vi LIST OF TABLES 1 Statistics for granular/intergranular areas. See text for details. . . . 55 vii LIST OF FIGURES Mean atmosphere for the initial state, averaged horizontally. The solar surface is at depth 2 = 0, where < T >2 1. The atmosphere spans five pressure scale heights above the surface and six below. It is stably stratified in the photosphere, above z=-0.1 Mm and unstably stratified in the convection zone below the visible surface ........ 30 Image of the emergent intensity with overlay of the horizontal flow field at the surface, for the 400 G run. Granules have a horizontally diverging flow ................................ 34 Image of the vertical velocity at the surface with contours of the magnetic field at the surface. Contours are in increments of 200 G. Same time and parameters as Figure 2. The magnetic field is swept into the intergranular lanes by the diverging granular flow and for this much flux nearly fills the lanes .................... 35 Image of the vertical velocity at z = 2.5 Mm with contours ofthe mag- netic field at z = 0. Contours show B > 1 kG in 100 G increments. The surface field is concentrated over the downflow boundaries of the underlying mesogranules .......................... 37 Distribution of vertical magnetic field for the 200 G run (top) and the 400 G run (bottom) at z = 0. Bin size is 10 G. Dotted line is best hyperbolic fit. ............................ 40 Distribution of horizontal magnetic field for the 200 G run (top) and the 400 G run (bottom) at z = 0. Bin size is 10 G. Horizontal fields produced by the stretching and twisting of the predominantly vertical fields have an exponential distribution. ................. 42 Distribution of vertical magnetic flux through individual grid points at at z = 0 for the 200 G run (top) and the 400 G run (bottom). . . . 44 Correlation of the magnetic field and density (top) and the magnetic field and temperature (bottom) at the visible surface for the 400 G run. Locations with strong field are cooler and evacuated. ...... 45 viii 10 ll 12 13 14 15 16 17 18 19 Correlation of the magnetic field and vertical velocity (top) and the magnetic field and horizontal velocity (bottom) at the visible surface for the 400 G run. Strong fields suppress both horizontal and vertical velocities. ................................. Correlation of the magnetic field and density (top) and the magnetic field and height (bottom) at unit optical depth for the 400 G run. . . Correlation of the magnetic field and vertical velocity (top) and the magnetic field and horizontal velocity (bottom) at unit optical depth for the 400 G run. Both horizontal and vertical velocities are sup- pressed by strong fields. ......................... Correlation of the magnetic field at z = 0 and emergent intensity for the 200 G run (top) and the 400 G run (bottom). High field locations are darker than average. ......................... Correlation of the magnetic field at z = 0 and its inclination, 7, for the 200 G run (top) and the 400 G run (bottom). Strong fields tend to be vertical. ............................... Typical emergent intensity images. The images have been replicated in each direction for illustration purposes. In the presence of large magnetic flux, granules become more rounded, with more diffuse edges and wider intergranular lanes .................... Images of the vertical velocity at z = 0 for the same times as the images in Figure 14. Downflows are bright. In the presence of large magnetic flux, fast downfiow lanes remain narrow ............ Distribution of the emergent intensity, summed over ten minute in- tervals. Large magnetic flux reduces both the brightness of typical granules and the darkness of the typical intergranular lanes. ..... Granulation size spectrum given by the emergent intensity power, averaged over ten minute intervals. The wavenumber has units of Mm‘1 .................................... Emergent intensity on a horizontal slice as a function of time. The granular intensity structures have longer lifetimes as the magnetic flux increases. ............................... Volume rendering of the magnetic field. Below the surface (layer 15) the field is collected into “flux tubes” by the diverging flow. Above the surface the field controls the flow and spreads out .......... ix 50 51 54 56 60 64 20 21 23 24 25 26 27 28 29 Background: filled contours of the magnetic field in 250 G intervals. Overlay: contours of the natural logarithm of the density ranging from -2.0 (topmost) to 4.0 in 0.5 intervals. The T = 1 depth is shown as the thick line around 2 = 0 Mm. The established, strong “flux tube” in the center has been evacuated. The smaller flux concentrations on either side are in the process of being evacuated, starting above the surface and piling up plasma below the surface. ............ Background: filled contours of the current density in 0.1 A/m2 inter- vals. Overlay: same as in Figure 20. A strong current sheet surrounds the “flux tube” in the center. ...................... Background: filled contours of the magnetic field in 250 G intervals. Overlay: contours of the temperature from 4000 K to 16000 K in 1000 K intervals. The 7' = 1 depth is shown as the thick line around 2 = 0 Mm. The “flux tubes” are significantly cooler than their surroundings at a given geometrical depth. ...................... Background: filled contours of the current density in 0.1 A / m2 inter- vals. Overlay: same as in Figure 22. .................. Background: filled contours of the magnetic field in 250 G intervals. Overlay: velocity field. The 7' = l depth is shown as the thick line around 2 = 0 Mm. The flow is suppressed in the interior of the strong established fiux tube in the center, but not the evacuating “flux tubes” on either side ......................... Background: filled contours of the current density in 0.1 A/m2 in- tervals. Overlay: same as in Figure 24. The downflows are confined to the current sheet at the boundary of the strong established “flux tube” in the center. ........................... Image of vertical velocity (red and yellow down, blue and green up) with contours of the magnetic field strength at 0.5 kG intervals at the surface. .................................. Image of vertical velocity (red and yellow down, blue and green up) with contours of the magnetic field strength at 0.5 kG intervals at the 1 Mm below the surface. ......................... Plasma beta along magnetic field lines. Solid lines indicate strong field regions, dashed lines indicate weaker field regions. ........ Horizontal slices of the magnetic field at 500 km intervals of depth. White/ yellow indicates strong field regions, black / blue indicates weak field regions ................................. 71 77 79 30 31 32 33 34 35 36 37 38 Radiative heating and cooling. Units are 1010 erg g"1 s“. The top three panels show the net radiative heating/ cooling and its relation to the temperature and magnetic fields. The bottom three panels show the contribution of vertical, inclined and nearly horizontal rays. The interior of the strong established “flux tube” is nearly in ra- diative equilibrium, with cooling by vertical rays nearly balanced by horizontal heating from the side walls. ................. 80 Total net heating and cooling by all processes (top two panels) and the radiative contribution (bottom three panels). Scale is the same as for Figure 30 ............................... 83 Total net heating and cooling by all processes (top two panels) and non-radiative contribution (bottom three panels). Scale is the same as for Figure 30 ............................... 84 Correlation of intensity with magnetic field strength as a function of time when a micropore forms. The lowest intensity first occurs at intermediate field strengths. Then the field strength in the micropore gradually increases to a maximum. ................... 88 Images of magnetic field (right panels), intensity (center panels) and mask showing only low intensity, strong field locations (left panels). Micropores form when magnetic flux from surrounding intergranular lanes is swept into locations where a small granule disappears. . . . . 89 Evolution of the central point of a micropore at different depths in the convection zone. The flow reverses direction at the surface at t = 0 minutes. Top panel: magnetic field. Middle panel: vertical velocity. Bottom panel: temperature fluctuations. See text for details. 92 Evolution of the central point of a micropore at different depths in the convection zone. The flow reverses direction at the surface at t = 0 minuntes. Top panel: density fluctuations. Bottom panel: gas pressure fluctuations. See text for details ................. 93 Evolution of the central point of a micropore. The flow reverses direc- tion at the surface at t = 0 minutes. Comparison at different depths in the convection zone of the total acceleration, arr, the acceleration due to the Lorentz force, a3, and the sum of the accelerations due to the pressure gradient force, ap, and the buoyancy force, ac. See text for details .................................. 94 Time — horizontal slice through the center of a micropore showing an image of the emergent intensity and magnetic field contours. ..... 97 xi 39 40 41 42 43 44 45 46 47 48 49 50 51 Two dimensional horizontal — time visualization of surface magnetic field (left) and emergent intensity (right) during formation of a mi- cropore. Vertical scale is time in half minute intervals. Size of the box is approximately 2x2 Mm ....................... Time — horizontal slice through the center of a micropore showing an image of the magnetic field strength for the same region as Figure 38. Amplitude errors relative to the analytic solution for first derivative schemes. Symbols represent sampled wavenumbers. .......... Gaussian profile used for horizontal error analysis ............ Error profile of the horizontal divergence of a Gaussian function for the fourth-order first derivative scheme .................. Error profile of the horizontal divergence of a Gaussian function for the sixth-order first derivative scheme. ................. Contour of the error profile of the horizontal divergence of a Gaussian function for the fourth-order first derivative scheme. Dotted lines are shown as aids in visualizing asymmetry .................. Contour of the error profile of the horizontal divergence of a Gaussian function for the sixth-order first derivative scheme. Dotted lines are shown as aids in visualizing asymmetry .................. Amplitude errors relative to analytic solution for second derivative schemes. Symbols represent sampled wavenumbers. .......... Error profile of the horizontal Laplacian of a Gaussian function for the second-order second derivative scheme ................ Error profile of the horizontal Laplacian of a Gaussian function for the sixth-order second derivative scheme ................. Contour of the error profile of the horizontal Laplacian of a Gaussian function for the second-order second derivative scheme. Dotted lines are shown as aids in visualizing asymmetry. .............. Contour of the error profile of the horizontal Laplacian of a Gaussian function for the sixth-order second derivative scheme. Dotted lines are shown as aids in visualizing asymmetry. .............. Images in this dissertation are presented in color. xii 98 99 108 110 111 113 114 119 121 122 Chapter 1 INTRODUCTION Magnetic fields play an integral role in many astrophysical processes, from the accel- eration of particles in galaxies, to the formation of stars and stellar activity, down to the evolution of planets. The Sun provides a unique laboratory to study astrophysical magnetohydrodynamics. It is the only stellar object close enough to resolve many of the features associated with magnetic activity, and it provides the only guide to pro- cesses that are not directly observable on other stars. For these reasons, the study of magnetic fields has been of vital interest in solar research ever since it was determined that sunspots were magnetic features, almost a century ago (Hale 1908). Ultimately, an understanding of the evolution of magnetic features is sought—linking the gener— ation of magnetic fields, to their appearance at the visible surface, to interactions in the upper atmosphere, and finally to their destruction. A variety of magnetic features are observed in and above the solar photosphere ranging in size from the cool, dark sunspots at 104 — 105 km down to bright magnetic elements which are on the order of 100 km in diameter. Pores, magnetic knots and micropores are intermediate-sized objects with diameters of 103 — 104 km, 103 km and 500 km, respectively. Regions of the solar surface are often classified into categories based on the total amount of magnetic flux they contain. Active regions are sites of emerging flux where sunspots are formed, and they also contain high mean—strength field regions known as plages. The area outside of active regions is collectively labeled the “quiet” Sun. While having a low average mean-field strength, these regions are by no means free of field; they are sites of a network of magnetic elements comprised of high flux concentrations. The abundance of magnetic features seen on the surface is a function of a 22 year activity cycle that reverses the polarity of the magnetic field every 11 years. The most favored theoretical paradigm is that the combination of convection and differential rotation drives a dynamo that generates magnetic field. This field is stored in a subadiabatic region at the base of the convection zone, where it gets amplified (Schiissler 1993, Riidiger 1994). When the magnetic field strength reaches a critical value, an instability causes a loop of field to enter the convection zone where it buoyantly rises and emerges at the surface (Moreno-Insertis 1994). Convective motions in the photosphere are responsible for compression of mag- netic “flux tubes” and play an integral part in their evolution. The buffeting of “flux tubes” by nearby granules generates transverse magnetohydrodynamic waves (Steiner, Knolker, & Schiissler 1994, Steiner et al. 1998). These waves dissipate their energy in the chromosphere and provide one possible source of heating in this region (Ofman, Klimchuk, & Davila 1998). Convective motions also shift the footpoints of coronal magnetic loops, twisting field lines and resulting in magnetic reconnection. This process heats the corona through resistive dissipation (van Ballegooijen et al. 1998, F urusawa & Sakai 2000). Furthermore, magnetic features may account for vari- ations in the Sun’s luminosity during the solar activity cycle. These changes in solar energy output could be responsible for climatic changes on Earth (Cubash & Voss 2000) It is possible to gain insight into solar subsurface and atmospheric processes by studying the interaction of magnetic fields and convection in the photosphere. Over the past decade, advances in instrumentation such as the Near-Infrared Magneto- graph (Rabin 1992), the Ziirich Imaging Stokes Polarimeter (Keller et al. 1992) and the Advanced Stokes Polarimeter (Elmore et a1. 1992) have provided vector magne- tographs, better spatial resolution and the opportunity to study magnetic configura- tions on scales smaller than supergranulation. These advances have allowed observers to not only study high flux active regions(Stolpe & Kneer 1997, Lites, Skumanich, & Martinez Pillet 1998), but also the quiet Sun network (Grossmann-Doerth, Keller, & Schiissler 1996, Sigwarth et a1. 1999, Grossmann-Doerth et a1. 2000, Sénchez Almeida & Lites 2000). Imaging techniques have been used with magnetograms to identify the existence of small magnetic elements by their associated bright points (Zhang et a1. 1998, Muller et al. 2000). Speckle reconstruction methods have been used to achieve the best spatial resolution ( 300 km) of magnetic features (Keller 1995, Koschinsky, Kneer, & Hirzberger 2001). Recent observations have shown evidence for a diffuse weak-field component to the surface magnetic field, in addition to the kilo— gauss “flux tubes” seen at supergranulation boundaries (Lin 1995, Meunier, Solanki, & Livingston 1998, Lin & Rimmele 1999, Stolpe & Kneer 2000). These findings could be evidence that a secondary turbulent dynamo is responsible for generating the weak—field component of the magnetic field seen in the photosphere (Cattaneo 1999). Unfortunately, observations cannot provide all the answers. Direct imaging of magnetic features is degraded by distortions introduced by the Earth’s atmosphere, and present—day telescopes are unable to fully resolve the smallest solar features. Therefore, spectropolametric methods are used as an indirect way to study magnetic fields. Difficulties arise with this technique also. The amount of polarization mea- sured is generally limited to several percent, resulting in the need for longer exposure times compared to unpolarized investigations, and limiting the minimum observable flux. Polarization spectra must furthermore undergo inversion procedures that rely on prescribed model atmospheres. The reliability of the derived quantities depends on how well these models represent the actual solar atmosphere. Numerical simulations play a complementary role to observations. Simulations can be used to study the inherently complex interaction between magnetic fields and convection on length scales or depths not accessible to direct observation. In turn, observational results are used to constrain and validate the simulations. Numeri- cal models of convection and magnetoconvection have provided us with much of our understanding of photospheric magnetic features. Early numerical studies of magne- toconvection were done in two dimensions using the Boussinesq approximation, where density fluctuations are only accounted for in the buoyancy force, making them nearly incompressible. Furthermore, these studies were kinematic, as the back reaction of the magnetic field on fluid motions was not included in the calculations. Convective motions were found to wind up the initially uniform magnetic field and concentrate it at the edges of the convective cell (Weiss 1964, 1966, Proctor & Weiss 1978), sim- ilar to the process of ‘flux expulsion’ described by Parker (1963). Three-dimensional studies have found the magnetic flux concentrated into an isolated tube by the sur- rounding convective flow (Galloway, Proctor, & Weiss 1978, Galloway & Moore 1979, Galloway & Proctor 1983). The maximum field strength found in these types of simulations is greater than the value expected from an equipartition of kinetic and magnetic energies. Dynamical simulations in the Boussinesq approximation, includ- ing the effects of the Lorentz force, show that the convective cells expel nearly all of the magnetic field into concentrated flux structures (Peckover & Weiss 1978, Proc- tor & Galloway 1979, Weiss 1981a, 1981b). The flux regions are separated from the convective regions by a thin current sheet, and fluid motions inside the flux regions are suppressed. These simulations have been useful in their ability to explore the incompressible magnetoconvective parameter space. The importance of compressibility in a stellar atmosphere, however, has been shown in simulations of non—magnetic convection. An asymmetry exists between upflows and downflows (Graham 1977, Nordlund 1984b, Porter & Woodward 1994) with the downflows being stronger as a result of pressure fluctuations on the density (Hurlburt, Toomre, & Massaguer 1984). Stratification and mass conservation lead to larger convective cell sizes with increasing depth below the surface (Nordlund 1985, Nordlund & Stein 1996). It has also been shown that solar convection is driven at the granulation scale by downflows as a result of radiative cooling in a thin boundary, rather than as an energy cascade from large to small eddies (Stein, & Nordlund 1989, Spruit, Nordlund, & Title 1990,Nordlund & Stein 1996, Spruit 1997,Stein & Nordlund 1998). As computing power increased, simulations of fully compressible magnetocon- 5 vection became feasible. A number of studies have used simplified physics to explore the parameter space. It was found that magnetic field is expelled from convective cells into concentrated flux sheets, which are partially evacuated. The motions inside the flux sheets are suppressed, while the sheets themselves are collared by downflows (Hurlburt & Toomre 1988). Other studies have investigated various parameters to gain a better understanding of strong field regions. Weiss et al. (1990) found that the nature of convective flow depends upon field strength and the ratio of magnetic to thermal diffusivity. At a given ratio, stronger field leads to steady convection while weaker field leads to oscillatory motions (Weiss et al. 1996). The structure of convective cells also depends on the angle of the imposed field (Hurlburt, Matthews, & Proctor 1996), as well as the aspect ratio of the computational domain (Blanch- flower, Rucklidge, & Weiss 1998, Tao et al. 1998). Simulations using more realistic physics (but at low spatial resolution due to the large computational effort required) have provided a comparison to non-magnetic convection(Nordlund 1984a, Nordlund & Stein 1990b, Stein, Brandenburg, & Nordlund 1992,Bercik et al. 1998). The magnetic field gets concentrated in the intergranular lanes, forming “flux tubes” or “flux sheets”. The magnetic field affects the convective structure in turn, altering the shapes of the granules and intergranular lanes. The present work is a higher resolution, larger domain extension of these simulations. Besides general properties of magnetoconvection, much effort has gone into the investigation of isolated magnetic structures. Early work in this area investigated the temperature structure inside a flux tube using a magnetostatic model (Spruit 1976). The tube radius and heat exchange are important factors in determining the tube’s structure. Other investigations concentrated on how radiative transfer affected flux tube structure, although these models had no convective energy transport (Ferrari et al. 1985, Kalkofen et al. 1986, 1989, Steiner & Stenflo 1989). Dynamical simulations in two dimensions have explored the temperature and magnetic field strength of flux structures (Deinzer et al. 1984a, 1984b, Kniilker, Schiissler, & Weisshaar 1988, Knolker & Schiissler 1988, Grossmann-Doerth et al. 1989, Knolker et al. 1991, Grossmann-Doerth et al. 1994). These models showed strong external downflows and an average continuum intensity signature close to the quiet Sun average. Due to large viscosities, these simulations could not develop non-stationary convection. More recently, two—dimensional simulations have been used to study the interaction of flux sheets and non-stationary convection (Steiner et al. 1996, Steiner et al. 1998). The convective flow buffets the flux structure causing it to sway back and forth. Strong flows occur outside the sheet, and shocks are seen to propagate upwards inside the sheet, leaving a noticeable signature in Stokes V profiles. The magnetic field is swept to the periphery of granules where a strong downflow partially evacuates the magnetic region, compressing the field to kilogauss field strengths (Grossmann- Doerth, Schiissler, & Steiner 1998). The dependence of the structure on the total flux for ideal axisymmetric flux tubes shows good agreement with observed pores (Hurlburt & Rucklidge 2000). The field becomes more inclined at the periphery of the tube as the total flux increases, and surface flows converge on the tube. This work will investigate the interactions of magnetic fields and convective motions in a near-surface layer of the Sun through a series of realistic numerical simulations. A description of the code and numerical method is given in §2. General properties of the magnetic field are discussed in §3, followed by the differences between convection and magneto-convection (§4.) The structure of flux concentrations is described in §5. In §6 the formation and evolution of micropores is discussed. A summary of the results and plans for future work are given in §7. Chapter 2 NUMERICAL METHODOLOGY The numerical approach is based on the hydrodynamic code of Nordlund and Stein (Stein & Nordlund 1998, Nordlund & Stein 1990a), modified to include magnetic fields. The physical region under investigation is a thin solar surface layer. The dimensions are 12 Mm x 12 Mm x 3 Mm, extending 500 km above the surface to the temperature minimum and extending 2.5 Mm below the surface into the convection zone. The simulations are run on a three-dimensional non-staggered Cartesian grid. The horizontal directions have uniform spacing with 125 grid points, while the vertical direction has non-uniform spacing with 63 grid points. This gives a grid spacing of 96 km in the horizontal directions and 35 — 87 km in the vertical direction. The calculations are highly computer intensive and only feasibly run on super- computers. The simulations were run on multiple processors on an Origin 2000 at the National Center for Supercomputing Applications (NCSA) and an Enterprise 6000 at Michigan State University. Even so, the higher magnetic field strength runs take about one day for each solar minute. The following sections describe the numerical techniques utilized in the code. 2.1 Equations First, consider the purely hydrodynamical case in the absence of magnetic flux. Fluid motions are governed by the conservation of mass, momentum and energy. The simulation domain is a highly stratified medium spanning eleven pressure scale heights and nine density scale heights. A nonconservative form of the conservation equations is used, i.e., the equations are not of the conservative form 8F 5+V-(Fu), (1) where F is the conserved quantity. The logarithms of the pressure and density are used, since these quantities have a nearly linear variation with depth in most of the domain. The equation of mass conservation is 81np_ at ——u-Vlnp—V-u, (2) where p is the mass density and u is the fluid velocity. Conservation of momentum is given by ('3 P 1 #:_u.vu+g——V1nP+—v.r, (3) at p p where P is the gas pressure, g,- = gdiz is the constant gravity and 1' is the viscous stress tensor. 10 The energy conservation equation is .85:- = —u - Ve — SV - u + Qrad + Qvisc, (4) where e is the internal energy per unit mass, Qrad is the radiative heating and Qvisc is the viscous dissipation. In a highly conducting medium like the solar atmosphere, fluid motions gener- ate electromagnetic fields. These fields evolve according to Maxwell’s equations and Ohm’s law. We consider these equations in the reduced form of the magnetohydro- dynamic approximation. Since velocities are small compared to the speed of light, all terms 0(u2/c2) are ignored. The equations of interest are then VXB = ”OJ, (5) 8B VXE = —0—t, (6) v B = o, (7) v E = 0, (8) J = a(E+u>(2ref) and serves to reduce the viscosity in compressions in the deeper layers where fluid motions are slower and there are no shocks. The reference height Zref is taken to be a height just below the visible surface. The term 2 Altai = Z [ui+Ai—l - ui+Ail+ (32) Ai=—l represents the sum of compressional (positive) velocity jumps over four zones. This acts to spread out the contribution over the entire shock, not just the steepest part. The last term is proportional to the fiducial sound speed (Equation 20) and stabilizes sound waves and weak shocks. The direction dependent diffusive scheme has the advantage that the diffusion is proportional to the grid spacing in each direction, thus allowing less damping with higher resolution grids. Another advantage is that the diffusion is increased only in the direction(s) necessary; e.g., a strong compression in one direction increases diffusion in that direction only, without increasing the diffusion in the perpendicular directions. The inclusion of viscosity, albeit strictly numerical in nature, leads to dissipation of energy. The effect of the viscosity is to prevent the buildup of mechanical energy on the smallest scales by converting it to heat. The energy dissipation that corresponds 19 to adopted viscosity scheme is 3 3 821,- V Bu.- _ V Bu,- {<—)1<—)<.;>1}, The magnetic field is stabilized by a hyperresistive scheme having the form 3A,- 1 . . (a?) = —§’I(0? + OZM' (17‘é J # k), (34) diffusive where 77 is the resistivity coefficient and a? is the hyperresistivity enhancement factor |A2J.-| 77 : - . . 0' scale height) on the boundary through derivatives, so it is only necessary that the values prescribed for this layer be reasonable. The mean density in the fiducial layer is adjusted with the net mass flux through the boundary to maintain hydrostatic equilibrium and satisfy the continuity equation, HE = 99 — Atf—ng—l (39) where Hp = (P)/((p)g) is the pressure scale height and brackets denote horizontal averages. The density fluctuations at the boundary are then scaled and added to the mean value in the fiducial layer keeping the vertical derivative hydrostatic. The vertical velocity derivative at the top is set to zero (stress free) by requiring 1.11 = 112. The internal energy is set to a constant value which is allowed to evolve slowly in 22 time through the extrapolation el = 0.95eold,1 + 0.05 [(62) — (M) (22 — 21)] , (40) 27—22 where (82) and (67) are averaged horizontally. The value of 61 is re-evaluated at 30 second intervals, allowing the model to reach a relaxed state. The condition on the vector potential at the top of the domain is a potential field (J = 0) extrapolation from the magnetic field values at the boundary. The field is taken to be decreasing at infinity; this condition means that the Fourier components of the vector potential are just scaled by the factor 6‘”?21 )k, where k = W is the horizontal wavenumber. The bottom boundary is more complicated. Since the motion is mostly laminar, the effect of the boundary on convection should still be small. In order to facilitate the study of acoustic waves (radial p—modes oscillations), a node in the mass flux is applied at the boundary so that there is a node in the radial p—mode waves. A horizontally uniform boundary pressure condition is enforced to stabilize the non- radial modes. The pressure at the boundary is determined by adiabatically adjusting the density and internal energy such that the total pressure, P101 2 P + 32/2110, is a constant equal to the horizontal mean, (PM), Inp —+ Imp—[P—(—%)] (as?!) . (41) 0 s 32 Be e ——+ e—[P—(—flo-)](5p) (42) P —> (PW—2%. (43) 23 To keep the influence of the boundary on the solution as small as possible, variable values are only modified at those points at which fluid is entering the domain; the variable values at outflow points are unchanged. The variables at the inflow points are evolved to their desired values with a time constant (44) This time constant is half the time it takes for the fastest incoming fluid to traverse the boundary zone. The entropy, S, of the incoming fluid is taken to be isentropic, and is evolved toward a constant entropy at constant pressure, e ——> e — £135 (32),), (45) Be Blnp lnp ——) lnp—n—AS(-a—S)P ( Be >P, (46) BS BS AS = (e — ebot) (52)!) + (p - pbot) (5)6 , (47) where At is the timestep defined in Equation 16 and pbot and ebot are constants in space and time that give the correct effective temperature at the surface. The incoming velocities are evolved toward flat profiles, At “:1: -_) Us: (1_ —) a (48) At 11,, —> uy (l - t_—), (49) At uz —> 11,. + t—'-( — us). (50) 24 The term (uzjn) is the horizontal average of incoming velocities. The vector potential is evolved toward a vertical magnetic field, (BAx/Bz, BAy/Bz) ———+ 0, AZ —> 0, At Axmz ) Axmz + (Axum—1 "‘ Ax,nz)T, (51) At Aymz —'* Aymz + (Aymz-l — Ay.n2)?‘a (52) t A zero mean mass flux through the bottom boundary is enforced by altering the pressure gradient, Ausznz-l —> Auz,l:nz-l — <—P‘_> aln(P)) _ aln

(54) (p) 82 (puz)nz:0 82 Uz,nz __> Um; _ At(P)nz (Bln(P)) _ Bln(P) , (55) pm Bz (Wand, B2 "2 where brackets denote horizontal averages. The final step in the bottom boundary routine is to ensure that total mass is conserved. The value of the total mass (f(p)dz) is compared to its original value, and if it is not within the specified tolerance level, then the density of all points in the boundary plane is adjusted iteratively until the following condition is satisfied, 93 < 3 x10“7. (56) (p) 25 2.6 Input Physics The applicability of the solutions that the simulations generate depends not only on the accuracy of the numerical scheme, but also on the realism of the input physics. It is desirable to be a realistic as possible in order to quantitatively compare the results with observations. However, realism comes at a heavy computational price. A balance must therefore be struck between the complexity of the physics and the feasibility of completing a simulation run in a reasonable amount of time. 2.6.1 Equation of State The system of equations given by Equations 13— 15 is incomplete. An equation of state (EOS) which relates the pressure to the density and internal energy is required to close the system. In the surface region under consideration, approximately two- thirds of the internal energy is in the form of ionization energy with the other one- third in the form of thermal energy (Stein & Nordlund 1998). The EOS in the code is generated using the Uppsala stellar atmosphere package (Gustafsson 1973) and includes ionization and excitation of hydrogen and other abundant elements and molecules present in the upper convection zone. The pressure, temperature and their derivatives are stored in tabular form as functions of the internal energy and the logarithm of the density. 2.6.2 Radiative Transfer The primary energy transport mechanism changes from convection below the photo— sphere to radiation above. It is the local deviations from radiative equilibrium that 26 cause the entropy fluctuations that drive the convection. This transition is at optical depth unity, and both the optically thin and optically thick (diffusion) approxima- tions to radiative energy exchange break down. In addition, highly evacuated “flux tubes” in the photosphere exchange heat with their surroundings mainly by radiative means. For these reasons, it is desirable to include detailed radiative transfer in the model. The three—dimensional treatment of radiative transfer in the simulations is as— sumed to be in local thermodynamic equilibrium (LTE), which only breaks down near the upper boundary at the temperature minimum. The radiative heating in the energy conservation equation (Equation 15) is the difference between absorption and emission, Qrad = 4Wp[\ [Qt-6,411.9 — SA)deA, (57) where It) is the monochromatic absorption coefficient, [1,9 is the monochromatic intensity and SA is the monochromatic source function (which is the Planck function, BA, for LTE). The intensity is related to the source function through the equation of radiative transfer, (ll/pg d7) = Im — 5,. (58) The monochromatic optical depth, 7'), is defined as 7') = pm sec de, (59) where 0 is the inclination from the surface normal. 27 It would be impossible to solve the transfer equation for the millions of spectral lines in the solar atmosphere, along a multitude of rays. Some approximations are therefore necessary to make the problem tractable. To take into account the effects of spectral lines without having to solve the transfer equation a prohibitive number of times, the wavelengths are sorted into four bins according to opacity (continuum, weak lines, intermediate lines, strong lines). The Planck function is then averaged in each bin. The Uppsala stellar atmosphere package is used to calculate tables con- taining the absorption coefficients and pseudo-Planck functions. The binned opacities used in these calculations are the opacity distribution functions of Gustafsson et al. (1975). All opacities are assumed to scale with height in the same way, which is reasonable except for the weak iron lines (/citeNord82, /citeNord+Stein90). The transfer equation is solved only for that part of the computational domain where T < 300. This subvolume is interpolated to a grid the same size as the original grid in order to increase the resolution of the radiation calculations. The radiative heating is found by solving a modified Feautrier equation (Nordlund 1982) along a vertical ray and four straight inclined rays. The azimuthal angle of the inclined rays is rotated 15° each time step to avoid a preferred direction of heating. Even with these approximations and computational tricks, the radiation routine still takes up to 50% of the total run time. For a detailed treatment of the radiation transfer scheme, see Nordlund 1982. 28 2.7 Initial Conditions The initial state for the simulations was created from a snapshot taken from a previous lower resolution, hydrodynamic simulation. The snapshot was then interpolated to the desired resolution, 125 x 125 x 63, and the desired physical size, 12 Mm x 12 Mm x 3 Mm. This new state was run with no magnetic field until it was relaxed thermally. The horizontally averaged mean atmosphere of the initial state is shown in Figure 1. The initial state was used to start three simulation runs, each differing by the magnitude of the initial vertical field strength. All imposed fields were unipolar. l. Bro = Byfl 2 32,0 2 0 at each grid point. This is a control run representing the purely hydrodynamical state. 2. 83,0 = Bw = 0, 32,0 2 200G at each grid point. This case represents a moderate field strength plage region. 29 Initial Atmosphere I I I I I I I I I I I I T I I I I I T I I I I I I I I I I 1000.00 E ,_..__-------""",'/' 1 - 20 l} ,/ :( >. :3 ./ A .2: *1. -' .0 (I) t: v-I : 1 J V cu 10.00 g: 15 a, Q E 1. .. $23 : 1‘ . 5 -. J “5 93 ' s a 1.00 .— 4 a. 8 E / ’/ _10 E L. : / —— Pressure a CL _/ ‘. ,’ '\ 5f - — - Density 0.10 :— '-. l E \ “i ----- Temperature '. I ., . ' 1'5 _.-._ E t . 1" n ropy 4 5 0.01 K-L’I 1\I l l 1 l l i l l l l i l l l l I l l L l I l J l 1 —0.5 0.0 0.5 1.0 1.5 2.0 2.5 Depth (Mm) Figure 1. Mean atmosphere for the initial state, averaged horizontally. The solar surface is at depth 2 = 0, where < 7‘ >= 1. The atmosphere spans five pressure scale heights above the surface and six below. It is stably stratified in the photosphere, above z=-0.1 Mm and unstably stratified in the convection zone below the visible surface. 30 3. B“, = 89.0 = 0, 32,0 2 400G at each grid point. This case represents a high field strength plage region. Each of the scenarios was run for approximately three hours of solar time. 31 Chapter 3 GENERAL PROPERTIES OF THE MAGNETIC FIELD The evolution of the magnetic field in a gas is determined by the relative contribution of the two terms in the induction equation (Equation 15). A measure of this is given by the magnetic Reynolds number, RM, which is the ratio of the advection term to the diffusion term. The high conductivity of the solar atmosphere leads to advection dominating over diffusion, RM >> 1. The magnetic field is effectively “frozen” into the fluid (Alfvénll950) and transported by convective motions. Current computational resources limit the feasible resolution in the simulations such that RM is much smaller than that of the Sun, where PM ~ 104 on the scale of granulation (Parker 1979). The simulations are in the realm RM > 1, however, and the results compare well to observations. In this chapter, we examine some of the general properties of the magnetic field; where the field is located, what its magnitude is, and how it relates to other physical quantities. These properties give insight into the magnetic features to be discussed in later chapters. 3.1 Location of the Magnetic Field As a consequence of mass conservation and the density stratification, the flow pat- tern in the convection zone is horizontally divergent in the upflows(Stein & Nordlund 32 1998). The result is that fluid, and therefore magnetic field, is swept to the gran- ular boundaries as first described in incompressible models (W'eiss 1964). In the intergranular lanes the strongest convergence of flow is at the vertices where several granules meet, so these vertices become the sites of the strongest field concentrations. Figure 2 shows the horizontal flow field at the surface. This process, known as flux expulsion (Parker 1963, Weiss 1981a, Weiss 1981b), is one mechanism responsible for concentrating the magnetic field. That the magnetic field is excluded from granules and is confined to the intergranular lanes is evident from Figure 3, an image of the vertical velocity overlayed with magnetic field contours representing the same time snapshot as Figure 2. On a larger scale, the flow converges toward stronger downflow regions, which penetrate deeper into the convection zone. These downflows are longer-lived than typical granular lifetimes of 10 — 15 minutes, and are associated with the largest- scale convective structures at the bottom of the computational domain (/citeStein89, /citeStein+Nord98). For convenience, this size scale will be referred to as mesogran- ulation; this is not to imply a distinct convective cell size, but rather as a means to differentiate between the horizontal size of the upflows at the bottom of the simula- tion and those at the surface. Mass conservation requires that most of the gas flowing up through a granule must exit through its sides over approximately a scale height. By this argument, the cell radius, r, can be estimated to be uh r ~ H—. (60) 33 Figure 2. Image of the emergent intensity with overlay of the horizontal flow field at the surface, for the 400 G run. Granules have a horizontally diverging flow. 34 2 Figure 3. Image of the vertical velocity at the surface with contours of the magnetic field at the surface. Contours are in increments of 200 G. Same time and parameters as Figure 2. The magnetic field is swept into the intergranular lanes by the diverging granular flow and for this much flux nearly fills the lanes. 35 At 2.5 Mm, the larger scale height, H, leads to larger cells. The magnetic field has more time to collect at these locations before getting convectively disrupted. The larger magnetic field structures in the simulations, such as micropores, are found at these sites. Figure 4 shows an image of the vertical velocity at z = 2.5 Mm, with magnetic field contours at the surface and reveals that the strong surface fields occur at the boundaries of the underlying mesogranules. However, flux expulsion cannot explain the peak field strength of approximately 2 kG at the visible surface. There is a back-reaction of the magnetic field on the velocity due to the Lorentz force when the magnetic energy density approaches the kinetic energy density. At this point, flow is suppressed, and further field compression diminishes. This flow suppression can be seen in the two dark micropores in Figure 2. The maximum field strength attained by this process is not expected to be appreciably greater than the equipartition field strength given by the condition 33, 1 Z _ 2 2,1 2,011 . (61) 3 and u = 2kms‘l Using typical values for the photosphere of p : 3 x 104g cm— implies that we should not expect field strengths much greater than 400 G by means of flux expulsion. One of the consequences of flux expulsion is that regions where B m Beq are no longer supplied with the convective heat necessary to balance radiative energy losses, and cool as a result. This cooling reduces the gas pressure in the region, initiating a downflow. The region is thermally isolated from the surroundings and the adiabatic 36 A Figure 4. Image of the vertical velocity at z = 2.5 Mm with contours of the magnetic field at z = 0. Contours show B > 1 kG in 100 G increments. The surface field is concentrated over the downflow boundaries of the underlying mesogranules. 37 downflow results in the region becoming cooler than the superadiabatic surroundings. This superadiabatic effect (Parker 1979) drives a convective instability that further accelerates the downflow and cooling, and a partial evacuation develops. A pressure balance is maintained by compressing the magnetic field until the magnetic pressure is sufficiently large to counter further contraction. The largest field strength possible by this convective collapse process occurs when a completely evacuated region is in horizontal pressure equilibrium with the surrounding gas pressure, or B2 —=Pas. 2 For the photosphere, this leads to a maximum field strength of approximately 2 kG, which is consistent with observations and numerical simulations. 3.2 Magnetic Field Distribution The vertical magnetic field is constrained to have a constant average field strength of either 200 G or 400 G, in order to be comparable to moderate to strong plage regions. This is a convenient way to set the amount of magnetic flux that passes through the domain. A consequence of this choice is that the number of points at a given field strength cannot be significantly altered without affecting the rest of the distribution. If the flux level is low enough that the evolution of the magnetic field is controlled predominantly by convective processes, the system is well represented by 38 a hyperbolic distribution of the form P(-.r) = (c0+c11:)'1. (63) This can be seen in the distribution for the weaker average field strength run (top plot, Figure 5). The distribution is truncated at low and high values of the magnetic field. At large B2, the field strength is limited to values less than that required to maintain a horizontal pressure balance. Flux in the direction opposite to the initial flux (negative B2) is generated as the field gets twisted around. Only relatively weak fields that are easily manipulated by fluid motions are affected in this manner, limiting the maximum field strength. The distribution of these opposite polarity fields decrease exponentially, as seen in other advected quantities. An increase in the amount of flux threading the surface increases the probability for high field strength points. A greater filling factor of high strength fields leads to filling more of the intergranular lanes with magnetic flux, to a clumping of the field and the formation of magnetic structures large enough to resist being destroyed by convective processes. By damping out fluid motions through the back reaction of the Lorentz force, these structures can survive for many convective turnover times. The net effect is an overabundance of high field strength points when compared to a hyperbolic distribution (bottom plot, Figure 5). The distribution is nearly flat for the for the majority of field strengths. The same conditions that lead to the cutoffs in the lower flux case still exist in the higher flux case. Horizontal field is generated as convection shuffles around vertical field, causing it 39 Frequency Frequency Figure 5. Distribution of vertical magnetic field for the 200 G run (top) and the 400 0.1000; __ I 1 0.0100; —_ E .. 3 0.0010; .. . _ .......... -= 000010.”... a- ‘ 0 500 1000 1500 B. (G) 400 G 0.1000;— -; 0.0100; — 0.0010; ------------------------------- 1: E 1 0.0001”....._- .. ‘ 0 500 1000 1500 B. (G) G run (bottom) at z = 0. Bin size is 10 G. Dotted line is best hyperbolic fit. 40 to bend and flex as it gets bumped by expanding granules and turbulent downdrafts. The horizontal field is not constrained directly, but through the constraint on the vertical field and zero divergence of the field. The weaker the magnitude of the magnetic field, the more easily it is influenced by fluid motions, and the greater the chance that it will be horizontal. The distribution of the magnitude of the horizontal field for both runs is shown in Figure 6. The distributions are exponential. The stronger average vertical field strength case has a larger maximum horizontal field, at the expense of low horizontal field strength points. Since each grid point in the simulations represents the same area, the distribution of vertical fields can be used to determine the distribution of vertical fluxes per grid point. The total flux is set by the initial conditions N T =< Bz>ZA=< Bz>NA, (64) i=1 where N is the total number of points and A is the area represented by each point. The fraction of the total flux in bin i is then given by (Di __ Bani/4 (PT _ < B; > NA, (65) where n,- is the number of points in bin i. The distributions of (Di/(1)7 are plotted in Figure 7. The consequences of having a hyperbolic distribution are evident here. A hyperbolic distribution has the property of being scale-free, so that smaller values contribute as much to the total as larger values. The regions in which a hyperbolic 41 I 1 0.1000 I I IIIIIII 111111 J 1 1111111 0.0100 I I I IIIIII Frequency 0.0010 1 1 1111111 1 11 I I I IIIIII 0.0001 0 200 400 600 800 r 0.1000 0.0100 I IIIIIITr I ITIIIIII Frequency 0.0010 I DI IIIIII 1 1 0.0001 0 200 400 600 800 Bh (G) Figure 6. Distribution of horizontal magnetic field for the 200 G run (top) and the 400 G run (bottom) at z = 0. Bin size is 10 G. Horizontal fields produced by the stretching and twisting of the predominantly vertical fields have an exponential distribution. 42 is a good fit will therefore have a flat flux distribution. In the 200 G run smaller flux values contribute more to the total flux than in the 400 G run, where the increase in higher field strength values leads to a dominance of higher flux values. 3.3 Magnetic Field Correlations The thermodynamic and kinetic structure at the surface (2 = 0) is weakly dependent on the magnetic field for low values of the field strength (Figure 8 and Figure 9). As the magnetic field surpasses the equipartition strength, its magnitude is controlled by pressure balance, and therefore correlates well the thermodynamic quantities density and temperature. The evacuation in strong magnetic regions leads to a density that is an order of magnitude less than that of typical intergranule regions, and about five times smaller than that of granules. These regions are also cooler than the typical intergranule lane by 1000 — 1500K. The majority of low field strength regions are upflows and regions with greater than equipartition strength are confined to downflows. As the Lorentz force starts to dominate the evolution of the dynamical structure the flow is inhibited, both vertically and horizontally, until it is under 1 km s‘1 for the highest field strengths. This is especially apparent in pores, where the flow inside is nearly non-existent. Now consider the surface of optical depth unity. All points on this surface have similar thermal structure since at any given time there is approximate radiative— convective equilibrium. As a result the spread in temperature and density will be smaller than for a surface at constant geometric height. As an example, the density at unit optical depth is shown in Figure 10. The surface of unit optical depth is 43 x 200 G 30.10005'”: r LL. : : '3‘ 1 1 *5 0.0100; g E— : 3 “—1 " .. 0 t 5:: 0.0010;— 2 .2 E 3 +9 _ O _ 20.0001 f“ 0 5 10 15 Flux (10m Mx) >< 400 G 3 0.10005 ' ' ' ' ' ' i I 'T I T LL. : 2 g 0 0.0100; —_ 9* 3 3 “5 t 1 :1 0.0010; 5 O : : 3:3 : : C) _ 1 20.0001 f“ 0 5 10 15 Flux (10“ Mx) Figure 7. Distribution of vertical magnetic flux through individual grid points at at z = 0 for the 200 G run (top) and the 400 G run (bottom). 44 .0 Lia‘uuhiiiiiiiiliiiiin 1111111111111111 'r (103 K) 1 DC a.Mllllllllllllllllllllllllllllll111111111llllllllllllllll Figure 8. Correlation of the magnetic field and density (top) and the magnetic field and temperature (bottom) at the visible surface for the 400 G run. Locations with strong field are cooler and evacuated. 45 1 .1 J i 4 q 4 -l 4 J J o C a. ‘5 a. 2' , \. -, i; :"c-' " . O 5 ' \ i. .0 ¢ 00 s -1 1 L1 1 14144 Uz (km s") 0 311111111]; —4 ..‘éi 7'” 0.5 1.0 1.5 B (kG) 2:0 0“ ‘ o.\: : '. '. u .5 ....... u p...q....e..:.' Uh (km s") a . 4 ldllllllllll llll11111111lllllllllllllllllll1lllllllllllll Figure 9. Correlation of the magnetic field and vertical velocity (top) and the magnetic field and horizontal velocity (bottom) at the visible surface for the 400 G run. Strong fields suppress both horizontal and vertical velocities. 46 4 -4 O % 1111111111111 p (10‘7g cm"°) I 1 -4 d q q —( 1 1 3.0 77" .q d O. ..— \ _. ... —4 D. d 2‘. q -1 d — d d 1 E A 2 I V -' Z N : d —1 d 4 q d d -01 r" ' :1 ‘ ....1....1....11L111141111111.1" 0.5 1.0 1.5 2.0 2.5 3.0 B (kG) Figure 10. Correlation of the magnetic field and density (top) and the magnetic field and height (bottom) at unit optical depth for the 400 G run. 47 corrugated in height; hotter regions (granules) have T = 1 above 2 = 0, whereas cooler regions (intergranular lanes) have 7' = 1 below 2: = 0. As can be seen in Figure 10, the height at which T = 1 is well correlated with the magnetic field strength. The evacuation and cooling in magnetic regions causes the level of unit optical depth to occur deeper in the atmosphere (Wilson depression). The depression reaches a maximum of approximately 300 km in the strongest field structures. The velocities at 7' = 1 are shown in Figure 11. The results are an extension of those at z = 0; the 7' = 1 surface goes deeper into the convection zone where the pressure is greater, and stronger magnetic regions can develop. The flow is further suppressed at larger field strengths. At 3kG, the horizontal flow is restricted to be less than a few hundred ms‘l, and the vertical flow is less than 1km s“. As shown in Figure 3, fields in excess of about 100 G are found predominantly in intergranule lanes. The emergent intensity of the intergranule lanes ranges from average brightness to 30% below average (Stein & Nordlund 1998). The emergent intensities for both the lower and higher flux runs are shown in Figure 12. Practically all points with moderate to high field strengths are darker than average, most indis- tinguishable from their intergranular surroundings. Those points that are of average brightness or brighter are distinguishable because of the 20% contrast with the in- tergranule lane. The points that comprise micropores are also readily apparent since they have I < .6 < I >. Magnetic structures with weak field points are more or less at the mercy of con- vection, and are easily buffeted by the expansion of granules. Weak fields are therefore found with inclinations ranging from vertical to horizontal (Figure 13). As the field 48 . 1 . O '1 1 1_1 1 1 l 1 l '1 1L11111 . 3.0 d d 4 —4 —( 1 J _1 4 d lllllllllllllllllll 1111llllllllllllllllllllllllllllllllllllll h ‘. 3.0 Figure 11. Correlation of the magnetic field and vertical velocity (top) and the magnetic field and horizontal velocity (bottom) at unit optical depth for the 400 G run. Both horizontal and vertical velocities are suppressed by strong fields. 49 11111411 I/ 111111111111 93 o 0.0 0.5 1.0 1.5 1.2 411111111 1.0 I/ «...-fiv-.." .ku—w‘q ‘0 a V ‘ 0.8 0.6 0.4 11111111111 .0 O O 01 1—1 0 y... 01 I“ C Figure 12. Correlation of the magnetic field at z = 0 and emergent intensity for the 200 G run (top) and the 400 G run (bottom). High field locations are darker than average. 50 1111;11111 7 (deg) 11111111 '9 O U! y—n O p—s 01 5" O 80 60 . . : .‘ c . ”k: .oz" .Oh. 2:. ..’¢ . o . J, . ., .§' .. . o . ‘ 7 (deg) L11111111111_1114111 0.0 0.5 1.0 1.5 B (kG) 93 c Figure 13. Correlation of the magnetic field at z = 0 and its inclination, 7, for the 200 G run (top) and the 400 G run (bottom). Strong fields tend to be vertical. 51 strength increases, “flux tubes” become more rigid and are able to resist being moved about by fluid motions. In addition, horizontal pressure balance causes the density to decrease as the field strength increases, making regions of larger field strength more buoyant. High strength fields are predominantly vertical, with maximum inclinations of 20° from the vertical. The moderate strength fields with appreciable inclinations tend to be located at the borders of magnetic features. From these correlations, we can predict that a strong-field region will be located in the intergranular lanes. It will have a relatively low density, will be cooler than the surrounding fluid and will have its velocity suppressed. Furthermore, it will be darker than average and nearly vertical. 52 Chapter 4 CONVECTION VS. MAGNETO—CONVECTION Convection plays an important role in concentrating magnetic field, but in the process undergoes changes itself. Observations of plage regions have shown granulation to have a lower contrast, smaller size scale and to evolve more slowly than in the quiet Sun (Title et al. 1992, Sobotka, Bonet, & Vazquez 1994). In this chapter we examine how magnetic field affects the convective structure, and check how well the results compare to observations. A comparison of typical emergent intensity images from the three simulation runs is shown in Figure 14. The most striking detail is the reduction in contrast when magnetic fields are present. Purely hydrodynamical granulation exhibits a sharp boundary between granules and intergranular lanes. The boundary contrast is softened with the introduction of moderate to strong magnetic fields for two reasons. First, lower intensity granules are present. Larger granules that become isolated and collared by magnetic field have decreased brightness. Inside of mesogranular boundaries, where clusters of granules can be virtually field free, normal granulation exists. Second, evacuated areas of moderate downflow containing magnetic field are brighter than typical down'flow regions in the non-magnetic case. The optical depth in these downflows is depressed allowing radiation from hotter gas deeper in the 53 Figure 14. Typical emergent intensity images. The images have been replicated in each direction for illustration purposes. In the presence of large magnetic flux, granules become more rounded, with more diffuse edges and wider intergranular lanes. 54 Table 1. Statistics for granular/intergranular areas. See text for details. 0G 200G 400G < bright area/dark area > 0.93 0.87 0.84 < downflow area/upflow area > 0.63 0.75 0.98 atmosphere to be seen. These downflows are adjacent to granular boundaries, and any increase in brightness serves to make the edges of the granules less distinct. Another visible difference concerns the shapes of granules. Ordinarily, gran- ules have a polygonal shape determined by the dynamics of their local environment. Isolated granules that are surrounded by magnetic field tend to be more rounded. The shape and size of the intergranular regions also exhibit changes with increasing amounts of magnetic field. As the magnetic field gets compressed in the intergranu- lar lanes, it essentially becomes a wall to incoming material. If convection is still to exist, outflows from granules that cool must still descend, and therefore a downflow forms between the magnetic region and the granular boundary. The result is that the intergranular region becomes wider, but the width of the downflow itself is similar to the non—magnetic case (Figure 15). Statistics related to the ratio of the total granular to total intergranular areas are listed in Table 1. The first row is the ratio of the number of points having an intensity greater than the mean intensity to the number of points having an in- tensity less than the mean value. The results are averages over ten minute intervals. There is a trend toward more dark points with increasing field strength. The second 55 Figure 15. Images of the vertical velocity at z = 0 for the same times as the images in Figure 14. Downflows are bright. In the presence of large magnetic flux, fast downflow lanes remain narrow. 56 row of the table gives the ratio of downflow points to upflow points, also averaged over ten minute intervals. Non-magnetic convection has granules occupying approx— imately 60% of the total area, whereas the presence of strong magnetic field makes the intergranular regions larger and the downflow area rivals the upflow area. The aforementioned reduction in contrast is evident in the distributions of emer- gent intensity (Figure 16). The distribution for the non-magnetic case is bimodal; the bright peak represents the distribution of granules, while the dark peak represents intergranular lanes. Magnetic fields have two effects on the distribution. The domi- nant effect is that there are many more intermediate brightness features that smooth the distinction between the bright and dark peaks. The growth in the number of average brightness points is at the expense of both bright and dark points. Less in- tense granules account for the loss of brighter points, and the edge brightening of the intergranular lanes is responsible for the loss of dark points. The other effect caused by magnetic fields is the addition of a tail to the dark end of the distribution. When there is enough flux to form extended magnetic features, these regions can become 20—30% darker than non-magnetic intergranular lanes. This second effect therefore varies more greatly with time than the first. Besides affecting the brightness of the granulation, observations of magneto- convective regions on the solar surface have shown that magnetic fields also affect the size distribution of the granulation (Title et al. 1992, Sobotka, Bonet, & Vazquez 1994. The power spectra of the intensity shows that regions of high flux tend to have a smaller-scale granulation pattern. The size spectrum (power spectrum of the emer- gent intensity) of the simulated granulation is shown in Figure 17. There is a small 57 4; T ' ' T I ' j 1.5><10 . ......... 0 G , _ __ 200 G 4 3 1.0x104” — Q _ . E _ . :1 r I Z 3’ ‘ 5.0x10 ~— 4 0” 1 . . . . i —04 —0.2 0.0 02 04 AVG) )- A 1 1 I T I q 1.5x104~ ' ......... 0 G P _ _ 400 G ( 5 1.0x104— “i Q _ _ E _ _ 2 _ - k . 5.0x103" T r- -i 1 - _ J h J O 1 1 1 1 A A 1 —0.4 —0.2 0.0 0.2 0.4 (AI/<1) Figure 16. Distribution of the emergent intensity, summed over ten minute intervals. Large magnetic flux reduces both the brightness of typical granules and the darkness of the typical intergranular lanes. 58 increase in the power at larger scales for higher average magnetic field. This might be a result of the stabilizing effect of long-lived field structures at mesogranular scales, but the power increase is too small to be conclusive. At granulation scales and smaller there is no conclusive evidence for smaller—scale granulation. The power at these scales is reduced with increasing magnetic field strength for all well-resolved scale lengths. The intergranular lanes increase in width with increasing total flux, resulting in less overall area for granules. Higher resolution simulations will be necessary to test for features on the scale of a few hundred kilometers or less. Finally, magnetic field affects the lifetime of convective cells. The magnetic field, especially when there is a large magnetic flux which has a high filling factor in the intergranular lanes, acts a stabilization mechanism for the granulation. The field in these magnetic “walls” provides rigidity to the convective pattern. This can result in the convective flows inside of granules to be channeled upward with larger than typical velocity. The downflows that collar larger magnetic structures likewise stabilize these structures against fragmentation. The net effect is that the lifetime of convective cells with the presence of magnetic field tends to be longer than without magnetic field. This phenomenon is apparent in time sequence movies of the emergent intensity. The effect can also be seen in images of a horizontal slice of the emergent intensity taken over time. This type of image is shown in Figure 18 for each of the three simulation runs. There are more long-lived features evident in the runs containing magnetic field. Some of these features exist for several typical granule lifetimes. 59 0.0100_1111[ I T I 11111! I 1‘ Power 0.0010; 0.0001 1 l 1 l L 1 1 1 1 1 1 1 1 1 1 1 Wavenumber Figure 17. Granulation size spectrum given by the emergent intensity power, aver- aged over ten minute intervals. The wavenumber has units of Mm‘l. 60 Chapter 5 MAGNETIC FLUX TUBES Having examined the general properties of the magnetic field—convection interaction, it is time to turn to the magnetic structures themselves, namely “flux tubes” and “flux sheets”. The larger, dark structures, micropores, will be detailed in the next chapter. This chapter will focus on the structure and energy balance of “flux tubes/ sheets”. 5.1 Structure It was discussed earlier (Section 3.1) that the magnetic field is confined to intergran- ular regions by the convective flow. The shapes of the magnetic structures should then be influenced by the shapes of the intergranular lanes. Sheets of flux are formed in the lanes between granules, whereas structures that resemble ideal “flux tubes” are formed at the vertices of the intergranular lanes. These “flux tubes” are not the sym- metric, cylindrical, isolated structures used in many mathematical models, however. Their cross sections near the surface are polygonal, the shapes determined by the sizes of the surrounding granules. At least near the surface, they either extend into the connecting intergranular lanes forming star-shaped patterns, or connect to nearby “flux sheets”, depending on how much flux is threading the surface. The magnetic field strength is visualized in Figure 19. Below the surface the field is concentrated in 62 the narrow downdrafts. Approaching the surface (layer 15) the magnetic field spreads out as it is no longer confined by the gas pressure and horizontally diverging flow. The structure of magnetic features can be examined in more detail by taking vertical slices through the features, Figures 20 — 21. The slice cuts through three flux structures, one strong and two weaker. Figures 20 and 21 show the contours of the natural log of the density overlayed on filled contours of the magnetic field and current density, respectively. The density contour near 2 = 0 shows the granulation structure, being less dense in the hotter, diverging granules and more dense in the cooler, converging intergranular lanes. The strong “flux tube” shown here has its upper layers evacuated, and has its optical depth unity height depressed 300 km below the surface. The horizontal density is steep along the outer edges of the “flux tube”, where the horizontal magnetic field gradient is also steep; namely where there is a current sheet that collars the tube and is a boundary layer between the material being pushed outward from granules and the evacuated interior of the “flux tube”. The density is high at the outer boundary of “flux tubes” as well, since there exists a downflow. The weaker structures have started to evacuate their uppermost levels, but not below the surface. Since these regions have not fully undergone convective collapse, there is a pile up of material in the layers below the surface, as material from above is evacuated by a strong downflow. These “flux tubes” are also in the process of having their magnetic fields concentrated by the evacuation, surpassing equipartition strength, but still at or below the 1 kG level. The temperature structure of the “flux tubes” (Figures 22 and 23) give some insight into their appearance. The temperature contours are depressed inside of 63 20.00 10.00 Figure 19. Volume rendering of the magnetic field. Below the surface (layer 15) the field is collected into “flux tubes” by the diverging flow. Above the surface the field controls the flow and spreads out. 64 10 Mm Figure 20. Background: filled contours of the magnetic field in 250 G intervals. Overlay: contours of the natural logarithm of the density ranging from —2.0 (topmost) to 4.0 in 0.5 intervals. The T = 1 depth is shown as the thick line around 2 = 0 Mm. The established, strong “flux tube” in the center has been evacuated. The smaller flux concentrations on either side are in the process of being evacuated, starting above the surface and piling up plasma below the surface. 65 10 Mm Figure 21. Background: filled contours of the current density in 0.1 A/m2 intervals. Overlay: same as in Figure 20. A strong current sheet surrounds the “flux tube” in the center. 66 10 Mm “1W Figure 22. Background: filled contours of the magnetic field in 250 G intervals. Overlay: contours of the temperature from 4000 K to 16000 K in 1000 K intervals. The T = 1 depth is shown as the thick line around 2 = 0 Mm. The “flux tubes” are significantly cooler than their surroundings at a given geometrical depth. 67 Figure 23. Background: filled contours of the current density in 0.1 A/m2 intervals. Overlay: same as in Figure 22. 68 “flux tubes”. This is partially due to the fact that the magnetic field is expelled from granules to cooler intergranular lanes, and partially due to energy losses during convective collapse. At any given geometric height near the surface and below, the “flux tube” interior will be cooler than its surroundings, typically by a few thousand degrees, but occasionally more (the strong “flux tube” is about 6000 degrees cooler that its surroundings at z = .3 Mm). The temperature that would normally be near the surface in a typical downflow is pushed deeper inside the “flux tube” due to the initial evacuation. The temperature gradient below the unit optical depth depression is then much greater than the temperature gradient of the surrounding gas. This super-adiabatic gradient causes the interior of the collapsing “flux tube” to cool relative to the exterior medium and assists its evacuation. It is the how the temperature structure at optical depth unity inside compares to the temperature structure outside that determines whether a “flux tube” appears bright or dark, however. Small “flux tubes” with strong downflows are evacuated and cool, lowering their temperature and opacity. The evacuation is rapid, and since the opacity is highly temperature dependent (H‘ opacity is proportional to T8), optical depth gets depressed faster than the temperature. This leads to radiation from deeper, hotter gas than the surroundings, and the tube will appear bright. After the magnetic field in the tube has been compressed, flow will be suppressed and the evacuation will slow. If the tube is large enough to be thermally isolated from outside heating, it will continue to cool through vertical radiation. This allows the temperature inside the tube to continue to drop faster than the optical depth such that the interior of the “flux tube” is cooler than the surroundings at unit optical 69 depth. As in the case of the density, near the surface the horizontal temperature gradient is large in the current sheath that surrounds the “flux tube”. The gas dynamics must also be considered in any discussion on the structure of “flux tubes”. The velocity field for the example “flux tubes” is shown in Figures 24 and 25. The difference between a strong, dark “flux tube”, and a weaker, bright or neutral one is clear in the dynamical structure. The velocity inside the dark “flux tube” is greatly diminished compared to its surroundings. The tube is devoid of convective energy transfer, and thus cooled as seen in the temperature structure. Downflow in the tube is limited to the current sheet, and the only means of heating the tube interior is through lateral irradiation. In the weaker “flux tubes”, there is downflow throughout the interior, causing evacuation and lateral compression of the field. In smaller “flux tubes”, the current sheet extends through most of the tube, allowing downflow but not much cooling. Only when the field is compressed to to high enough strength that there is a back reaction from the Lorentz force will there be any plateau in magnetic field inside the tube, a suppression of the fluid flow and a chance for the tube to cool to the point where it appears dark. The flow pattern in weaker “flux tubes” and the resulting temperature and density structure imply that at least some bright points may be part of the evolutionary process associated with convective collapse. Lest it be thought that “flux tubes” are nice circular structures with downflows surrounding them, Figures 26 and 27 show images of the fluid velocity with contours of the magnetic field strength at the surface and 1.5 Mm below. The downflows occur in a few fast downdrafts on the periphery of the irregularly shaped “flux tube”. 70 Figure 24. Background: filled contours of the magnetic field in 250 G intervals. Overlay: velocity field. The T = 1 depth is shown as the thick line around 2 = 0 Mm. The flow is suppressed in the interior of the strong established flux tube in the center, but not the evacuating “flux tubes” on either side. 71 Figure 25. Background: filled contours of the current density in 0.1 A/m2 intervals. Overlay: same as in Figure 24. The downflows are confined to the current sheet at the boundary of the strong established “flux tube” in the center. 72 0.0 0.5 1.0 1.5 2.0 2.5 3.0 Mm surface, z=O Mm Figure 26. Image of vertical velocity (red and yellow down, blue and green up) with contours of the magnetic field strength at 0.5 kG intervals at the surface. 73 Figure 27. Image of vertical velocity (red and yellow down, blue and green up) with contours of the magnetic field strength at 0.5 kG intervals at the 1 Mm below the surface. 74 The difference in the flow field between weak and strong tubes is qualitatively similar to the kinematic/ dynamic regimes found in two—dimensional Boussinesq sim- ulations (Weiss 1964) and later in two-dimensional fully compressible simulations (Hurlburt & Toomre 1988). In the kinematic regime, the effect of the Lorentz force is small, so the magnetic field is swept to the convective cell boundaries and concen- trated into flux sheets. There is flow inside the magnetic regions, as seen in the weaker flux tubes in Figure 24. In the dynamic regime, the Lorentz force is strong enough to affect convective flows and suppress motion inside of magnetic regions. This behavior is analogous to the stronger tube in Figure 24. Quantitative comparisons are diffi- cult because the two—dimensional simulations use impenetrable boundary conditions resulting in convective rolls that are not seen in the current simulations. More recent simulations of this type use axisymmetric flux tubes to model pores and sunspots (Hurlburt & Rucklidge 2000). It was found that flows converge on the magnetic structure at the surface with downflows along the outer edge of the tube. The gen- eral properties of these tubes compare well with the results presented here, but there are differences in the details of the structure due differences in boundary conditions and input physics. Two—dimensional simulations of the dynamics of flux tubes with more realistic physics (Steiner et al. 1996, Steiner et al. 1998) are in closer agreement with the work presented in this section. However, there are several notable differences. The two-dimensional simulations show more swirls and rolls in the velocity field. This is most likely due to the choice of boundary conditions and the constraints of a two- dimensional geometry. The tubes also get buffeted more forcefully (have a greater 75 horizontal displacement), and have stronger downflows surrounding the tube in the two—dimensional simulations. This could be a result of higher resolution in the two- dimensional case, or differences in the structure due to geometry. The consequence of having different velocity and temperature profiles is that spectral line profiles will differ. ’ can be gained by examining how More insight into the structure of “flux tubes’ the relationship between the gas pressure and the magnetic pressure changes with depth. The local plasma beta is the ratio of the local gas pressure to the local magnetic pressure. Figure 28 shows the local B with depth, traced along several field lines. Field lines of higher magnetic field strength will have higher magnetic pressures and hence a lower B than low strength field lines. The stronger field regions are also less affected by gas motions and tend to have a tighter B distribution than weak field lines that are easily twisted around by convective motions. For strong—field regions, there is a change in the gradient of B around 2 = 0.5 Mm. Below this level, B changes by a factor of four from 2 = 1.0 Mm to z = 2.0 Mm, while the gradient is steeper above this level where B changes by more than two orders of magnitude over a megameter in height. The value of B decreases with height because the magnetic field scale height decreases more slowly than the gas pressure scale height. The structure in strong-field regions contributes to the sharp change around 2 = 0.5 Mm. Strong tubes are partially evacuated above the optical depth unity level, which is depressed 300 — 400 km below the z = 0 Mm surface. The lower density and cooler temperatures reduce the gas pressure and thus reduce the value of B. The local B is useful to determine the relative importance of the gas and 76 _O 5 I I I I IIIIHIVIIWI IV I\IIII I I I IIIITfi Ifi I TfTIII I I I I IIII o ' (A \ \ \ _‘ )_ .. f ‘ ‘ \ " )— \ —1 l. \ / _. _ ’m‘ q / / _." 0.5 P— , I I _ \ A " \ u / \ E 7 X? \ - \ _ Vh\ \ \ .. L \ \ \\ \ \ \ k \ \ x \— \ \ \\ \ ‘\\ '- \ \ ‘ \ V ._ \ \ \ I \ \ _ 1 5 _ \ \ I _\. . \ / ’— \ \ \ \ \/ _. \ \ |\ ’ \ \ \\ ' .. \ \ -‘ \ \ \ \ I \ \ 2.0 1 1 l 1 11111 1 14111111 M 1 1 111111 1 1 lullll 1 1u1111\11\ \ 1 1 11111 0.001 0.010 0.100 1.000 10.000 1000001000000 13 Figure 28. Plasma beta along magnetic field lines. Solid lines indicate strong field regions, dashed lines indicate weaker field regions. 77 magnetic pressures on the dynamics at a given point, but cannot be used to explain how the size of magnetic features changes with height; the pressures interior and exterior to the magnetic region must then be considered. The cross—sectional area of magnetic features increases above the surface. This can be seen in horizontal slices of the magnetic field at different depths (Figure 29). The magnetic region expands until the magnetic field weakens to the point where the exterior gas and dynamic pressures can balance the interior magnetic pressure. 5.2 Energy Budget The energy budget in a “flux tube” controls its temperature, pressure scale height, evacuation and appearance. The suppression of convective heat transport is believed to be the reason that larger magnetic structures are dark. The dominant mode of energy transport in magnetic features will therefore have consequences on the their evolution. This section investigates which terms in the energy equation (Equation 15) control the structure of “flux tubes” and their immediate surroundings. The example “flux tubes” from the previous section are utilized again to show the connection between heating, cooling and structure. Radiative transfer plays a vital role in the development of “flux tubes”. Figure 30 shows details of radiative heating/cooling for the same vertical slice through three “flux tubes” previously studied (Figures 20 — 25). The background color masks contributions between 3% of the minimum and 3% of the maximum. The top three panels show the net radiative heating/cooling. Contours of the temperature and magnetic field are given to mark the location of the magnetic structures. Except 78 0.516 Mm 2.010 Mm Z 2 1.521 Mm z=-0.010 Mm 2 0.500 Mm 0.992 Mm 1:— z: Figure 29. Horizontal slices of the magnetic field at 500 km intervals of depth. White/yellow indicates strong field regions, black/blue indicates weak field regions. 79 Net Radiative Heating/Cooling —0.5 _ I -1__ A °-°|.- JP‘““ ' —=—- _-‘*«-=a. m 2:?" g ‘1. if" V 0.5f ‘ Lo‘mfi...“ 1... .2 ... ..1. . ... 0 1 2 3 4 5 6 7 B 9 10 11 12 (Mm) Net Radiative Heating/Cooling with 1000 K Temperature Contours —0.5 \f« V - o . w 6 (him) Net Radiative Heating/Cooling with 250 G B-Field Contours -0.5 1 ’ ’ \J - E 0.0 . ~ As.w%;:-..ba ‘9qu m _~:‘ 5 0.5 < 1.0 .. .1”... .01. .. _ . 0 1 2 3 4 5 6 7 8 9 10 11 12 (Km) Radiative Heating/Cooling ,u.=:t l E "$.33 -“:‘“'-_"-' :tagum—Zw __—.~ ‘9": '2. "33.5.3521“;- 'n’ w 5 m ' g l 2 3 4 5 8 7 0 9 10 ll 12 (Jim) Radiative Heating/Cooling 11:1 0.5 -0.5 e — — i r- ' ——-. A o_o W13- .1_ 52‘ E P. «If ‘3‘ r J“ a—IIIr—I. " a.” I ,, v 0.5- ‘; 1n . . 0 1 2 3 4 5 8 7 8 9 10 11 12 (um) Radiative Heating/Cooling 11:1 0.05 _0.5 3 _. _ .N, . , E 0.0? ‘ . r.- < V 0.51 .l 1.0 ‘ 0 1 2 3 4 5 B 9 10 ll 12 8 (Dim) Figure 30. Radiative heating and cooling. Units are 1010 erg g'1 s“. The top three panels show the net radiative heating/cooling and its relation to the temperature and magnetic fields. The bottom three panels show the contribution of vertical, inclined and nearly horizontal rays. The interior of the strong established “flux tube” is nearly in radiative equilibrium, with cooling by vertical rays nearly balanced by horizontal heating from the side walls. 80 near their boundaries, “flux tubes” tend to be in radiative equilibrium with a balance between radiative heating and cooling in their interiors. The established tube in the center shows significant radiative energy transfer from material flowing down along the sides of the tube, cooling as it laterally heats the material near the edges of the tube. The other “flux tubes” in this slice show little net (above 3%) energy transfer through radiation, although there is some net radiative cooling as they are being evacuated. The bottom three panels of Figure 30 show the contributions to the total ra- diative energy transfer from rays at different angles. The first of these panels shows the results from vertical rays. The granulation structure is apparent. There is large cooling of gas in the granules at the level where the optical depth is unity. The radiation from this thin cooling layer illuminates and heats the gas above it in the photosphere. The overturning granule gas continues to cool radiatively as it is swept to and descends into the downdrafts. The two smaller “flux tubes” in this example have large downflows through their interiors and cooling that is similar to downdrafts devoid of any magnetic field. Near unit optical depth the large, central “flux tube” cools vertically and is heated horizontally from the sides. If we now consider rays traveling 60° from the vertical, the situation changes. The granules still Show cooling, but to a lesser extent, and there is much less heating of the gas above them. Even the interiors of the downflows show less overall cooling; this includes the weaker “flux tubes”. Cooling in the larger tube is now restricted to the surrounding downflow region. The more oblique radiation is able to heat the interior walls of the tube to a greater depth, however. For nearly horizontal rays the 81 radiative energy exchange is drastically different than the vertical case. Horizontal radiation heats rather than cools. There is a lateral heating above downdrafts as converging material near the surface radiatively cools. The established tube heats preferentially along the sides of the tubes. The level of heating/ cooling occurring inside the “flux tube” is an order of magnitude less than the radiative cooling from the top of granules. The question remains as to how large a role radiation plays in the overall energy budget of “flux tubes”. Contributions to the total heating and cooling are shown in Figures 31 and 32. The top two panels of Figures 31 and 32 give the net energy transfer in the example “flux tubes”. The contributions due to radiation are shown in the bot- tom three panels of Figure 31. There is a small degree of radiative cooling above downdrafts and weak “flux tubes” from overturning gas. This cooling is balanced by non-radiative heating and the regions above the intergranular lanes experience no net energy transfer, except near the level of optical depth unity. The established central tube is in radiative equilibrium in its interior, but still cools by vertical radiation from its walls. The net effects of lateral irradiation are minimal, as energy transferred from the exterior downflow is balanced by vertical radiative losses. A comparison between the net energy exchange and the energy transfer due to radiation reveals that other heating/ cooling mechanisms are at work. The contribu- tions from non-radiative processes are shown in the bottom three panels of Figure 32. The non-radiative contributions to the energy are from the advection of internal en- ergy (u-Ve) and work from expansion or compression (P(V-u)). Hot gas is advected 82 (Hm) (Hm) (Hm) (Mm) Net Heating/Cooling —0.5 a. A 0.0 _——~--:—., ..b—‘fil‘l mar—~1- .... 4" “ .3 1.0 5:. .._J1......L. I. 5.. . .. “.11., . 1 2 3 4 5 6 7 B 9 10 11 12 (lm) Net Heating/Cooling with 250 G B—Field Contours - ,- ' _— A 6 7 B 9 10 11 12 (Mm) Net Radiative Heating/Cooling . _ . 1 . _1 H1. fl, fir-21“ -l - 5“?“— _— 1-1 ‘I ‘ 5 ' -i 3 4 B 7 B 9 10 11 12 (Km) 0 Radiative Heating/Cooling 1.0 2112 0.5 —? ‘ ‘ .'_...,~-_Q ’ r'- '—-' 3’ .—-~ 1 "' —‘-‘ rah-at "i l“- -==- ...- "- ~ , —5 -1 3 4 7 8 9 10 11 12 0 (him) Radiative Heating/Cooling 0.5 2,112 0.0 Mfr—g. r; .._—. .-—-I Figure 31. radiative contribution (bottom three panels). Scale is the same as for Figure 30. Total net heating and cooling by all processes (top two panels) and the 83 (Mm) (Mm) (Mm) (Mm) (Mm) Figure 32. Net Heating/Cooling -0.5 0.0 _ ____ _ F:‘fi m7“ *m‘ _ __ 0.5 _ .’ fl 1.0 II.‘.: I. I r O 1 2 3 7 8 9 10 11 12 (Mm) Net Heating/Cooling with 250 G B-Field Contours 11 12 a (Mm) Total Non—Radiative Heating/Cooling -0.5 I. f 0.0 ~__ __ _ ...—" ‘ ' F -- a? . - 0.5 _d- ' ‘ .5 l N 1.0 $— '1” ,. V , ,K_‘ 0 l 2 3 4 5 8 7 a 9 10 11 12 (um) u-Ve Heating/Cooling _O.5 ; A-i,.“. 7‘ ;a. . .... , .i v‘ .7 , ...i.‘ z .z . F 0.0 x -. 0.5 _. {,— a' 1.0 . i 0 l 2 3 4 (Mm) P(V-u) Heating/Cooling Total net heating and cooling by all processes (top two panels) and non-radiative contribution (bottom three panels). Scale is the same as for Figure 30. 84 upwards in granules, heating these regions. However, the flows are divergent, and cool by expansion. These two effects are balanced except for the thin boundary layer at the tops of granules, where the steep temperature gradient make the advection process dominant and results in a net heating. The situation is reversed in the inter- granular lanes. Gas cools and is transported downward in converging flows. Neither of these mechanisms has much effect in strong “flux tubes” due to the suppression of velocities. The interiors of the weaker tubes show a net cooling as a result of down- flow. The outer walls of these tubes heat as a consequence of subsurface granular flows transferring energy into the surrounding downflow. 85 Chapter 6 MICROPORES The previous chapter detailed the structure of magnetic flux concentrations. Larger flux structures can become appreciably darker than their surroundings and are called micropores and pores. The question remains as to what conditions are necessary in the local environment for these structures to form. This chapter will discuss the formation and evolution of micropores, the extended, dark, intense magnetic structures found to develop in the simulation runs. A typical micropore has an area given by its intensity signature (I < 0.8 < I >) in the range 1 — 2 Mmz, although its magnetic size is slightly larger. The edges of magnetic structures are brighter where material flows down along the current sheet. The micropores have a fluxes of 1.5 — 3.0 X 1019 Mx. They do have observational analogs, and exhibit sizes, intensities and fluxes consistent with the smallest observed pores (Zirin & Wang 1992, Bonet, Sobotka, & Vazquez 1995, Soltau 1997, Siitterlin 1998, Leka & Skumanich 1998, Keil et al. 1999). The micropores in the simulations occupy between 1 — 3% of the surface, similar to the 2% coverage found for pores and micropores in observations (Berger et al. 1995). 86 6.1 Formation The defining characteristic of micropores is that they are so much darker than other magnetic structures (on the scale of granulation). To what do they owe this distinc- tion? Some clue lies in the formation process. Figure 33 shows a time sequence of relative intensity versus magnetic field strength during the formation of a micropore. The first thing to note is that the length of the sequence is only four minutes; the formation is necessarily shorter than the convective time scale to have any reasonable chance at success. As time progresses, a tail breaks loose in the distribution of inten- sity and magnetic field and becomes distinct from those points located elsewhere in the intergranular lanes. Intermediate strength field is darkest first, with higher field strength becoming darker as the tail gets swept (magnetic field gets concentrated) to larger values. The mechanism responsible for this process cannot be flux expulsion or convective collapse, or at least it cannot be solely these mechanisms at work. Even though it is apparent that the field is being compressed, neither of these can explain why lower strength field gets darker first. Also, if either of these mechanisms was responsible, we would expect to see dark magnetic field regions with a much higher frequency than is seen. The answer lies in a convective mechanism. Another time series covering the formation of a micropore is shown in Figure 34. The left panel shows a mask including high field strength points (B > 1.5 kG) and low emergent intensity points (I < 0.8 < I >). A high-field, dark ring forms first, with the interior subsequently decreasing in intensity and gaining in field strength. The middle panel shows the emergent 87 t=56.5 min t=55.5 min Y r V“) V“) 1/(1) 0.0 0.5 1.0 1.5 3 (k6) t=57.5 min V") 1/(l) l/(l) t=59.0 min VV’VYY 'fi w—fi l/(l) V“) V“) Figure 33. Correlation of intensity with magnetic field strength as a function of time when a micropore forms. The lowest intensity first occurs at intermediate field strengths. Then the field strength in the micropore gradually increases to a maximum. 88 . t: 5 i 59.0 min 57.0 min Q :55 5 mm a o 2 f 53.0 min 2 (In) Magnetic Field 1 Emergent Intensity N (In) Figure 34. Images of magnetic field (right panels), intensity (center panels) and mask showing only low intensity, strong field locations (left panels). Micropores form when magnetic flux from surrounding intergranular lanes is swept into locations where a small granule disappears. 89 intensity for the same times. A small granule is in the process of disappearing, leaving behind a dark hole. The magnetic field is shown in the right panel. The disappearing granule is surrounded by magnetic field that gets pushed into the space the granule leaves behind. There is a gradient in the magnetic field strength at the interface between a granule and intergranular lane, with weaker field closest to the granule and increasing into the lane. That explains why weaker field is darker first, it enters the dark region first. It is then compressed by the gas pressure gradient to high field strength. Magnetic features are generally thought to be formed when a “flux tube” emerges at the surface after having risen from below. Due to the initial conditions, we do not simulate this process. The merger process described here is an alternative method of magnetic feature formation. This type of formation process is much like that seen in observations of network bright points, which has been named formation by granule compression (Muller & Roudier 1992). The network bright points appear to form in large, dark spaces between granules as the magnetic field gets compressed by converging granules. The formation time was found to be rapid, occurring in approximately four minutes. Further observations confirmed the process in the network (Roudier et al. 1994) and also in plage regions (Berger & Title 1996). There are currently no observational studies detailing the formation of micropores, due to the difficulty in achieving the spatial and temporal resolution necessary to follow their evolution. There are no reliable proxies to find the micropores as there is with bright points. The micropore formation process may be just the high flux extension of the process observed in bright points. A large magnetic filling factor in the intergranular 90 lanes is a necessary prerequisite. Large amounts of high—strength field surrounding the disappearing granule ensures that the inflow of heat by convective means is reduced and that the magnetic feature will have a large enough radius that lateral radiation cannot heat the interior. This reduction in heat flow and the resultant cooling is similar to the process thought to occur in sunspots. Collapsing granules also appear to be sources of transient acoustic waves (Skartlien, Stein, & Nordlund 2000). If granule compression is a fundamental means to form small to mid-sized magnetic structures, then the formation process itself may be important in transferring heat to the chromosphere. 6.2 Evolution Sunspots are observed to have lifetimes ranging from several hours up to months, while pores are observed to have lifetimes ranging from under an hour up to days. Since micropores are smaller structures, it is reasonable to assume that their lifetimes are shorter than pores. It is difficult to track pores in time observationally, so the best recourse is to study them through simulations. This section examines the properties of micropores over the course of their lifetimes. The evolution of micropores can be studied by following in time the properties of the central point of the collapsing granule (and the subsequent magnetic feature) at different depths in the convection zone (Figures 35-37). The granule collapses because its decreasing heat flux cannot balance the radia- tive losses at the surface, and the gas cools and sinks. The energy flux decreases as a result of a diminished upflow. This deceleration starts at t = —3 min at the surface, 91 11111 I I I I 77 I I I I I I I I I Y I I Y I rfi I 77 I I I I "l Ifi 17" I T I 4’— ‘ ‘ ,_ / C * z=00 Mm / a h ‘— ‘ I. —I b- I —: ~ ......... 2:0 5 Mm / -. F- - .. / -4 ~ z=10 Mm / ~ — — — — ",/ u p— a‘ — 3_ u / _ >— / — i— -I >— / — A '— q P— " / —‘ u - . . , _ .3 _. " ’ ‘ n— ' / —: v 2>—— ‘. / — I— , \ —1 m r— ' // -1 >— -1 r- -4 p— —‘ i- d _ _ F— — .— 1— n— —: r— .1 F- — i.- d r — — — —- ~ — —— r;- — ... b- —4 )— ...... -— g— ......... -4 u- -I O TTJ I 11111111 LJ 111111111 l ......... l 111111111 l 1 +1 fT vvvvvvvvv 1 VVVVVVVVV l rrrrrrrrr I VVVVVVVV l _ F - II— I‘ - q — u- .— n- u— .- 1 ._ .— .— ._ — ,_ __.. —. - . - — I“ h ‘1 A _ — — _ u ' O.._ _ m _ — r- q r— —4 — d _ - .— —I .54 .. _ v " .. u— - - _ _ D K _ )- —1 — _ p— —I y— — r— — >— ....... — I.- _A -2h-U— - — _ p- d —‘ ... _. _4 —A A l I A l L J. L A L l A A A l AAAAAAAAA l LLLLLLLLL l A A A l L L A A A J_‘_J_ . bIn T—(T) (10’ K) )- c- -4 P— h—n-t P- — L L A 4L 1 A L L L l l L l A A A A L4 L4 LLLLLLLL 1 J_ LLJ l L -10 0 10 20 30 Time (min) Figure 35. Evolution of the central point of a micropore at different depths in the convection zone. The flow reverses direction at the surface at t = 0 minutes. Top panel: magnetic field. Middle panel: vertical velocity. Bottom panel: temperature fluctuations. See text for details. 92 77777 IV V Y I V V V I 1' T—Y [V V V V "V V 77 1‘77 Y—Y V V Y Y T Ifrfi" 1" '1 r- —< L. a Z i A i. .4 .— -4 is o: : )—- o .,,¢ an E“ 1/ N. 'r E O u I v-1 _ _. A : W.“ /: s t v : l C ' 2 <1 : ‘-: +— ‘.-4 -2_— "E i— d r- -4 ” 1 ........ 1 ......... 1 14" I .................. 1 j J —1 A q 'u" 0 \ 4 \ s. \ — “U - \ - 2 2:0.0 Mm \ ._ \ V b ......... z=0.5 Mm ‘~-.\ ‘ A _ o. 4”— .___z=l.OMm \ F v .. \ _ 7‘ I \ \ .. D.- r \ _ '1 \ .. \\ 7 _6.— \ x— . 1 1 ......... 1 ......... 1 ......... 1 ......... 1 . 1 . -10 O 10 20 30 Time (min) Figure 36. Evolution of the central point of a micropore at different depths in the convection zone. The flow reverses direction at the surface at t = 0 minuntes. Top panel: density fluctuations. Bottom panel: gas pressure fluctuations. See text for details. 93 VVVVV IV V V V V V V r1 I V V V‘rV Y V V V I V V V V U f7 V V ‘jfifiYY Y [TY 1'—‘ 2:— , 1 : 2:00 Mm a . : >— T -1 : """"" an : : _ _. _ 3p + ac ., 1*— _. A Z 3 ’1‘ ~ . .1 m : r: E I / 3 O 0”— ------------ ; -------- l ------------------ _ V r _1 S : ‘ ’ \ 3 w ~ \/_1 I \ / Z -l:' ‘ ’ 1 : \ I Z : \_ / : "1 41 1 11111111 L1 11 1 111111111 111w 111111 1111‘ ‘VYV I VVVVVVVV rr j—fiv TV ‘ VVVVVVVVV I T‘Vfi VVVVVV ‘l ” z—0.5 Mm ‘ 0.05~ _ b I‘ .1 :— !--\/ q _/ ‘ H . ... f ‘ §‘ 0.00 ,, .3 . 1.». .-, a. 1 1 .— m \ I I '. /\ ~\ I i: E \\A I \ J‘ I \ I \ /\ '1‘— o r \/ “ \ \ x . -0.05— \ ’ \' \ — o _ f / g /\/“\ /1 v _ I ‘ _4 —0.10- 1 ' / — \ ‘1 - / —O.15— ..... L AAAAAAAAA l 111111111 1+ 11111111 1 111111111 1; A . . . ., ......... , ......... I ......... 1 ......... 1 . . . ’ z=l.0 Mm ‘ >- I —1 0.05 \ _ l / ‘ I \ 5" — / I l U} ._ / I" E 0.00 ‘7'" ;./ ................. / 0 it x \ I f \ ~ - b “N I“ \/ \\_ \ " d L- \/ \ \ \ \ — v __ I \/ l \ _ 005 “ ‘ \ ’\ H o F- l‘ / \I _‘ >— \ l J ..1 V F #111111 11E1111 111111111 1 111111111 111111.11111 -lO 0 10 20 30 Time (min) Figure 37. Evolution of the central point of a micropore. The flow reverses direction at the surface at t = 0 minutes. Comparison at different depths in the convection zone of the total acceleration, aT, the acceleration due to the Lorentz force, a3, and the sum of the accelerations due to the pressure gradient force, ap, and the buoyancy force, ac. See text for details. 94 and several minutes earlier at the deeper levels. The deceleration of the upflow leads to rapid cooling with the cool gas piling up on top of the still rising gas. The gas pressure starts to decrease as the temperature falls, allowing the magnetic field to get compressed and increase its field strength. When the flow reverses at the surface (t = 0 minutes) the cool gas is evacuated, resulting in a further drop in the gas pres- sure and increase in magnetic field strength. The compression process takes about six minutes from the time the granule starts sinking until maximum field strength is achieved. The downflow is suppressed as the field strength increases, and remains small until the micropore disappears 15 minutes later. This halts further evacuation, so the gas pressure and density stay nearly constant during this time period as well. Prior to the granule sinking at the surface, there is a decrease in the upflow at the deeper levels. A positive pressure fluctuation at z = 1.0 Mm causes an upward acceleration due to the steepened pressure gradient, and the small downflow prior to t = 0 minutes becomes an upflow by the time the flow reverses at the surface. A similar effect occurs at z = 0.5 Mm a few minutes later. The gas that is evacuated near the surface can be seen later as density enhancements in deeper layers. As the material piles up the buoyancy force increases and starts to counter the upward acceleration due to the pressure gradient. The gas is subsequently evacuated when the subsurface upflows become downflows. These layers continue to cool, decreasing the pressure and increasing the magnetic field strength even after the intensity sig- nature of the micropore has disappeared from the surface. This particular micropore brightens as the result of a deformation that allows it to get heated by an inflow of gas. The density increases rapidly and the temperature and pressure increase, 95 causing the magnetic field to expand and decrease in strength. A downflow starts as the Lorentz force increases to prevent the decay of the magnetic field, but ulti- mately is not enough. A weaker magnetic region continues to exist, but the increased temperature results in the feature becoming indistinguishable from the intergranular lanes. The evolution of the magnetic structures is greatly affected by convective mo- tions. Flux tubes are continually buffeted about by the expansion of granules and flows along the intergranular lanes. Even the more rigid micropores are no match for granular flow. The shape of a micropore is constantly changing, from being elongated to being squeezed into nearby intergranular lanes, the “flux tube’s” life is a dynamic one. Magnetic structures undergo many mergers and splittings during the course of their lifetime (Figure 39 and Figure 40). One other facet of a “flux tube’s” evolution is apparent in Figure 39. Both small diameter and large diameter tubes can lose and regain their intensity signatures of the course of their lifetimes. Flux structures can become elongated to the point where they fragment, or are squeezed to such a narrow width that heating from the ambient environment becomes important. This does not invalidate the use of photometric proxies to locate magnetic field structures. It does however, make it more difficult to draw conclusions about the lifetimes of magnetic features based solely on those proxies. This fact can be seen clearly in three-dimensional renderings of a magnetic flux region (Figure 39). The left figure is the magnetic field strength at the surface (red is high, blue is low) and the right figure is the corresponding emergent intensity (red is darker, blue is brighter). The total length of time represented is three hours. 96 Figure 38. Time — horizontal slice through the center of a micropore showing an image of the emergent intensity and magnetic field contours. 97 Figure 39. Two dimensional horizontal — time visualization of surface magnetic field (left) and emergent intensity (right) during formation of a micropore. Vertical scale is time in half minute intervals. Size of the box is approximately 2x2 Mm. 98 Minutes Figure 40. Time — horizontal slice through the center of a micropore showing an image of the magnetic field strength for the same region as Figure 38. 99 Not only are there many mergers in the region, but the pre-darkened state of a micropore is evident around snapshot 225. It is also apparent that not only does the brightness of the field change often, so does the magnetic field strength. Therefore, magnetic structures can lose their intensity signature on a convective time scale. Only those flux concentrations that reside in stable downdrafts (i.e., those at the vertices of mesogranules with deeply penetrating downflows) really have any chance to survive more than a few convective turnover times without being significantly altered. 100 Chapter 7 CONCLUSIONS 7.1 Summary A series of numerical simulations was used to investigate magnetic fields in a convec- tively unstable atmosphere in order to get a better understanding of the interactions and features seen in solar observations. The general properties of the magnetic field agree well with observations. Magnetic field is transported from granules to inter- granular lanes by convective flow. The field gathers in the intergranular lanes and is compressed. The field compression is further amplified by a convective instability, and the magnetic regions cool and become evacuated. Magnetic features form pref- erentially at the vertices of intergranular lanes, and in particular at the vertices of mesogranular lanes. Convection is responsible for the compression and formation of magnetic features, but magnetic fields, in turn, affect fluid motions. Low to moderate mean field strength regions tend to have granules of smaller size, while moderate to high mean field strength regions show clumping of clusters of field—free granules. Granule boundaries get more rounded, and more low-brightness granules appear as the flux levels increase. The magnetic sheets in the intergranular lanes provides rigidity to the convective cell 101 pattern, and granules have longer average lifetimes in the presence of field. The interior of “flux tubes” and “flux sheets” is partially evacuated, lowering the opacity, and causing a depression of the optical depth scale. The largest structures have a depression of the unit optical depth level of approximately 300 km. The interiors also cool since the temperature gradient is steeper than in the environment outside the tube. The “flux tubes” are collared by a downflow, which tends to be brighter than the average intergranular lane intensity. Inside of an established “flux tube”, the flow is suppressed by the Lorentz force. This greatly reduces the amount of convective energy transfer into the tube; radiative energy transport thus becomes important. The interior of the “flux tubes” lose energy by vertical radiation, but the sides of the tube heat as the gas in the downflow outside of the tube radiates laterally inward. Dark micropores form at the vertices of mesogranules where a small granule submerges. Magnetic field that surrounded the granule is then pushed into the dark region left behind by the sinking granule as neighboring granules expand toward the hole. The micropores are large enough to survive for more than a few convective turnover time scales. Even though these structures survive, they are constantly be- ing deformed and squeezed into the intergranular lanes by changes in the granulation pattern. Over the course of their lifetimes, larger magnetic structures undergo many splittings and mergers. This dynamic evolution can cause the magnetic features to lose their intensity signature while still being a coherent magnetic structure, suggest- ing that their evolution cannot be traced (at least in its entirety) by these signatures. The Solar-B mission planned for launch in 2005 will provide photospheric mag- 102 netic field data with enough spatial resolution to study small—scale magnetic features. The satellite will contain an optical vector magnetograph to study photospheric mag- netic fields and coordinated X—ray/XUV imaging telescopes to get simultaneous data in the corona. In addition, Solar—B will be able to accurately measure all three com- ponents of the photospheric magnetic field. Comparisons of the simulations to the Solar—B data will be a valuable aid in constraining future calculations, and determin- ing how the input physics can be refined. 7.2 Future Work This work has just scratched the surface of possible knowledge that can be gleaned from this type of numerical simulation. As with observational studies, there is always the desire to get better spatial resolution. There always seems to be substructure when we gain the ability to look at smaller scales. The computational cost to increase the resolution with the current code is large, and work is already being done to better utilize massively parallel computing resources. The scope of investigations can also be broadened. Only one magnetic field configuration was studied here, and a mainly unipolar field is not generally applicable to the Sun as a whole. A more realistic scenario might be to include horizontal field that is advected in through the bottom of the computational domain. Another option would be to feed in “flux tubes” formed in deep convection simulations to determine if there are any important differences in structure as the tubes emerge at the surface. More work can be done to make comparisons to observations. Two dimensional simulations have been used to develop diagnostics for observations (e.g., Steiner et 103 al. 1998). Now we have data to make Stokes profiles for three dimensional magnetic features. It would also be useful to do the radiative calculations in the G-band, where observations have shown that magnetic bright point have a greater contrast compared to the continuum (Muller & Roudier 1984, Berger et al. 1995). 104 APPENDIX A Derivative Schemes This appendix details the formalism and accuracy of the derivative schemes used in the code. Prior versions of the code used tridiagonal finite difference schemes, so the goal of the update was to keep a tridiagonal scheme, but to increase the formal accuracy by increasing the number of stencil points. A. 1 First Derivatives A.1.1 Horizontal Direction: Uniform, Periodic Grid Compact finite difference schemes can be used to approach the global dependence of spectral methods without being as computationally intensive (Collatz 1966, Lele 1992). For a uniform grid with grid spacing, h, a tridiagonal scheme with a five-point stencil is given by I I I a b 0‘ i+l + fi + a i—l = fifffil — fi-1)+ 217241242 — fi—2)- (66) The values of a, a, and b can be determined by matching the coefficients of the Taylor expansions of the left and right sides of Equation 66. The first unmatched coefficient then gives a measure of the formal truncation error of the scheme. Due 105 to the symmetric nature of Equation 66, no information can be gained from the odd—order accurate equations. Matching coefficients for the even—order terms gives 2nd Order (ff): 1 + 20 = a + b (67) m , 20 _ 1 2b MhOMm(fl) §T_fi(m+2) ma 2a 1 6th Order (12(5)); I = ‘5‘! (a + 2%). (69) The old derivative scheme uses the fourth-order accurate, one-parameter (a) family. Solving Equation 67 and Equation 68 gives the dependence of a and b on a, 2 1 The truncation error for this family is given by the coefficients of the ff”) terms, but since these coefficients vanish, the formal error is given by ezéwa—nmfifl (n) The value of a used in the code is chosen to give a cubic spline scheme with the spline condition that second derivatives are continuous at the end points of each interval. Symmetry adds an extra order of accuracy, making the scheme equivalent to the fourth-order classical Padé scheme with a = 1/4. With this scheme, the other coefficients are 3 1 a = 3, b = 0, C = —g?h4fi(5). (72) 106 Note that this scheme only requires a three—point stencil. The new first derivative scheme uses the sixth-order accurate scheme given by solving Equations 67—69. Solving this system of equations leads to the coefficients 6 = —h6f.‘7’. (73) l. The accuracy of these schemes is analyzed by investigating the amplitude errors relative to an analytical solution. A test function of the form f(.r) = expfikiv) (74) is used where k = 27r/A is the wavenumber and /\ is the wavelength of the oscillation. The results are shown in Figure 41 for a sampling of wavenumbers. The maximum fractional amplitude error, 6, is given by _ ffd 7 v (75) 6=max where the subscript fd represents the finite difference derivative. At small wavenum- bers, the sixth-order scheme is two orders of magnitude more accurate than the fourth order scheme, as might be expected. The sixth-order scheme is also more accurate at shorter wavelengths. The exception is for k = 71’, where the oscillation is a saw- tooth function, and both schemes give zero for the derivative. Since both schemes are central difference schemes, the errors are dispersive in nature, not dissipative. The symmetry of the schemes can be investigated by taking the horizontal di- 107 Y I Y I I Y I T T r T Y T I I I I I I l I I I I I I I I Y I Y Y 100 — x —« l— + -1 + x e 10—2 — + BK _ LU + Q) - —4 8 + ”i '71 10‘4 e + 3" — E + 6 ” + 3" ‘ s 6 1 1.. f, 10 — — 8 + L‘“ ~ X + Compoct 4th — :E, —8 9K 3K Compact 6th .E 10 — .1 X 0 3K 2 a .. 10‘1O — _. 3K '- 1 1 1 1 l 1 1 1 4 l 1 1 1 1 I 1 1 1 1 l 1 1 1 1 l 1 1 1 1 l 1 T 0.0 0.5 1.0 1.5 2.0 2.5 3.0 Wovenumber Figure 41. Amplitude errors relative to the analytic solution for first derivative schemes. Symbols represent sampled wavenumbers. 108 vergence of a two-dimensional Gaussian function and comparing the results to the analytical solution. The error, 5, is given by a (9 a a 5: [31: + 5;] “51719) — “8—5)“ + (fig—>111] “$131), (76) where f(:r,y) is the Gaussian profile shown in Figure 42. Figure 43 and Figure 44 show the profiles ofé for the old and new schemes, respectively. The asymmetries of these schemes can been seen in the contours of 5. As seen in Figures 45 and 46, both schemes lack circular symmetry. A.1.2 Vertical Direction: Non-Uniform, Non-Periodic Grid Since the vertical direction uses a non-uniform grid, we must generalize the previ— ous derivative formalism for the interior regions and derive schemes for the top and bottom boundary zones. Interior First we will consider the interior domain. A tridiagonal scheme with a three—point stencil is given by a 1’11 + f.-’ + fif.’-1 = af1+1 + bf1+ cf.-_1. (77) It will be convenient to define the following grid spacings: h+ : ZH-l — Zia h- : Zi _ Zi-l- (78) 109 1.0- Figure 42. Gaussian profile used for horizontal error analysis. 110 Figure 43. Error profile of the horizontal divergence of a Gaussian function for the fourth-order first derivative scheme. 111 5x10’ ..—......0—.—--—.—.———.-————--.-.- -—---——-—---o--—---o--- .......... -----—-—---—-C’---—------—--—--¢ o—-——----—--—--——-------- ....... .---—-------------—---—--— ....... -————----—--v o-—-——----—--- ----------.-- -—--------- ‘-- --.------ ------—---- -- ------—-——— ---------— ---------’- -n---———--' cu---— -— o...—-. .-.—9- ------ --———- ---—- ----o ---‘- - ---- ---“.. .. ~ ‘ \\ ‘0‘“ \\“‘ 1 \x s .1 “\\‘1 ,filr I. O‘ ’l' ‘ I 11“ I ‘ \ ‘ a. \\\ \\“ ‘: 6 ”at I ---- ----- ' ...-C. ..—C‘ ----u ’-00 ..- -‘CD- -----. ‘.~ -—--- on--- - ------‘OD- cc--- -- unnu- ------v- ...-------c--- ---—--—-—-—-—--- --‘—-’----- ----—---- --—----- ---- -..—-— ..-- ------- --—-——-- ...... c-‘----- ........... Figure 44 . ' Error SlXth- Profile of th . order first derivatiVe SchgrfllOFlZontal divergence of a G 6' aUSSian funC . tion for th e 112 T .llIIIIIITTITIIIIlllllIIIIIIIIIIIIIIIIIIIIIIIIIIIYIIIIIIIIIIII. 60 50 4O llllllllilllllllllilllyllIll 3O 20 III]IIIIIIIITT—FTTIIIIIIIIIIIIIITIIIIIIIIIIIIIIIIIII i 10 111111111lllllllllllllllJlJJl IIIIIIIII O"lllllllliLlLllllllliLLllllllilllllllllilllLllllllllllllJlllIfl 1O 20 30 4O 50 6O 0 Figure 45. Contour of the error profile of the horizontal divergence of a Gaussian function for the fourth-order first derivative scheme. Dotted lines are shown as aids in visualizing asymmetry. 113 . IIIIITIIIIIIIIIIIITIIIIIIIITIIIIIIITIIIIIIIIIIIIIITIIIITTTII, 50 a I 50 4O lLlJLIJliJlllILlllilllll]llli 3O 20 IIIIIIIIIIllIIIIIIllllllIlllllllllTIITITrTTIITrllll lllllllllllllllllllllllllllll lllllllll O ’llllllllllllllllLllILllllllllillllllilliLlJLllllIilIlllllll[I'— 10 20 30 4O 50 6O 0 Figure 46. Contour of the error profile of the horizontal divergence of a Gaussian function for the sixth-order first derivative scheme. Dotted lines are shown as aids in visualizing asymmetry. 114 The system of equations defining a given order scheme are found by matching the the coeflicients of the Taylor expansion of Equation 77, 0th Order (f1): 0 = a + b + c (79) 1st Order ()7); a+0+1= h+a —h_c (80) 2nd Order (fg'): 1210 — 11-5 = 21701111 + hie) (81) 3rd Order (f,.’”): hia + 123/3 = 3‘; (Ilia -— h3_c) . (82) The current scheme uses the third-order accurate, one parameter family given by the solution of Equations 79—82. The parameter is chosen to give a cubic spline scheme such that the second derivatives are continuous across adjacent zones. This solution gives 1 h_ 3011—11-) = — , b : ______, " 2h++h_ 2 1.111- _ 1 h+ 3 121 ”i _ 2h++h_’ C ‘ _§h_(h++h_)’ (83) a — § h‘ — 1 h h h h (4) 2h+(h++h—)’ 6 ‘ '2.41+ ‘( +" —)f.-. The truncation error of the scheme, 6, obviously depends on the details of the grid. Top Boundary We are restricted to using one-sided derivative formulations at the vertical boundaries of the computational domain. At the top boundary we have a fiducial layer, so the derivative scheme does not need to be highly accurate. Therefore we use an explicit 115 first-order scheme of the form fi =af1+bf2 (84) The coefficients of this scheme are a : —b = —1/(22 — 21). Bottom Boundary At the bottom boundary we use the cubic spline condition, continuity of third deriva- tives. The third derivative in the last two zones can be expressed as 1,1,2: : afnz + bfnz-l + Cfnz—Z + afTIIZ + (Bffiz-l (85) f.’.’:_1 -- afm. +Bfm_1+éfm_2 +c‘rf,’.. +Bfi..-1, (86) where n2 is the number of points in the vertical direction. Then, for "’ = ’” nz nz-l’ (0 - 6tlf'z + (5 - Elfin—1 = (fl - a)fnz + (5 - blan + (5 - C)fnz—2- (87) The solution to this equation is determined by matching the coefficients of the Taylor expansion of Equation 87. This yields l a-“ z h_h2_’ 13 b _ 3h_—2h2_ - 1 _ _ hZ(h_ — h2_)2’ 3"fl = h_(h2_—h_), _ _ h_ (88) &_a = h_+2h2_ “—6 ’ —h§_(h_—h2_)2’ hih§_ ’ where h- = znz — znz_1 and h2_ = znz — znz_2. 116 A.2 Second Derivatives A.2.l Horizontal Direction: Uniform, Periodic Grid In an effort to reduce the computational load of the calculations, the old second derivative scheme is an explicit forward difference scheme. The scheme takes the derivative of the spline derivative scheme detailed in A.1.1 and is of the form [3 f1” = 5; (L11 — f.-) + 52ft. + Eff, (89) where h is the uniform grid spacing. The system of equations governing the accuracy of this scheme are determined by matching the coefficients of the orders of the Taylor expansion of Equation 89: 0th Order (ff): 0 = a + a + B (90) lst Order (f,-”): 2 = a + 20 (91) 2nd Order (f,-'”): 0 = a + 30. (92) Solution of Equations 90—92 gives 2 2 (4) a=6,a=—2,fl=—4,e=—-4—'hf, , (93) where c is the formal truncation error of the scheme. The new horizontal second derivative scheme is a compact tridiagonal scheme similar to the compact first derivative scheme. The increase in computing speed has 117 made this feasible, without a large increase in run time. The new scheme utilizes a five-point stencil of the form II II II a b a 1.1+ f, + a .-. = h—,(f.+1 — 2f.- + f:_1)+ 4712—0.” — 212+ 1-2). (94) Matching coeflicients of the Taylor expansions gives 2nd Order (ff’): 1 + 20 = a + b (95) 4th Order (fig) 2 — i (a + 22b) (96) ’ ' 2! — 4! a 1 6th Order (f,(6)): 4—! = 6—1 (a + 24b) . (97) The new scheme uses the sixth-order accurate scheme determined by the solution of Equation 95 through Equation 97. The parameters of this scheme are 2 12 3 184 a 11’ a 11’ 11’ 6 11-8! f’ (98) The sixth-order compact scheme is a centered scheme while the second-order explicit scheme is a forward difference scheme. This has some effect on the symmetry of the schemes, which will be seen in the discussion on errors. The errors will be analyzed in the same manner as for the first derivative schemes. The accuracy of the second derivatives is determined by investigating the ampli- tude errors of the schemes as defined in Equation 75. The values of 6 for various wavenumbers are plotted in Figure 47. The new scheme is much more accurate at 118 I I I V I j 7 f Tj T I T Y I T Y Y 1 I l' I T T T 1 I I I I 7 Y 100~ __ + + 3‘ + + c + 9" 8 10‘2— + — 11‘] ++ ”i Q) >— + —. -: —4 Q 10 3K '* E < B ’ x ‘ 8 6 45, 10 — 9K — 9 LL )- 316 + Explicit 2nd ‘ E 8 3K 3K Compoct 6th .E 10- — ‘— X E _ x 1 10‘10— _. x *— l l l l l l l 1 l l l l l 1 l l l 1 l L 4L 1 1 1 J. l l l l l 1 J¥d 0.0 0.5 1.0 1.5 2.0 2.5 3.0 Wovenumber Figure 47. Amplitude errors relative to analytic solution for second derivative schemes. Symbols represent sampled wavenumbers. 119 small wavenumbers, and results in improved values for the current density, J. The schemes converge at k = 7r, but unlike the first derivative schemes that pass two-zone oscillations, the second derivative schemes have some damping at this wavelength. The symmetry of the second derivative schemes can be seen by looking at the errors in the horizontal Laplacian of a Gaussian function relative to the analytic solution. The Gaussian profile is the same used for the error analysis of the first derivative schemes, and is shown in Figure 42. The error is given by (‘32 82 02 02 . 1——11(—)(——>1 where f(.r,y) is the Gaussian profile. The profiles of 6 for the second derivative schemes are shown in Figure 48 and Figure 49. The contours of these profiles are shown in Figure 50 and Figure 51. Both schemes are symmetric about the horizontal directions and the 45° lines. The explicit scheme is not symmetric about the center of the Gaussian because it is a forward difference scheme, but rather it is symmetric about a point shifted by half a zone in both the horizontal directions. A.2.2 Vertical Direction: Non-Uniform, Non-Periodic Grid Interior The second derivative scheme in the vertical direction is of the form 1” = affi+1— fi) + b(f: - fi-1)+ 0(ff+1+ 2f:) + “ff—1 + Zfil- (100) 120 Figure 48. Error profile of the horizontal Laplacian of a Gaussian function for the second-order second derivative scheme. 121 —' '10' 333 ||I IIIIIIII I" 11111 C“.I“.1‘ "Wit """""" "11")" """""" "I.‘ ==== 11"11‘". 111111111 11111111111 ''''''''''' ............ '1'". iiii I" .| '1" ...-I ---—¢ ---—- ..-—-— 731;: \ ..| 1‘. I .| 'I I I ...- f— -_ .. ‘-- cu--- ““_-- ——--- ._---.:':.‘.-.'.'.-. Ilununcncn-nu-unc. ..-—--- ...-.--ucuuu --‘------- ---—--- -..— ..... - -------- ‘-—‘-’-— —-——--- ---- — ’---——----—-—------— ---—------uuoucnc- —-—--¢--—-——---——-- ------—--C‘—-‘-‘—‘------- ---‘—---—------—-----—-‘ -----——---‘---——-——-—--- ........ -----------—-—----——-—---——o— ———————— gun-------u-a—----_-c—nuco----—n----_—— —-—¢-——--———--—-— ....................... A0 50 Figure 49 . ° Errol- Slxt - Profile of t . h order second derivative :fhhorlzontal LaPleaCian of a G eme, aUSSian quC - 122 —< .IIIIIIIIIIlITTTTIIIIITITTWIIIIIIIIIIIIIIIIllllllllllllllllll bl 60 50 4O 1111111111111111111lllllllllLl 3O 20 [TIIIIIIIIIIFIIIIIIII[IlllIIITIIFTITTIIIITTTTTITTII 'lllllllll111111111111111111111 10_ OTIIIIIIJIII11111lllllllllLJlJLliILlLilllliJlllllllllllllllllli‘l 0 10 20 30 40 50 60 Figure 50. Contour of the error profile of the horizontal Laplacian of a Gaussian function for the second-order second derivative scheme. Dotted lines are shown as aids in visualizing asymmetry. 123 I .TTTTTTIIIIIIIIIIIIIIIIITITIITIITIIIIIIIIIIIIIIII]111111111111. 60 50 4O 1111]lllillllLllllillLLLiLLJil 3O 20 IIIIIITTIIIIIIIIIIIIIIflTIITTIIIIIITIIIIITIIIITPTWI 11111111111111llllllJlllllllll ITIITTIII O ‘i111111ll1111111111111111441111111111111111111111111111111111 10 20 30 4O 50 60 0 Figure 51. Contour of the error profile of the horizontal Laplacian of a Gaussian function for the sixth-order second derivative scheme. Dotted lines are shown as aids in visualizing asymmetry. 124 Then equating the coefficients of the Taylor expansion of Equation 100 and letting h+ = z,“ — z,- and h- = z,- — 2,-_1 gives 0th Order (ff): 0 : h+a+h_b+3a+3fl (101) h2 2 lst Order (f,-"): 1 = —1'-"a — 2—Tb + h+a — 12-6 (102) hi hi hi hi 2nd Order (fz ): 0 = -3Ta + Eb + Ea + 5TB. (103) The scheme used in the code is a member of the one-parameter family represented by Equation 101—Equation 103, with the parameter chosen to make the scheme con- sistent with the first derivative cubic spline spline scheme detailed in Section A.1.2. This leads to the following values of the coefficients: a _ 71’ b- 73’ a _ 71’ 91— —h_, c _ ——4! (11+ +h_)f,- . (104) Boundaries The second derivative at the top boundary uses a one—sided forward differencing scheme of the form f{’ = a(f2 - f1) + afi + fifé- (105) Matching coefficients and solving the subsequent equations yields a second-order boundary scheme with (106) 125 where h+ = 22 — 21. The bottom boundary scheme is derived in a similar manner, using a backward one-sided differencing scheme. The form is similar to Equation 105, frlz’z =a(fnz_fnz—l)+af7,1z+fif7’13—1' (107) The coefficients of this second—order system are 6 —E’ 4 2 (I — K, )6 = E, (108) a: where h- = 2m —- 2,124. 126 REFERENCES Alfvén, H. 1950, Cosmical Electrodynamics, Claredon Press, Oxford Bercik, D. J., Basu, S., Georgobiani, D.,Nordlund, A., & Stein, R. F. 1998, in ASP Conf. Ser. 154, The Tenth Cambridge Workshop on Cool Stars, Stellar Systems and the Sun, ed. R. A. Donahue & J. A. Bookbinder, (San Francisco: ASP), CD-568. Berger, T. E., Schrijver, C. J., Shine, R. A., Tarbell, T. D., Title, A. M., & Scharmer, G. 1995, ApJ, 454, 531 Berger, T. E., & Title, A. M. 1996, ApJ, 463, 365 Blanchflower, S. M., Rucklidge, A. M., & Weiss, N. O. 1998, M.N.R.A.S., 301, 593 Bonet, J. A., Sobotka, M. & Vézquez, M. 1995, A&A, 296, 241 Cattaneo, F. 1999, ApJ, 515, L39 Collatz, L. 1966, The Numerical Treatment of Differential Equations, Springer- Verlag, New York, 538 Cubasch, U., & Voss, R. 2000, Space Science Reviews, 94, 185 Deinzer, W., Hensler, G., Schiissler, M., & Weisshaar, E. 1984a, A&A, 139, 426 Deinzer, W., Hensler, G., Schiissler, M., & Weisshaar, E. 1984b, A&A, 139, 435 Elmore, D. F., Lites, B. W., Tomczyk, S., Skumanich, A. P., Dunn, R. B., Schuenke, J. A., Streander, K. V., Leach, T. W., Chambellan, C. W., Hull, H. K., & Lacey, L. B. 1992, Proc. SPIEE 1746, 22 Ferrari, A., Massaglia, S., Kalkofen, W., Rosner, R., & Bodo, G. 1985, ApJ, 298, 181 Furusawa, K., & Sakai, J.-I. 2000, ApJ, 540, 1156 Galloway, D. E., Proctor, M. R. E., & Weiss, N. O. 1978, J. Fluid Mech., 87, 243 Galloway, D. E., & Moore, D. R. 1979, Geophys. Ap. Fluid Dyn., 12, 73 Galloway, D. E., & Proctor, M. R. E. 1983, Geophys. Ap. Fluid Dyn., 24, 109 Graham, E. 1977, in Lecture Notes in Physics, Vol. 71, Problems of Stellar Convection (Berlin: Springer), 151 127 Grossmann-Doerth, U., Knéilker, M., Schiissler, M., & Weisshaar, E. 1989, in Solar and Stellar Granulation, ed. R. J. Rutten & G. Severino (Dordrecht: Kluwer), 191 Grossmann-Doerth, U., Knolker, M., Schiissler, M., & Solanki, S. K. 1994, A&A, 285, 648 Grossmann-Doerth, U., Keller, C. U., & Schiissler, M. 1996, A&A, 315, 610 Grossmann-Doerth, U., Schiissler, M., Sigwarth, M., & Steiner, O. 1998, A&A, 337, 928 Grossmann—Doerth, U., Schiissler, M., Sigwarth, M., & Steiner, O. 2000, A&A, 357, 351 Gustafsson, B. 1973, Upps. Astr. Obs. Ann. 5(6) Gustafsson, 13., Bell, 11., Eriksson, K., & Nordlund, A. 1975, A&A, 42, 407 Hale, G. E. 1908, ApJ, 28, 315 Hurlburt, N. E., Toomre, J., & Massaguer, J. M. 1984, ApJ, 282, 557 Hurlburt, N. E., & Toomre, J. 1988, ApJ, 327, 920 Hurlburt, N. E., Matthews, P. C., & Proctor, M. R. E. 1996, ApJ, 457, 933 Hurlburt, N. E. & Rucklidge, A. M. 2000, MNRAS, 314, 793 Hyman, J. 1979, in R. Vichnevetsky, R. S. Stepleman (eds.), Adv. in Comp. Meth. for PDE’s—III, 313 Kalkofen, W., Rosner, R., Ferrari, A., & Massaglia, S. 1986, ApJ, 304, 519 Kalkofen, W., Bodo, G., Massaglia, S., & Rossi, P. 1989, in Solar and Stellar Gran- ulation, ed. R. J. Rutten & G. Severino (Dordrecht: Kluwer), 571 Keil, S. L., Balasubramaniam, K. S., Smaldone, L. A., & Reger, B. 1999, ApJ, 510, 422 Keller, C. U., Aebersold, F., Egger, U., Povel, H. P., Steiner, P., & Stenflo, J. O. 1992, LEST Technical Report, 53 Keller, C. U. 1995, Reviews of Modern Astronomy, 8, 27 Knélker, M., Schiissler, M., & Weisshaar, E. 1988, A&A, 194, 257 Knolker, M., & Schiissler, M. 1988, A&A, 202, 275 Knélker, M., Grossmann-Doerth, U., Schuessler, M., & Weisshaar, E. 1991, Advances in Space Research, 11, 285 128 Koschinsky, M., Kneer, F., & Hirzberger, J. 2001, A&A, 365, 588 Lele, S. K. 1992, J. Comput. Phys., 103, 16 Leka, K. D., & Skumanich, A. 1998, ApJ, 507, 454 Lin, H. 1995, ApJ, 446, 421 Lin, H., & Rimmele, T. 1999, ApJ, 514, 448 Lites, B. W., Skumanich, A., & Martinez Pillet, V. 1998, A&A, 333, 1053 Meunier, N., Solanki, S. K., & Livingston, W. C. 1998, A&A, 331, 771 Moreno—Insertis, F. 1994, in Solar Magnetic Fields, ed. M. Schiissler & W. Schmidt (Cambridge Univ. Press), 117 Muller, R., & Roudier, Th. 1984, Sol. Phys., 94, 33 Muller, R., & Roudier, Th. 1992, Sol. Phys., 141, 27 Muller, R., Dollfus, A., Montagne, M., Moity, J., & Vigneau, J. 2000, A&A, 359, 373 Nordlund, A. 1982, A&A, 107, 1 Nordlund, A. 1984a, in The Hydromagnetics of the Sun, ed. T. D. Guyenne (ESA SP-220; Noordwijkerhout: BSA), 37 Nordlund, A. 1984b, in Small-Scale Dynamical Processes in Quiet Stellar Atmo- spheres, ed. S. L. Keil (Sunspot, NM: Sacramento Peak Observatory), 181 Nordlund, A. 1985, s01. Phys., 100, 209 Nordlund, A., & Stein, R. F. 1990a, Computer Phys. Communications, 59, 119 Nordlund, A., & Stein, R. F. 1990b, in Solar Photosphere: Structure, Convection, Magnetic Fields, IAU Symp. 138, ed. J. O. Stenflo, 191 Nordlund, A., & Stein, R. F. 1996, in Proc. 32d Liege Colloq., ed. A. Noels et al. (Liege: Univ. Liege, Inst. d’Astrophys.), 75 Ofman, L., Klimchuk, J. A. & Davila, J. M. 1998, ApJ, 493, 474 Parker, E. N. 1963, ApJ, 138, 552 Parker, E. N. 1979, Cosmical Magnetic Fields, Claredon Press, Oxford Peckover, R. S., & Weiss, N. O. 1978, M.N.R.A.S., 182, 189 Porter, D. H., & Woodward, P. R. 1994, ApJS, 93, 309 Proctor, M. R. E., & Weiss, N. O. 1978, in Rotating Fluids in Geophysics, ed. P. H. Roberts & A. M. Soward (London: Academic Press), 389 129 Proctor, M. R. E., & Galloway, D. J. 1979, J. Fluid Mech., 90, 273 Rabin D. 1992, ApJ, 390, L103 Roudier, Th., Espagnet, O., Muller, R., & Vigneau, J. 1994, A&A, 287, 982 Riidiger, G. 1994, in Solar Magnetic Fields, ed. M. Schiissler & W. Schmidt (Cam- bridge Univ. Press), 77 Sénchez Almeida, J., & Lites, B. W. 2000, ApJ, 532, 1215 Schiissler, M. 1993, in IAU Symposium 157, The Cosmic Dynamo, ed. F. Krause, K.-H. Radler & G. Riidiger (Dordrecht: Kluwer), 27 Sigwarth, M., Balasubramaniam, K. S., Knéilker, M., & Schmidt, W. 1999, A&A, 349, 941 Skartlien, R., Stein, R. F., & Nordlund, A. 2000, ApJ, 541, 468 Sobotka, M., Bonet, J. A., & Vazquez, M. 1994, ApJ, 426, 404 Soltau, D. 1997, A&A, 317, 586 Spruit, H. C. 1976, S01. Phys., 50, 269 Spruit, H. o, Nordlund, A., & Title, A. 1990, ARA&A, 28, 263 Spruit, H. C. 1997, Mem. Soc. Astron. Italiana, 68, 397 Stein, R. F., & Nordlund, A. 1989, ApJ, 342, L95 Stein, R. F., Brandenburg, A., & Nordlund, A. 1992, in ASP Conf. Ser. 26, The Seventh Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, ed. M. S. Giampapa & J. A. Bookbinder, (San Francisco: ASP), 148 Stein, R. F., & Nordlund, A. 1998, ApJ, 499, 914 Steiner, O. & Stenflo, J. O. 1989, in Solar Photosphere: Structure, Convection, Magnetic Fields, IAU Symp. 138, ed. J. O. Stenflo, 181 Steiner, O., Kniilker, M. & Schiissler, M. 1994, in NATO ASI Ser. 433, Solar Surface Magnetism, ed. R. J. Rutten & C. J. Schrijver (DordrechtzKluwer), 441 Steiner, O., Grossmann-Doerth, U., Schussler, M., & Knolker, M. 1996, S01. Phys., 164, 223 Steiner, O., Grossmann-Doerth, U., anilker, M. & Schiissler, M. 1998, ApJ, 495, 468 Stolpe, F., & Kneer, F. 1997, A&A, 317, 942 Stolpe, F., & Kneer, F. 2000, A&A, 353, 1094 130 Siitterlin, P. 1998, A&A, 333, 305 Tao, L., Weiss, N. O., Brownjohn, D. P., & Proctor, M. R. E. 1998, ApJ, L39 Title, A. M., Topka, K. P., Tarbell, T. D., Schmidt, W., Balke, C., & Scharmer, G. 1992, ApJ, 393, 782 van Ballegooijen, A. A., Nisenson, P., Noyes, R. W., Liifdahl, M. G., Stein, R. F., Nordlund, A., & Krishnakumar, V. 1998, ApJ, 509, 435 Weiss, N. O. 1964, Phil. Trans. Roy. Soc. London, A., 256, 99 Weiss, N. O. 1966, Proc. Roy. Soc. London A, 293, 310 Weiss, N. O. 1981, J. Fluid Mech., 108, 247 Weiss, N. O. 1981, J. Fluid Mech., 108, 273 Weiss, N. O., Brownjohn, D. P., Hurlburt, N. E., & Proctor, M. R. E. 1990, M.N.R.A.S., 245, 434 Weiss, N. O., Brownjohn, D. P., Matthews, P. C., & Proctor, M. R. E. 1996, M.N.R.A.S., 283, 1153 Zhang, H., Scharmer, G., Lofdahl, M., & Yi, Z. 1998, Sol. Phys., 183, 283 Zirin, H. & Wang, H. 1992, ApJ, 385, L27 131 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIII ‘WWWWW