Exploring Convection and Dynamos in the Cores andEnvelopes of Stars
by
Nicholas Andrew Featherstone
B.S., North Carolina State University, 2003
M.S., University of Colorado, 2005

A thesis submitted to the
Faculty of the Graduate School of the
University of Colorado in partial fulfillmentof the requirements for the degree ofDoctor of Philosophy
Department of Astrophysical and Planetary Sciences
2011This thesis entitled:
Exploring Convection and Dynamos in the Cores and Envelopes of Stars
written by Nicholas Andrew Featherstone
has been approved for the Department of Astrophysical and Planetary Sciences

Juri Toomre

Bradley W. Hindman

Allan Sacha BrunMark P. Rast

Philip J. Armitage

Keith Julien

Date

The final copy of this thesis has been examined by the signatories, and we find that both thecontent and the form meet acceptable presentation standards of scholarly work in the above
mentioned discipline.

iii
Featherstone, Nicholas Andrew (Ph.D., Astrophysical and Planetary Sciences)
Exploring Convection and Dynamos in the Cores and Envelopes of Stars
Thesis directed by Profs. Juri Toomre and Bradley W. Hindman

We present theoretical studies in two complementary areas dealing with convection, rotation,and dynamos in A-type stars, and with local helioseismology in the Sun. Our studies begin withthe main-sequence A stars (stars of about 2 M⊙ ) that possess a radiative envelope overlying a
convective core. Using 3-D simulations with the Anelastic Spherical Harmonic (ASH) code to
study full spherical domains, we examine the effects of a primordial magnetic field on the dynamo
action realized in the turbulent core. Dynamo activity realized in the presence of such a field issignificantly more efficient than in its absence, yielding magnetic energies that are roughly tenfold
those of the kinetic energy associated with the convective motions. Both convective motions andmagnetic fields assume a decidedly global-scale topology in this regime, with convective downdrafts
from one side of the core streaming freely across the rotation axis, advecting and stretching magnetic
fields across distant portions of the core in the process. We examine the topology of these strong
magnetic fields and aspects of their generation in this super-equipartition dynamo.
We next develop a 3-D inversion method for helioseismic measurements of horizontal flows
obtained using ring-diagram analysis. Helioseismology uses the broad range of acoustic oscillations
observed at the solar surface to study properties deep within the Sun. Our inversion method (called
ARRDI) incorporates measurements of the wavefield made at multiple horizontal resolutions todiscern the subsurface structure of horizontal flows within the star. We adopt a regularized leastsquares (RLS) approach for these inversions and develop a novel iterative extension to the RLS
scheme wherein the flow field across the entire solar disk may be efficiently recovered. We havecalculated the set of 3-D sensitivity kernels necessary for the application of our inversion techniqueto MDI data. We explore the horizontal- and depth-averaging properties of these sensitivity kernels,
and find they differ substantially between measurements made at different horizontal resolutions.

iv
After characterizing the errors and averaging properties of our inversion algorithm, we examine
the subsurface flows around sunspots. We find that sunspots possess outflows which extend to adepth of 10 Mm. These outflows possess a noticeable two-component structure, characterized bya near-surface moat outflow and another deeper outflow at 5 Mm. Our 3-D inversion procedure
should be very useful in interpreting the vast helioseismic data sets now becoming available.

Dedication

For Manissa. For the songs we have sung together as well as those we’ve yet to sing.

vi

Acknowledgements

Throughout my tenure as a graduate student, Juri Toomre and Brad Hindman have beena consistent source of guidance and encouragement. Our many scientific discussions have helped
to shape the themes pursued in this thesis, and their attention to detail has been invaluable incrafting its final presentation. Juri has fostered a research setting within his group wherein the
pursuit of such disparate topics as those presented within this thesis can not only be undertaken,
but can flourish. Working with Brad, I was introduced to the novel notion of peering beneath thesurface of the Sun, and it was under Juri’s guidance that I first gained footing in the theoretical
exploration of magnetism in the deep interiors of stars. Together, these two have played a pivotal
role in my development as a scientist, and words cannot begin to express the debt of gratitude I
feel toward them for this.
I am also grateful to Sacha Brun and Matthew Browning who journeyed to Boulder several
times each year to help see the A-star work to completion. Their thoughtful suggestions were crucial
in deciding upon which avenues of that work to emphasize and which to abandon all together. I
thank Deborah Haber for her invaluable assistance as I learned the essentials of ring analysis. I
would also like to acknowledge Michael Thompson for his many useful suggestions and contributions
to the development of ARRDI, and I thank Aaron Birch for his kernel calculation routines. I
gratefully acknowledge the many helpful conversations with Mark Miesch, Benjamin Brown, andKyle Augustson concerning convection and dynamo theory that have occured over these past few
years. Finally, I thank Gwen Dickinson, Pete Ruprecht, and Jim McKown. Their tireless work
ensured that my own energies were focused on the science presented in this thesis.

viiThis research has been supported by NASA GSRP awards NNX07AP34H and NNX06AI29H,
by NASA Heliophysics Theory Program grants NNX08AI57G and NASA NNG05G124G, by NASA
Living with a Star grants NNG05GM83G, NAG5-13520, and NAG5-12491, by NASA Heliophysics
Guest Investigator grant NNX08AJ08G, and finally, by NSF grant ATM-0219581. The simulations
presented herein were carried out with NSF PACI support of PSC, NCSA and TACC, and by
NASA HEC support at Project Columbia.

viii

Contents

Chapter
1 Overview of Convection and Magnetism in Stars

1

1.1

Stellar Magnetic Field Measurements . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

1.2

Magnetism in Low Mass Stars . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

1.3

Magnetism in Massive Stars . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6

1.4

Seeking to Address Magnetism in A-type Stars . . . . . . . . . . . . . . . . . . . . .

8

1.5

Magnetism in the Sun . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

1.6

Acoustic Oscillations of the Solar Interior . . . . . . . . . . . . . . . . . . . . . . . . 12

1.7

Probing Solar Magnetism through Local Helioseismology

1.8

1.7.1

Local Helioseismic Techniques . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

1.7.2

Meridional Circulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19

1.7.3

Flows Around Active Regions . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

1.7.4

Active Region Influence on Global-Scale Flows . . . . . . . . . . . . . . . . . 23

1.7.5

Flows Around Sunspots . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25

Achieving 3-D Inversions of Ring-Analysis Measurements . . . . . . . . . . . . . . . 26

2 Primordial Magnetic Fields with Core Dynamos in A-type Stars
2.1. . . . . . . . . . . . . . . 15

Ap Stars and Magnetism

28. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

2.1.1

Observations of Magnetic Fields . . . . . . . . . . . . . . . . . . . . . . . . . 29

2.1.2

Possible Primordial Origin of Magnetism. . . . . . . . . . . . . . . . . . . . 30

ix

2.2

2.3

2.1.3

Dynamic Origins for Fields . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

2.1.4

Interaction of Core Dynamo with Fossil Field . . . . . . . . . . . . . . . . . . 33

Simulation Elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.2.1

Anelastic MHD Approach

2.2.2

Implementation of ASH . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37

2.2.3

Modeling the A-type star . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38

2.2.4

Imposing a Fossil Magnetic Field . . . . . . . . . . . . . . . . . . . . . . . . . 39

Core Dynamo in Progenitor. . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

3 Surprising Implications of a Primordial Magnetic Field
3.1

Modified Dynamo Action . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42
3.1.1

3.2

42

Temporal Evolution of Energies in Case A . . . . . . . . . . . . . . . . . . . . 46

Dynamics of Super-Equipartition State. . . . . . . . . . . . . . . . . . . . . . . . . 47

3.2.1

Overall Properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

3.2.2

Complex Interplay of Flows and Field Structures . . . . . . . . . . . . . . . . 50

3.2.3

Statistical Properties of Super-Equipartition State . . . . . . . . . . . . . . . 51

3.3

Growth of Axisymmetric Helical Fields . . . . . . . . . . . . . . . . . . . . . . . . . . 53

3.4

Nature of Differential Rotation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55

3.5

Distribution of Magnetic Energy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57

3.6

Growth of Large-Scale Magnetic Structure . . . . . . . . . . . . . . . . . . . . . . . . 58

4 Maintaining the Strong Dynamo
4.1

62

Sustaining a Super-Equipartition State. . . . . . . . . . . . . . . . . . . . . . . . . 62

4.1.1

Surviving Lorentz Feedbacks . . . . . . . . . . . . . . . . . . . . . . . . . . . 64

4.1.2

Global Connectivity of Structures

4.1.3

Details of Magnetic Field Generation . . . . . . . . . . . . . . . . . . . . . . . 69

. . . . . . . . . . . . . . . . . . . . . . . . 67

4.2

Multipolar External Fields. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72

4.3

Conclusions and Perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74

x
4.3.1

Sustaining the Strong Field Dynamo . . . . . . . . . . . . . . . . . . . . . . . 77

4.3.2

Possible Field Emergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80

5 Helioseismic Ring Analyses with 3-D Inversions

82

5.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82

5.2

Ring-Analysis Technique. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85

5.2.1

Obtaining Doppler Measurements. . . . . . . . . . . . . . . . . . . . . . . . 85

5.2.2

Measuring Flows with Ring Analysis . . . . . . . . . . . . . . . . . . . . . . . 86

5.2.3

Interpretation of the Fitted Velocities . . . . . . . . . . . . . . . . . . . . . . 88

5.2.4

Characterization of the Data and its Errors . . . . . . . . . . . . . . . . . . . 91

6 Properties of the Sensitivity Kernels
6.1

97

Development of the Sensitivity Kernels. . . . . . . . . . . . . . . . . . . . . . . . . 97

6.1.1

Modeling the Solar Wavefield . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

6.1.2

Assembly of the Sensitivity Kernels . . . . . . . . . . . . . . . . . . . . . . . . 102

6.1.3

Structure of the Sensitivity Kernels. . . . . . . . . . . . . . . . . . . . . . . 103

7 Attaining 3-D Inversions of Ring Analysis Data

112

7.1

7.2

Elements of the Inversion Procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . 112
7.1.1

Regularized Least Squares Approach . . . . . . . . . . . . . . . . . . . . . . . 113

7.1.2

Stitched Inversions. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116

Averaging Properties of the Inversion Algorithm . . . . . . . . . . . . . . . . . . . . 119
7.2.1

Data Combination Properties of the Inversion Algorithm

7.2.2

Averaging Kernels . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122

7.2.3

Effects of Regularization. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125

8 Exploring Convection and Magnetism with ARRDI Inversions
8.1

Comparison of ARRDI Results with Prior Assessments
8.1.1. . . . . . . . . . . 119

129. . . . . . . . . . . . . . . . 129

Surface Comparisons with the 2◦ data . . . . . . . . . . . . . . . . . . . . . . 131

xi

8.2

8.3

8.1.2

Comparisons at Depth: Flowmaps and Their Relation to the 4◦ and 16◦ Data.135

8.1.3

Establishment of Large-Scale Means . . . . . . . . . . . . . . . . . . . . . . . 137

Sunspot Outflows . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140
8.2.1

Ubiquitous Sunspot Outflows at Depth. . . . . . . . . . . . . . . . . . . . . 140

8.2.2

Temporal Evolution of Sunspot Outflows . . . . . . . . . . . . . . . . . . . . 144

8.2.3

Sunspot Outflows as Revealed through ARRDI and Time-Distance Inversions 146

Further Reflections on 3-D Inversions with ARRDI. . . . . . . . . . . . . . . . . . 150

9 Reflections and Future Directions
9.1

9.2

9.3

153

Convective Cores and Primordial Magnetic Fields in A stars . . . . . . . . . . . . . . 153
9.1.1

Possible Magnetic Buoyancy Instabilities . . . . . . . . . . . . . . . . . . . . . 153

9.1.2

Convective Overshooting With and Without Primordial Fields . . . . . . . . 155

9.1.3

Primordial Magnetism and Dynamos in Pre-Main Sequence Stars . . . . . . . 156

Possible Further Improvements to 3-D Inversions with ARRDI . . . . . . . . . . . . 158
9.2.1

Obtaining Optimal Horizontal and Vertical Resolution . . . . . . . . . . . . . 159

9.2.2

Further Possibilities for Improvements to the Sensitivity Kernels . . . . . . . 160

Future Prospects for Ring Analysis with 3-D Inversions . . . . . . . . . . . . . . . . 161
9.3.1

Detailed Probing of Sunspot Outflows . . . . . . . . . . . . . . . . . . . . . . 162

9.3.2

Probing Supergranulation in the Near-Surface Shear Layer

9.3.3

Discerning the Surface Signal of the Solar Giant Cells . . . . . . . . . . . . . 163

Bibliography. . . . . . . . . . 162

165

xii

TablesTable
3.1

Time and volume averaged magnetic energies compared to kinetic energies. . . . . . 59

8.1

Target depths and horizontal resolution of the inversion grid for the 3-D inversionsof the April 24, 2001 data.

8.2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131

Target depths and horizontal resolution of the inversion grid for the 3-D inversionsof sunspots from CR 1985.

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143

xiii

FiguresFigure1.1

Zeeman Doppler Imaging of the dynamo cycle in the F star τ Boo. . . . . . . . . . .

4

1.2

Magnetic Doppler Imaging of the rapidly rotating Ap star HR3831. . . . . . . . . . .

7

1.3

Simulations with ASH of core convection in a rotating A-type star . . . . . . . . . .

8

1.4

The many aspects of magnetism in the Sun. . . . . . . . . . . . . . . . . . . . . . . . 11

1.5

Probing solar interior rotation with global helioseismology. . . . . . . . . . . . . . . . 13

1.6

Subsurface flows as viewed with ring-analysis helioseismology . . . . . . . . . . . . . 16

1.7

Solar power spectrum of a localized patch of the Sun as viewed when cut at three
temporal frequencies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

1.8

Example of a time-distance cross covariance function. . . . . . . . . . . . . . . . . . 18

1.9

Meridional circulations as measured using ring-analysis helioseismology . . . . . . . . 20

1.10 Surface inflows into active regions transition to outflows at depth. . . . . . . . . . . . 22
1.11 Torsional oscillations and meridional flow variations over the solar cycle. . . . . . . . 24
2.1

Flow and magnetic field properties of the core convection dynamo in the progenitor
case C4m. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

3.1

Superimposing a fossil field into the progenitor case to initiate case A. . . . . . . . . 43

3.2

Separating the role of flux and helicity. . . . . . . . . . . . . . . . . . . . . . . . . . . 45

3.3

Temporal evolution of overall energy densities . . . . . . . . . . . . . . . . . . . . . . 46

3.4

Snapshot of radial velocity vr and the two magnetic field components Br and Bφ . . 48

xiv
3.5

Volume rendering of magnetic stuctures (with rotation axis vertical) near day 15,000
accompanying Fig 3.4. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51

3.6

Sequence of views showing the evolution of magnetic energy density . . . . . . . . . 52

3.7

Time-averaged pdf’s for case A and case C4m. . . . . . . . . . . . . . . . . . . . . . 53

3.8

Temporal evolution of the mean toroidal and poloidal magnetic fields for case A. . . 54

3.9

Differential rotation realized in Case A. . . . . . . . . . . . . . . . . . . . . . . . . . 56

3.10 Variation with proportional radius of horizontally-averaged kinetic and magnetic
energy densities for case A. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
3.11 Magnetic and kinetic energy spectra sampled at mid-core in case A for early andlate times. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
3.12 Volume rendering of magnetic energy density viewing the entire convection zone . . 61
4.1

Dual rendering of kinetic energy and magnetic energy near day 15,000. . . . . . . . . 63

4.2

Comparing partitioning of fast flows and strong fields. . . . . . . . . . . . . . . . . . 64

4.3

Probability distribution functions (pdfs) of the relative angle θ between v and B. . . 65

4.4

Columnar convection for case A visualized using instantaneous streamlines near day
15,000. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66

4.5

Sampling the evolving flow streamlines and accompanying magnetic energy density . 68

4.6

Generation and dissipation of fluctuating magnetic energy for case A. . . . . . . . . 71

4.7

Multipolar cases. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73

4.8

Mean field evolution for the quadrupolar ℓ = 2 case. . . . . . . . . . . . . . . . . . . 75

4.9

As in Fig. 4.8 showing the ℓ = 4 multipolar case. . . . . . . . . . . . . . . . . . . . . 76

5.1

Structure of a typical solar power spectrum . . . . . . . . . . . . . . . . . . . . . . . 87

5.2

The effect of magnetic field strength on ring-analysis errors . . . . . . . . . . . . . . 92

5.3

Average effects of magnetic field strength on measurements made using different
analysis-region sizes. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93

5.4

Average ring-analysis errors as realized with MDI data. . . . . . . . . . . . . . . . . 94

xv
5.5

Data filling factor over the solar disk for ring-analysis measurements . . . . . . . . . 96

6.1

Sample sensitivity kernel for a 16◦ measurement. . . . . . . . . . . . . . . . . . . . . 104

6.2

Variation of sensitivity kernel structure with horizontal wavenumber. . . . . . . . . . 106

6.3

Effects of a finite observing region on the observed power spectrum. . . . . . . . . . 108

6.4

Sample sensitivity kernel for a 2◦ measurement . . . . . . . . . . . . . . . . . . . . . 109

6.5

Variation of the vertical structure of a 2◦ sensitivity kernel . . . . . . . . . . . . . . . 110

7.1

Illustrating the need for a modified RLS inversion to bring adjacent inversion regions
into agreement. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115

7.2

Convergence properties from a full-disk inversion using the iterative stitching method.117

7.3

Comparision between an iterative stitched inversion and its companion full inversion 118

7.4

Depth-averaging properties of the inversion scheme.

7.5

Horizontal structure of the averaging kernels. . . . . . . . . . . . . . . . . . . . . . . 122

7.6

Depth properties of the averaging kernels. . . . . . . . . . . . . . . . . . . . . . . . . 124

7.7

Variation of the tile-size contribution parameter B(T, z) . . . . . . . . . . . . . . . . 125

7.8

Variation of 3-D inversion errors with changing horizontal and vertical regularization 126

7.9

Sensitivity of the averaging kernels to variations in the regularization. . . . . . . . . 128

8.1

Comparison of the ARRDI results to 2◦ f -mode fits . . . . . . . . . . . . . . . . . . 130

8.2

Closeups of the regions of interest in Figure 8.1 . . . . . . . . . . . . . . . . . . . . . 133

8.3

Comparison of ARRDI solution to fitted data a depth of 4 Mm. . . . . . . . . . . . . 134

8.4

Closeups of the regions of interest indicated in Figure 8.3 . . . . . . . . . . . . . . . 136

8.5

Flowmaps at a depth of 14 Mm resulting from the ARRDI algorithm and 1-D OLA. . . . . . . . . . . . . . . . . . 121

inversion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137
8.6

Establishment of the large-scale mean flows. . . . . . . . . . . . . . . . . . . . . . . . 138

8.7

Synoptic map of Carrington rotation 1985 at 5 Mm as measured by ARRDI . . . . . 141

8.8

The depth variation of a sunspot outflow . . . . . . . . . . . . . . . . . . . . . . . . . 142

xvi
8.9

Sampling sunspot outflows in CR 1985 . . . . . . . . . . . . . . . . . . . . . . . . . . 145

8.10 The slow evolution of a persistent sunspot . . . . . . . . . . . . . . . . . . . . . . . . 147
8.11 Emergence and disruption of a sunspot. . . . . . . . . . . . . . . . . . . . . . . . . . 148
8.12 Comparison of outflows beneath NOAA AR 9787 as revealed through time-distance
helioseismology and ring analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149
8.13 Sunspot outflow as realized through 3-D numerical simulations . . . . . . . . . . . . 151

Chapter 1

Overview of Convection and Magnetism in Stars

It has been known for over a century now that stars are magnetic entities. Magnetism was
first detected in the Sun by Hale (1908), with the first nonsolar detection by Babcock (1947)
following a few decades later. Magnetic fields are now known to be ubiquitous across the main
sequence and indeed across all stages of stellar evolution. They play an active role in stellar
formation, affecting molecular cloud collapse and modulating accretion onto protostars throughtheir interactions with the circumstellar disk. On the main sequence, magnetic fields regulate thespin down of stars through their interaction with stellar winds, later serving as seeds for the strong
magnetism observed in white dwarfs and neutron stars.
Within main sequence stars, magnetic fields rarely exist in isolation. Whether in the coresof massive stars or in the envelopes of stars more like our Sun, magnetism is inextricably linked
to the convective motions of an ionized gas. The interaction of convection and magnetism canlead to some interesting and surprising effects, such as magnetic cycles and explosive events. The
latter is of obvious importance to modern society, where solar phenomena such as coronal massejections can severely impede communications satellites and disrupt power grids. However, throughthe periodic emergence of sunspots at the solar surface, the Sun’s magnetism has perhaps exterted
a more subtle influence on society as well, fascinating ancient sungazers, notably the Chinese, long
before the existence of any theory of electromagnetism (e.g., Wittman & Xu 1987).
The methods for understanding magnetism and its generation in stars run the full gamutof theoretical to observational. This thesis delves into both of these areas, beginning with three-2dimensional (3-D) numerical simulations of convection and dynamo action in the cores of A-typestars. The convective cores of these stars lie buried beneath deep radiative envelopes that have long
been suspected of harboring intense and large-scale magnetic fields. We explore what implications
the presence of such strong, organized fields might have on the core dynamos within these stars
before turning to convection in the Sun.
The proximity of the Sun allows us to probe the interplay between its convection and magnetism in great detail through helioseismology. Following our A-star explorations, we develop a local
helioseismic technique useful for probing flows around magnetic features in the shallow near-surface
layers of the solar convection zone. This technique combines velocity measurements obtained withvariable horizontal resolution using ring-analysis helioseismology to yield a self-consistent description of the large and small-scales of convection present in the upper 20 Mm of the solar convection
zone. Before beginning, let us proceed with a general discussion of magnetism in stars, followed bya more focused discussion of what helioseismology has revealed about the interaction of magnetismand convection in our own star.

1.1

Stellar Magnetic Field Measurements
Magnetic fields in stars are measured by observing the Zeeman splitting of magneticallysensitive spectral lines. In the simplest case, a magnetic field will split a spectral line into what is
known as a Zeeman triplet, comprised of a linearly polarized component at the center of the unsplit
line, and two circularly polarized components on either side of line center. The degree of splitting
realized by the circularly polarized components is proportional to the strength of the magnetic
field. Unfortunately such triplets are difficult to resolve in all but the slowest rotators owing to
the competing effect of Doppler broadening of the line resulting from rotation of the star. Some
resolution of Zeeman triplets has been achieved using near infrared lines in cool stars such as the
M and K dwarfs (e.g., Saar & Linsky 1985; Valenti et al. 1995) and in some slowly rotating Apstars (Mathys 2001).
When the Zeeman splitting is unresolved, a broadening of the line is observed. Fitting

3of the observed line to a model accounting for Doppler broadening, magnetic fields, and otherline-formation effects allows the average (i.e., disk-integrated) magnetic flux to be estimated (e.g.,
Robinson et al. 1980a; Saar 1988). This technique is insensitive to any small-scale tangling present
in the magnetic field and has found wide use as a general diagnostic of the field strength in stars,
particularly in the cool dwarf stars (e.g., Robinson et al. 1980b) and protostars such as the classicalT Tauri stars (Johns-Krull 1997). While measurements of the splitting or broadening of the line
provide only an estimate of the magnetic field modulus, more sophisticated techniques enable therecovery of the magnetic field geometry.
The magnetic field orientation, it turns out, is encoded in the polarization properties ofthe spectral line. The splitting of the circularly polarized components of the line can be used todeduce the intensity-weighted disk-average of the line-of-sight magnetic field. Similarly, the linear
polarization allows one to infer the intensity-weighted averaged of the transverse field strength.
These simple measures of the field orientation are good for large-scale fields with little variation
across the stellar disk and have been applied widely to the persistent fields of the chemically peculiar
Ap stars (e.g., Borra & Landstreet 1980; Leroy 1995).
Initial attempts to discern the geometry of stellar magnetic fields began with the large-scale
magnetic fields of the Ap stars. The dipolar component of the magnetic field was obtained byobserving a changing spectral line over several rotation periods. These variations were fitted to amodel of what might be expected from a dipole inclined with respect to the rotation axis (e.g.,
Preston 1967). Zeeman Doppler Imaging (ZDI; Semel 1989; Brown et al. 1991) and the related
technique of Magnetic Doppler Imaging (MDI; Piskunov & Kochukhov 2002) have expanded on
this general idea through the use of higher-order multipoles for both the poloidal and toroidal fields
in the fitting procedure. They have yielded a much more detailed understanding of the surface
toplogies of magnetic fields in a wide variety of stars.
The advantages of surface-field imaging techniques such as ZDI are readily apparent in Figure
1.1 where an F7 star has been imaged three times over the course of three years. The poloidal field isclearly non-axisymmetric in nature and is shown to flip over the course of the three years observed,

4

Figure 1.1: Zeeman Doppler Imaging of the dynamo cycle in the F star τ Boo. Polar view ofmagnetic fields imaged in (a) June 2006 (Catala et al. 2007), (b) June 2008 (Donati et al. 2008)
and (c) June 2009 (Donati et al. 2009). Concentric circles indicate lines of latitude. The stellar
equator is indicatd by the thick black line and a latitude of −30◦ by the outer solid line. Bothpoloidal and toroidal components of the field are seen to reverse between subsequent years, hinting
at the presence of a dynamo with a cycle roughly ten times faster than that found in the Sun.
(Adapted from Donati & Landstreet 2009)

as does the more axisymmetric toroidal field. The ability to discern the topology of magnetic fields
on stars through techniques such as ZDI and MDI holds great promise for observing stellar dynamo
behavior. These techniques have also been used to analyze the magnetic fields of the more massive
stars, such as the chemically peculiar Ap and Bp stars (e.g., Kochukhov et al. 2004) as well as
magnetic fields in the pre-main-sequence (PMS) T Tauri stars (e.g., Donati et al. 2008, 2010;

5
Hussain et al. 2008).

1.2

Magnetism in Low Mass Stars
Stars across the H-R diagram fall into two broad categories: those with convective cores andthose with convective envelopes. Stars with masses less than about 1.5 M⊙ (spectral type F andlater) have convective envelopes that deepen with decreasing mass, ultimately encompassing theentire stellar radius in the case of the fully convective M-dwarfs. The outer convection zones ofthese stars exhibit a wide range of solar-like phenomena. Dark star spots that persist on timescales of days have been observed in numerous studies (e.g., Radick 1982; Dorren & Guinan 1982;
Richards 1990; Barnes et al. 2004). The absorption signatures of prominence-like clouds of cool
plasma, locked into corotation with the star by coronal magnetic fields, have also been observed(Collier Cameron et al. 1989; Donati et al. 2000), suggesting that a magnetic corona may bepresent in these stars. Cyclic chromospheric emission in solar-like stars observed through calcium
H and K has hinted at the presence of dynamo cycles in cool stars (e.g., Wilson 1978; Baliunas etal. 1998). More recently, as shown in Figure 1.1, ZDI has demonstrated the occurrence of a fieldpolarity reversal in an F star star (Catala et al. 2007; Donati et al. 2008, 2009).
The surface magnetism observed in the low-mass stars is thought to result from the presenceof a rotationally influenced convective dynamo, as opposed to a primordial field. Indeed, withinstars of a given mass, tracers of magnetic activity, such as X-ray flaring and chromospheric HK
emission, increase with the stellar rotation rate (e.g., Noyes et al. 1984). This trend of increased
magnetic activity ceases for all low mass stars at some mass-dependent rotation rate for all of thelow-mass stars, such that increased rotation does not further strengthen the magnetic fields (e.g.,
Pizzolato et al. 2003; Reiners et al. 2009). Presumably the fields become so strong as to feed backon the convection through some as yet unknown dynamo quenching mechanism.

6

1.3

Magnetism in Massive Stars
Stars more massive than about 1.5M⊙ possess deep outer envelopes that are primarily qui-escent and radiative in nature. While stars in the 1.5M⊙ -3M⊙ range are thought to possess thinconvection zones near their surface, driven by the ionization of hydrogen and helium, strong convection is mainly driven deep within the cores of these stars where vigorous burning of hydrogenis accomplished largely through the CNO cycle. Within the stars in this mass range (O,B,A, andearly F) magnetism has been observed primarily in the chemically peculiar Ap and Bp stars, butsome recent observations have begun to demonstrate the presence of magnetism in the normal O
and B stars (e.g., Donati et al. 2006; Bouret et al. 2008). Searches for magnetism in the normal A
and early F stars have been both extensive and fruitless (e.g., Shorlin et al. 2002; Bagnulo et al.
2006).
The chemically peculiar Ap and Bp stars are so named for the strong surface abundances
they exhibit in Si and many other metals (e.g., Hg, Cr, Eu; Mestel 1999). These stars constitute a
small percentage of the total A and B stars, peaking at 10% for the 3 M⊙ stars and tapering off tonearly zero at 1.5 M⊙ (Power et al. 2008). All known Ap stars of the Si and Sr-Cr-Eu peculiarity
classes are observed to be magnetic, with field strengths ranging from 300 G up to 30,000 G (e.g.,
Babcock 1947; Auriere et al. 2007). However, searches for magnetism in the other pecularity classes(such as the Hg and Mn peculiarity classes) have failed to find magnetic fields (Shorlin et al. 2002).
Imaging of these stars through MDI and ZDI has shown that the different peculiar elemental
abundances tend to concentrate in different regions of the star and do not always reflect the geometry of the magnetic field (see Figure 1.2). The observed surface fields are typically large-scale
and predominantly poloidal in nature, and when variable, they do so only at the stellar rotation
rate (e.g., Mestel 1999; Donati & Landstreet 2009). The large-scale and persistent nature of theobserved magnetic fields in the Ap and Bp stars, uncharacteristic of dynamo driven fields, suggeststhat these fields may be of primordial origin, deriving from a remnant magnetic field amplified
during the collapse of the star from its originating molecular cloud.

7

Figure 1.2: Magnetic Doppler Imaging of the rapidly rotating Ap star HR3831. (upper row ) Radial
component of the magnetic field imaged at five intervals spanning one rotation period. Many Apstars exhibit a relatively simple surface magnetic field topology. Here the field is well approximatedwith a dipole inclined by 90◦ with respect to the rotation axis. (lower rows) Surface maps ofthe relative abundances of different elements. ε(M ) ≡ log(NM /Ntotal ). Abundance maps differ
between elements. The Li abundance shows a strong correlation with the dipolar field, while otherabundances, such as Si and Mg, do not. (Adapted from Kochukhov et al. 2004)

Magnetic fields have been observed in more massive stars with some limited success. Detection of magnetic fields in the normal O and B stars are difficult owing to the intrinsic width andweakness of the magnetically sensitive spectral lines in these stars. Only a handful of B stars (six)have been observed to have magnetic fields (Henrichs et al. 2000; Neiner et al. 2003; Donati etal. 2006c; Petit et al. 2008b), and an even smaller number of O stars, namely three (Donati et al.

8Figure 1.3: Simulations with ASH of core convection in a rotating A-type star, with the coresurrounded by stable radiative envelope. (a) Volume rendering of the strong radial magnetic fieldachieved by dynamo action in the core, and (b) its longitudinal component. (c) Global view ofradial velocity at fixed radius near the middle of the core (downflows in blue tones). (d) Azimuthal
velocity in cross-section, with dashed lines showing the core-envelope interface and the equator.
Dark tones indicate central column of slow rotation. (e) Angular velocity variations with radius.
(Adapted from Brun, Browning & Toomre 2005)

2002, 2006b; Bouret et al. 2008). Much like the peculiar Ap and Bp stars, these few normal O andB stars are found to possess persistent, largely poloidal magnetic topologies with field strengths
in the kG range. Some notable exceptions do occur, such as with the early B star τ Sco, which
was deduced to have a rather complicated twisted toroidal topology thought to originate from theinteraction between a primordial surface field and the stellar wind (Donati et al. 2006). The static
nature of the fields in the O stars is also suspect since observations by Bouret et al. (2008) hint at
the presence of a dynamo driven field in the O star ξ Ori A.

1.4

Seeking to Address Magnetism in A-type Stars
The strong magnetic fields observed on the surfaces of Ap stars are intriguing in part for

their persistence in the apparent absence of any regenerative dynamo mechanism operating in theirradiative envelopes. Convection in the cores of these stars has long been suspected of harboringstrong dynamos, and recent 3-D numerical modeling has enabled a detailed exploration of the
dynamo states accessible to these stars (e.g., Browning, Brun & Toomre 2004; Brun, Browning &Toomre 2005). These simulations have revealed that core convection in A-type stars involves highly
time-dependent flows with complex, vigorous convection that can penetrate substantially into the9overlying radiative zone (Fig. 1.3). Such motions often span multiple scale heights and serve tocouple widely separated portions of the core. Convection tends to drive a prominent column of slow
rotation maintained throughout the cores of these stars, though this strong differential rotation isdiminished greatly in the presence of magnetic fields. When dynamo action is admitted, strong
magnetic fields are generated whose associated magnetic energy is in approximate equipartition withthat of the convective motions. Such magnetic fields were found to be largely non-axisymmetric innature, with the axisymmetric poloidal and toroidal fields contributing little to the overall energybalance.
It is thus entirely likely that intense magnetic fields are generated in the cores of A-type stars
but there are difficulties associated with transporting these fields to the stellar surface. Meridional
circulations are likely too weak to carry magnetic flux from the core to the surface on time scales
shorter than the stellar lifetime (e.g. Charbonneau & MacGregor 2001). Moreover, the transport ofdynamo-generated fields through magnetic buoyancy seems unlikely if magnetic fields attain only
equipartition levels in the cores of these stars (MacDonald & Mullan 2004). Indeed, no buoyant
magnetic fields were realized in the simulations of Brun, Browning & Toomre (2005).
It seems more plausible that the strong magnetic fields of the Ap stars originate instead
within the molecular clouds from which these stars have formed. Typical ohmic diffusion times
across the radiative envelope are comparable or longer than the stellar lifetime (e.g., Cowling 1945)and so strong magnetic fields amplified during the collapse phase may well persist throughout themain-sequence lifetime of an A star (about 109 years). An interesting aspect of the primordial fieldstory involves the 3-D simulations of Braithewaite & Spruit (2004), which found that a randomly
organized magnetic field, left to its own devices in a quiescent radiative envelope, will eventuallyrelax into a coherent large-scale field, with a twisted toroidal geometry.
Large-scale magnetic field strengths in the 10 kG range are commonly observed on the surfacesof the Ap stars and, if these fields truly are indicative of an organized subsurface field, equally strong
magnetism may well pervade the deep radiative interior. Magnetic field strengths in this range areapproaching equipartition levels with respect to the expected convective motions in the cores of10these stars (e.g., Browning, Brun & Toomre 2004). It is thus natural to ask what role a primordial
magnetic field, with some threading into the convective core, might play in shaping the behaviorof the core dynamo.
We explore this question through 3-D numerical modeling of convection in the interiors of
A stars, finding that dynamo action in the presence of a primordial field in the radiative envelope
is very different from that operating in its absence. Both convective motions and magnetic fields
assume a decidedly global-scale topology in this regime, with convective downdrafts from one
side of the star streaming freely across the rotation axis, advecting and stretching magnetic fields
across distant portions of the core in the process. Magnetic energies, which were equipartition in
the absence of the primordial magnetic field, attain levels that are tenfold that of the associatedconvective motions. We examine these super-equipartition magnetic states in detail in Chapters
2-4 before then proceeding to helioseismic studies of flows around magnetic features in the Sun.

1.5

Magnetism in the Sun
Our proximity to the Sun of course enables much more detailed observations of its magnetic

fields than is achievable in other stars. Magnetism in the Sun is most readily apparent throughthe appearance of dark sunspots on the solar surface (Figure 1.4a). The dark umbra of sunspots
are now known to harbor kG-strength fields and are typically cooler than the surrounding plasma
by about 1000-2000 K due to the inhibition of convective heat transport in the presence of suchstrong magnetism. Hosting typical diameters of a few tens of Mm, sunspots can live for weeks tomonths following their emergence before decaying and forming the small-scale flux which comprises
the active regions in which they reside (e.g., Solanki 2003).
Sunspots emerge in concentrated latitudinal bands, and it is the migration of these bandsof emergence that is the most easily recognizable signature of the 22-year solar cycle. Emergenceoccurs at mid-latitudes of ±35◦ initially and gradually migrates toward the equator over the courseof 11 years. After 11 years, the magnetic dipole of the Sun has reversed, and emergence begins
anew at the mid-latitudes. A variety of established symmetries and ordering rules, suggestive of an

11

Figure 1.4: The many aspects of magnetism in the Sun. (a) Continuum image of a sunspot groupnear disk center taken with the Swedish Solar Telescope (credit: G. Scharmer & K. Langhans).
Strong magnetic fields in the dark umbra of sunspots inhibit local convection (here prominently
visible on granular scales). (b) Magnetogram of the Sun obtained on 22 Oct 2010 showing magnetic
field along the line of sight in white and the opposite polarity in black (obtained with the Helioseismic and Magnetic Imager aboard the Solar Dynamics Observatory [SDO]). A clear preference
for polarity of the leading and trailing edges of active regions is evident within each hemisphereand opposite between hemispheres. Magnetism in the leading portion of the active region tends toobey Joy’s law, emerging at lower latitudes than in the trailing portion. (c) Image of the Sun from
6 Dec 2010 showing the eruption of a large filament, visible above the disk for over two week prior
to eruption (obtained with the Atmospheric Imaging Assembly (AIA) aboard SDO).

underlying large-scale dynamo, characterize sunspot emergence. Throughout the cycle, sunspots
emerge in pairs of opposite polarity, with the same polarity for all leading spots within the same
hemisphere (and the opposite sense in the other hemisphere). This law, known as the Hale polarity
law (Hale et al. 1919; Hale & Nicholson 1925), hints at the presence of a large-scale, subsurface
toroidal magnetic field generated by such a dynamo. Other phenomena, such as the tilting ofsunspots described by Joy’s law (Hale et al. 1919), whereby the leading polarity emerges nearer tothe equator than the trailing polarity (Figure 1.4b), and the equatorward migration of the active
region belts over the solar cycle, are similarly suggestive of a large-scale dynamo.
The prevailing picture is that solar active regions originate from subsurface magnetic flux
tubes that rise to the solar surface. The advent of helioseismology and the subsequent discovery
by it of the tachocline (see §1.6) has suggested that the base of the convection zone is the likely

12seat of the solar dynamo. Magnetic flux emerging from the tachocline must then first traverse thesolar convection zone, a region of highly turbulent flow, on its way to the solar surface. The effectsof convection on emerging magnetism are most visible in the upper convection zone where highly
vigorous convection occurs on many horizontal scales, ranging from the solar granulation (1 Mmscale) to the more global-scale flows, often termed solar subsurface weather (SSW), characterized bymeandering jet-stream-like patterns. As we discuss in detail in §1.7, flows at all scales are strongly
linked with surface magnetic features, and so it is thought that the twisting of magnetic fields bythese fluid motions, both in the photosphere and in the subsurface layers, is likely to be the driving
force behind solar eruptions such as flaring and coronal mass ejections (Figure 1.4c).

1.6

Acoustic Oscillations of the Solar Interior
Understanding the origin of the solar cycle and the associated energetic magnetic events areoutstanding challenges to modern solar physics. Serious progress is now being made in studying
the solar cycle through 3-D theoretical simulations of convection and magnetism in solar-like stars,
with some of these now beginning to exhibit dynamo cycles (e.g., Brown et al. 2010, 2011; Ghizaru,
Charbonneau & Smolarkiewicz 2010). On the observational side, helioseismic probing of the solar interior has proven to be a powerful tool for examining magnetism the upper reaches of theconvection zone, and has shown potential as a predictive tool for magnetic events.
The solar oscillations manifest as displacements of the solar surface that occur with periodsof about 5 min (3 mHz) and were observed by Leighton et al. (1962) using Mount Wilson Doppler
measurements. The prevailing interpretation now is that these oscillations arise from standing
acoustic modes trapped beneath the photosphere (Ulrich 1970; Leibacher & Stein 1971) and are
excited by near-surface turbulent convection (e.g., Goldreich & Keeley 1977; Libbrecht et al. 1986;
Goode et al. 1992). When power spectra (in horizontal wavenumber and frequency) are constructedfrom the observed Doppler data, modal power is found to fall along discrete ridges in frequencywavenumber space, with power focused at specific frequencies for each horizontal wavenumber(Figure 1.5a).

13

Figure 1.5: Probing solar interior rotation with global helioseismology. (a) Solar power spectrum
observed with MDI as a function of frequency and spherical harmonic degree ℓ. High acoustic
power is indicated in red tones, and low power by blue tones. Power is concentrated at discrete
frequencies specific to each ℓ. (b) Differential rotation of the Sun (angular velocity Ω) as realized
through inversion of global-mode frequency splittings. The fast equator is found to rotate once
every 25 days, and the slow pole completes one revolution in about 33 days. The differential
rotation transitions to solid body rotation in the radiative interior as illustrated more clearly in (c)where radial cuts at various latitudes are shown (Thompson et al. 2003).

14
A wide variety of techniques have been devised by which the properties of the solar wavefield
may be exploited to infer the internal structure and flows within the Sun. They broadly fall
into two categories. Global helioseismology employs large-scale and long-lived acoustic modes,
capable of circumnavigating the the Sun many times, to deduce properties of the solar interior as afunction of radius and latitude (e.g., Deubner & Gough 1984; Gough & Toomre 1991; ChristensenDalsgaard 2002; Thompson et al. 2003). Local helioseismic techniques, on the other hand, analyze
the wavefield over some localized portion of the solar surface to infer the local variation of somequantity (i.e., flow or temperature) beneath that region of the Sun (e.g. Braun et al. 1987; Hill
1988; Duvall et al. 1993; Lindsey & Braun 1997; Gizon et al. 2010).
Global helioseismology has been particularly useful in inferring the interior nature of thesolar differential rotation. Inversions of global-mode frequency splittings revealed that the surface
differential rotation (expressed in terms of the angular velocity Ω) imprints through to the baseof the convection zone. The differential rotation ultimately terminates in a thin layer of shear at
the base of the convection zone where a transition to solid body rotation in the radiative interior
occurs (e.g., Thompson et al. 1996, 2003; Schou et al. 1998; Figure 1.5b, c). This shearing layer,termed the tachocline (Spiegal & Zhan 1992), has become a common ingredient in solar dynamo
modeling due to its likely role in the generation of strong toroidal magnetic fields (e.g., Parker
1993; Charbonneau & MacGregor 1997; Browning et al. 2006). Global helioseismology has also
been instrumental in probing other aspects of solar interior structure, such as establishing thebase of the solar convection zone and determining the chemical composition of the solar interior.
Perhaps most notably, accurate determinations of the temperature in the solar core played a keyrole in encouraging the intensive search that revealed the existence of neutrino oscillations (see e.g.,
Christensen-Dalsgaard 2002 for a review).
The nature of global-helioseismic modes means that solar interior properties are averaged inlongitude and over both the northern and sourthern hemispheres. Details concerning the subsurface
structure of solar magnetic fields and near-surface convection, which manifest on small-scales andwhich are not axisymmetric or hemispherically symmetric, are thus inaccessible to global techniques.

15
Local helioseismic techniques do not perform such averaging and allow for full three-dimensional
mapping of regions of the Sun, making them much more useful for probing the fine-scale propertiesof the solar convection zone.

1.7

Probing Solar Magnetism through Local Helioseismology

1.7.1

Local Helioseismic Techniques

Modern local helioseismic studies typically probe flows in the upper convection zone using one
of two widely applied techniques. Ring diagram analysis infers the presence of a subsurface flow insome localized region of the Sun by fitting for Doppler shifts induced in waves propagating throughthat region. An example of the flow mapping achievable with ring-diagram analysis is shown in
Figure 1.6 where large-scale, meandering flows are seen to pervade the subsurface layers of the Sun.
Time-distance distance methods measure such flows by analyzing the wavefield in physical space,inferring the flow by measuring travel times for wave packets to propagate between two sites on
the solar surface. These two methods yield comparable results when compared (e.g., Hindman etal. 2004). However time-distance methods afford higher horizontal resolution, particularly below
the solar surface, than do current ring analysis implementations (see §1.8). We discuss these two
methods here briefly and provide a detailed desription of the ring analysis method in Chapter 5.
When viewed at a constant frequency, the 3-D power spectrum (in termporal frequency andhorizontal wavenumber) associated with a localized region of the solar surface appears as a seriesof concentric rings (Figure 1.7). Advection of waves by an underlying flow induce a Doppler shift
in these rings. Ring-analysis techniques (e.g., Hill 1988; Schou and Bogart 1998, Haber et al. 2002)assess subsurface flows by measuring anisotropies in the shifting of these rings and sound speed
variation through their isotropic shifts (e.g., Basu et al. 2004). Measurements of the frequency
shifts associated with different horizontal wavenumbers and modal orders average the subsurface
flow differently and may be combined through linear inversion to yield a map of the flow bothhorizontally and with depth.

16

Figure 1.6: Subsurface flows as viewed with ring-analysis helioseismology during (a) the magnetically active period spanning 3 April 1999 through 29 April 1999 and (b) the relatively inactive
period spanning 21 May 1997 through 16 June 1997. Blue arrows indicate subsurface flows at
a depth of 7 Mm beneath below the photosphere. Flow vectors at each latitude and Carrington
longitude represent a seven-day average over the measured flow at that coordinate. An averagemagnetogram for this time period is underlaid, with red and green indicating opposite polaritiesof strong line-of-sight magnetic field, and yellow denoting regions of weak magnetic field. Large,
meandering jet stream-like motions are visible which possess a clear poleward component withineach hemisphere. The meandering patterns of these flows are much less organized in the 1999 dataand seem to be influenced by the strong magnetic features present. (Adapted from Haber et al.
2002)

17

Figure 1.7: Solar power spectrum of a localized patch of the Sun as viewed when cut at three temporal frequencies. Power manifests as concentric rings that expand in horizontal wavenumber (kxand ky ) with increasing frequency ν. Ring-analysis methods deduce subsurface flows by measuring
Doppler shifts induced in these rings when an underlying flow is present. (adapted from Haber et
al. 2002)

Time-distance methods make direct use of the surface Doppler signal (rather than its power
spectrum) and determine flow velocities by measuring the time required for wave packets to travela fixed distance across the solar surface (e.g., Duvall et al. 1993; Kosovichev 1996). A backgroundflow produces an anistropy in the travel times for waves propagating in opposite directions, andthermal heterogeneities, which affect the sound speed, produce isotropic travel time shifts. Travel
time measurements are typically obtained by constructing the cross-covariance function and fitting
it with a Gaussian wavelet whose phase is related to the time lag between the two wave signals (e.g.,
Kosovichev et al. 1997; Figure 1.8). By Fourier filtering the Dopplergrams to isolate waves withspecific radial orders and phase speeds, travel-time maps of waves with varying degrees of depth
penetration can be constructed. Point-to-point correlations are often noisy and so travel times are
typically constructed using some degree of spatial averaging, either between a point on the surface
and a surrounding annulus (Duvall et al. 1996, 1997), or by averaging over all points that have thesame horizontal separation lying within some common region of the Sun (Featherstone et al. 2004).
These travel times may then be inverted to discern the subsurface flows (e.g., Gizon & Birch 2002;
Birch et al. 2004). A detailed review of the time-distance technique is provided by Gizon & Birch

18

Figure 1.8: Example of the time-distance cross covariance function. (a) Measured point-to-pointcross covariance as a function of spatial lag (shown here as angular separation in degrees) andtime lag. Broad, ridged bands of power are apparent. The slope of these bands indicates the timerequired for a wave packet to travel between two points separated by the indicated angular distance.
Different bands are denoted by colored and numbered lines and correspond to waves undergoing
different numbers of skips when traveling between points. (b) Geometry of the ray paths taken by
waves when creating the different bands of power in the cross covariance function. (Adapted from
Gizon et al. 2010)(2005).
We note that ring analysis and time-distance methods are by no means the only local helioseismic methods. One of the early applications of local helioseismology employed a techniqueknown as Hankel analysis (Braun et al. 1987), which examines differences between waves propagating into and out of some point on the solar surface. By looking at the power in incoming andoutgoing waves, Braun et al. (1987) found that sunspots were strong absorbers of acoustic power.
Hankel analysis has found wide application in the measurement of phase shifts and mode mixing
induced by sunspots (e.g., Braun, et al. 1992, Braun 1995). Another form of local helioseismology
is acoustic holography (Lindsey & Braun 1997), which uses the wavefield at the solar surface toestimate properties of the wavefield in the interior and has also been applied to the measurement

19of phase shifts in waves induced by magnetic fields (e.g., Lindsey & Braun 2005a,b) The mostnovel application of acoustic holography (and indeed local helioseismology) is its application to farside imaging of active complexes in the Sun (e.g., Lindsey & Braun 2000; Braun & Lindsey 2001;
González Hernández et al. 2007), making it a powerful predictive tool for active region emergence.

1.7.2

Meridional Circulations

One of the early applications of local helioseismic techniques to global-scale motions wasto the study of meridional circulations in the upper convection zone. Known through Doppler
measurements to possess a surface velocity of 20 m s−1 (Hathaway et al. 1996), these poleward
circulations are thought to be a key mediator of angular momentum transport in the solar convection
zone. There, along with Reynolds stresses resulting from the turbulent convective motions, themeridional circulations likely play a pivotal role in the establishment and maintenance of the solardifferential rotation (e.g. Brun & Toomre 2002). Rotational shearing in the tachocline, associatedwith the transition from the differentially-rotating convection zone to the solid body rotation of theradiative interior, is thought to play central role in the development of strong toroidal magnetic
fields by the solar dynamo. Measurements of meridional circulations are thus useful both for
informing 3-D numerical simulations of dynamo action in the convection zone (e.g., Brun, Miesch& Toomre 2004). Moreover, such measurements also provide important input to flux transport
dynamo models (e.g., Dikpati & Charbonneau 1999).
Initial measurements of the meridional circulations using local helioseismic techniques wereshown through subsequent time-distance measurements to remain relatively constant in strength
down to a depth of 26 Mm (Giles et al. 1997). Further studies using a variety of local helioseismic
approaches have confirmed this finding (e.g., González, Hernández et al. 1999, 2000; Duvall &Gizon 2000; Zhao & Kosovichev 2004). Temporal variation of the meridional circulations has beenassessed through time-distance (Giles 1999) and ring-analysis studies (Haber et al. 2002). As shownin Figure 1.9, these studies revealed the presence of a deep countercell of meridional circulationin the northern hemisphere. This submerged cell lacked a southern counterpart and varied in20

Figure 1.9: Meridional circulations as measured using ring-analysis helioseismology. Flows havebeen averaged in longitude and in time for the two periods spanning (lower panel) 14 April 1997
through 12 July 1997 and (upper panel) 7 March 1999 through 26 May 1999. Meridional circulations
are poleward within each hemisphere, and the strength of the circulations is relatively constantin depth during 1997. Variability of the circulations is visible in 1999 where the presence of asubmerged countercell is apparent in the northern hemisphere. (Adapted from Haber et al. 2002)

strength and sign with respect to the overlying shallow circulation over the course of the years
1998-2001. Subsequent studies (González Hernández et al. 2006) using more deeply penetrating
modes also found a submerged countercell, and reported that similar cells appeared periodicallywith the changing B0 -angle of the Sun. The existence of such multiple circulatory cells within asingle hemisphere might thus be an observational artifact (see also Gizon & Birch 2005).
Meridional circulations also vary over the course of the solar cycle. Upon subtraction of thetemporal mean from the meridional circulation profile, a surface inflow into the active latitudes
becomes visible. This 10 m s−1 inflow pattern migrates equatorward, along with the magnetic
activity belts, over the solar cycle (Haber et al. 2002; Gizon 2003; Zhao & Kosovichev 2004;
González Hernández et al. 2008; Basu & Antia 2010). Time-distance observations extending togreater depths suggest that these inflows transition to similarly migrating diverging flows below
about 18 Mm in depth (Chou & Dai 2001; Beck et al. 2002). The correlation between these21
variations in meridional circulations and the active region belts has led to investigations into the
influence that active-region related flows might have on the meridional circulations and differential
rotation.

1.7.3

Flows Around Active Regions

Local helioseismology has been extensively applied to the probing of the flows into andaround active regions, since the shearing of magnetic field lines by the flow is likely responsible fordriving the magnetic reconnection associated with energetic events such as flaring and coronal massejections. Both time-distance and ring analysis have shown that on large-scales (30◦ in extent),near-surface flows tend to converge on active regions with speeds of about 50 m s−1 (Gizon etal. 2001; Haber et al. 2001). When depth inversions are applied, these inflows are shown totransition to outflows of comparable strength around a depth of 10 Mm as in Figure 1.10 (Haber et
al. 2002; Zhao & Kosovichev 2004; Komm et al. 2005). Komm et al. (2004) have suggested thatthe divergence of large-scale flows around active regions can be coupled with mass conservation tocompute the vertical velocity of the flow and find 1 m s−1 downflows down to a depth of about 10
Mm beneath active regions and upflows below this depth.
Shearing flows tend to pervade active region complexes. Using ring analysis to study flows
beneath a quiescent filament, Hindman et al. (2006) reported that near-surface convection (onscales of about 40 Mm) underlying the filament tended to generate an average shear about theneutral line between its magnetic footpoints. On larger scales, Coriolis forces induce vorticity inthe surface inflows (cyclonic in the northern hemisphere) that transition to weakly anticyclonic flowin the surrounding quiet Sun (Komm et al. 2007; Hindman et al. 2009). Within an active complex,
the outflows associated with individual sunspots tend to be anticyclonic as well (Hindman et al.
2009). Moreover, the active regions themselves revolve faster about the rotation axis of the Sunthan their quiet Sun surroundings by about 20 m s−1 (e.g., Zhao et al. 2004; Hindman et al. 2009;
Komm et al. 2009).
Observations of these shearing flows in active regions have elicited much interest in using

22

Figure 1.10: Surface inflows into active regions transition to outflows at depth. OLA inversionof 16◦ ring measurements made around active region NOAA 9433 on 23 April 2001 are shown at
depths of (a) 7 Mm and (b) 14 Mm. Underlaid is the magnetogram for that day, with opposite fieldpolarities indicated in red and green. (c) Horizontal flows into and away from the active region as
a function of depth and latitude, averaged over the longitudes 142.5◦ -157.5◦ and the time period
23 - 27 April 2001. Inflows transition to outflows around a depth of 10 Mm. (Adapted from Haber
et al. 2004)

these measurements as a predictive tool for flaring and CMEs. The vertical velocity of active region
flows, along with its curl, may be used to construct the associated kinetic helicity. Subsequent ringanalysis studies have sought to use the vertical velocity and kinetic helicity to develop predictivemodels for magnetic events. The observations of Komm et al. (2008, 2009) hint that the appearanceof organized upflows might precede the emergence of large-scale magnetic flux. The measured the
kinetic helicity has also been shown to correlate well with both the level of a magnetic activity (i.e.
levels of high activity have high kinetic helicity), and the intensity of the X-ray flares arising fromthe active region (Komm et al. 2005, 2007; Mason et al. 2006)

23
1.7.4

Active Region Influence on Global-Scale Flows

Given the large-scale inflows into and circulations around active regions, it is natural to ask
whether the observed variations in the meridional circulation (Figure 1.11a) over the solar cycle are
due to migrating active region inflows or are they due to the presence of additional meridional cells(as posited by Zhao & Kosovichev 2004)? Moreover, what effect might active region flows haveon the so-called solar “torsional oscillations?” These latitudinal variations in the solar rotation
rate are characterized by alternating bands of 10 m s−1 slow and fast rotation (Figure 1.11b) that
migrate equatorward in similar manner to the bands of solar activity. The torsional oscillations
were first observed by Howard & Labonte (1980) through Doppler measurements of the surface
rotation rate, and examined in recent detail by Ulrich (2010). They have been studied throughinversion of global-mode frequency splittings (e.g. Kosovichev & Schou 1997; Schou 1999; Howe etal. 2000) and through local helioseismic studies using time-distance (Giles et al. 1998; Beck et al.
2002; Zhao & Kosovichev 2004) and ring analysis (Basu et al. 1999, Basu & Antia 2000; Haber et
al. 2000, 2002). In the vicinity of the active regions, these torsional oscillations manifest as regionsof fast rotation equatorward of the sunspot latitudes and slow rotation on the poleward side.
The driving of torsional oscillations has been attributed to a variety of phenomena, including
Lorentz forces from migrating dynamo waves, changes in the balance of Reynolds and Maxwell
stresses, and the suppression of turbulent viscosity by active regions (see e.g., Shibahashi 2004 fora review). One particularly interesting model is that of Spruit (2003) which proposes that thetorsional oscillations originate through downdrafts driven by radiative cooling associated with thesmall-scale magnetic fields of active regions. Coriolis forces acting on these downdrafts produce a
cyclonicity in agreement with observations which (in this model) is then responsible for the fasterrotation equatorward of the active belts and the slower rotation on the poleward side. Similarly,
the variation of meridional circulation over the solar cycle is appropriately reproduced.
Gizon & Rempel (2008) demonstrated that incorporating such a surface cooling term into aflux transport model can yield meridional and zonal flows in rough agreement with those measured

24

Figure 1.11: Torsional oscillations and meridional flow variations over the solar cycle. Shown are
the residuals for (a) meridional circulation and (b) zonal flows obtained after subtracting the timeaveraged component of the flow. Meridional flows were measured for the years 1996-2006, and zonal
flows during 1996-2002. Flows are extrapolated into the future to better illustrate the patterns
present. Thick black curves indicate the mean latitudes of magnetic activity. The migration of
the zonal flows and meridional circulations displays a clear linkage to the migration of the active
latitude belts. (adapted from Gizon & Rempel 2008; Gizon et al. 2010)

through time-distance. Other approaches to understanding this problem have involved analyses oflongitudinal means obtained which incorporate the active regions and those obtained when regionsof magnetic activity are masked out. Using f -mode time-distance, Gizon (2003) found that active
region flows alone do not fully account for the torsional oscillations (i.e. torsional oscillations occur
in quiet Sun as well), but can account for the variations observed in meridional circulations overthe solar cycle. However, ring-analysis measurements of quiet Sun meridional circulation obtainedusing similar magnetic masking indicate that active region inflows cannot fully account for thevariation in meridional circulation (González Hernández et al. 2006). The disentanglement of the
global-scale flows from active-region flows is thus tricky and the effects of the former on the latterremain unclear, particularly with respect to the meridional circulations.

25
1.7.5

Flows Around Sunspots

The emergence and migration of sunspots at the solar surface, the tell-tale signature of the
22-year solar cycle, has engendered much interest into their structure and maintenance. Of primaryinterest is whether sunspots are monolithic flux structures extending deep into the convection zone,
or whether they are coherent only in the shallow layers, transitioning to a much more diffuse fibril
structure below a depth of 1-2 Mm (e.g., Parker 1979, Schüssler 2005). Moreover, do plasma flows
in the vicinity of sunspots possibly help stabilize these magnetic structures, allowing them to persist
for days and even weeks? Theory and simulations have long suggested that sunspots must possessa stabilizing, convergent downdraft (e.g., Parker 1979, 1992; Hurlburt & Rucklidge 2000; Hurlburt& DeRosa 2008), but recent 3-D simulations by Rempel (2009) display systematic outflows withdepth.
Early measurements revealed the presence of surface outflows around sunspots, now termedmoat flows, visible in the movement of surface magnetic features and later through Doppler measurement at the solar surface (Sheeley 1969, 1972; Brickhouse & Labonte 1988). These moat flows
were found to be on the order of 1 km s−1 and are thought to be driven by pressure gradients resulting from heat buildup below the sunspot where convection is suppressed (e.g., Nye et al. 1988;
Solanki 2003). The advent of local helioseismic techniques has enabled further probing of the subsurface nature of these moat flows. Using f -mode time-distance methods, Gizon et al. (2000) wereable to confirm the presence of the moat flow in the upper 2 Mm of the convection zone. Deeper
probing using the solar p-modes has enabled the persistence of the moat flow with depth to beprobed using a variety of helioseismic techniques (e.g., Lindsey et al. 1996; Sun et al. 1997; Braun& Lindsey 2003; Gizon et al. 2009; Featherstone et al. 2011). There are differences concerning thestrength and depth to which the moat flows extend, but generally these approaches concur thata radial outflow on the order of 100-200 m s−1 persists around a sunspots to a depth of about 10
Mm.
A notable exception to this trend are the p-mode time-distance results of Zhao & Kosovichev

26(2001, 2003) and Zhao et al. (2010) which suggest that sunspots possess convergent downdrafts inthe upper 3 Mm of the convection zone. The disagreement between such results and other helioseismic studies, particularly the p-mode measurements of Gizon et al. (2009), is puzzling. However,some of the earlier time-distance studies using p-modes did indicate downdrafts in sunspots on theorder of 1-2 km s−1 in the upper 1 Mm (Duvall et al. 1996; Kosovichev et al. 1996). As noted in
Gizon et al. (2009), time-distance measurements are particularly sensitive to the filtering process
wherebye waves with a certain radial order and phase speed are isolated before measuring travel
times. Different filtering will undoubtedly affect the depth-averaging properties of the different
measurements, possibly leading to inconsistent results.

1.8

Achieving 3-D Inversions of Ring-Analysis Measurements
A number of the insights into the dynamics of the upper convection zone have been informedby ring-analysis methods. One goal of this thesis is to improve the capability of ring analyses
to assess small-scale and intermediate-scales of subsurface flows. Owing to the lack of fine-scale
resolution available at depth in their current inversion methods, ring-analysis studies yield little
information on the nature of the deep flows around sunspots. Ring analyses have typically been
performed on square regions of the Sun that are either 16◦ or 2◦ (in heliographic angle) in size. In thecase of 16◦ measurements, the data are inverted to yield the underlying flow field with depth, anda horizontally uniform flow is assumed, effectively making these 1-D inversions. Ring analyses on2◦ regions afford a much finer resolution, but have been restricted to fitting f -mode measurements
only, with no inversion performed. These higher resolution measurements are averaged in horizontal
wavenumber to yield a representative average over the upper 2 Mm of the convection zone.
Applications of 2◦ surface measurements to the study of sunspots within active regions andshearing flows beneath filaments (Hindman et al. 2006, 2009) hint at an unexploited potential
for ring analyses to probe these flows in similar detail at depth. In principle, such benefits may
be gained by performing 1-D inversions on 2◦ p-mode measurements, similar in nature to the16◦ inversions. In practice, noisy power spectra make the fitting of modes penetrating more deeply

27
than 5 Mm difficult for 2◦ analysis regions. However, when analysis regions with a somewhat
coarser horizontal resolution of 4◦ are employed, Doppler shifts for modes reaching down to about
10 Mm in depth can be reliably measured.
Ideally a flow map of the subsurface flow should be constructed by self-consistently combining
measurements made using different analysis region sizes. The region of the solar disk covered by asingle 16◦ measurement will also be covered by numerous higher-resolution measurements. These2◦ and 4◦ measurements can provide a much more detailed description of the near-surface layers,
but say nothing about the flow at depths beyond about 10 Mm. Combining measurements fromdifferent analysis regions can be accomplished through a 3-D inversion of the data but it requires
the development of 3-D kernels that describe both the vertical and the horizontal sensitivies of thevarious measurements to a flow.
We have developed a 3-D inversion method capable of combining multiple analysis region
sizes in such a fashion and have calculated the 3-D sensitivity kernels necessary for the inversion ofexisting MDI data. This inversion package, while tested using MDI data, will be implemented into
the HMI data pipeline in spring 2011. Collectively termed the ARRDI (Adaptably-Resolved Ring
Diagram Inversion) package, these programs enable the self-consistent inversion of ring-analysis
fits achieved using an arbitrary combination of modes and analysis regions sizes. We will detailthe development of ARRDI beginning in Chapter 5 and turn to the first applications of ARRDI
to flows around sunspots and active regions in Chapter 8. We begin now with an exploration of
convective dynamos in the deep cores of A stars.

Chapter 2

Primordial Magnetic Fields with Core Dynamos in A-type Stars

We begin with an exploration of magnetism and convection in the cores of A-type stars. Observations of magnetism in the peculiar A-type (Ap) stars indicate the presence of concentrationsof magnetic field that extend over large portions of the stellar surface. These patches of magnetic
field are thought to coincide with locations where a primordial magnetic field, originating in theformative collapse of the star, extends through the stellar surface. The following three chapters
discuss how such strong, organized magnetic fields may influence the dynamo achieved deep within
the cores of these stars and are largely a restatement of the work contained in Featherstone et al.
(2009)1 . As the principal architect of the text, I was also responsible for conducting the numerical simulations and performing their subsequent analyses. Throughout this study, my coauthors
provided guidance concerning the exploration of the strong dynamo states realized through thesesimulations.
We devote this chapter to a discussion of the observations of surface magnetic fields in the
Ap stars and to theories describing their origin. We close with a discussion of our methodology for
examining the interaction between a primordial magnetic field and a core dynamo. In Chapter 3,
we focus on how the overall magnetic and kinetic energy balance is modified in the dynamo staterealized when a primordial magnetic field and core dynamo are modeled together. The various
elements of the flows and magnetic field structures that characterize the super-equipartition statethat is realized, in which the magnetic energy can be tenfold greater than the convective kinetic
1

Featherstone, N.A., Browning, M.K., Brun, A.S. & Toomre, J., 2009, “Effects of Fossil Magnetic Fields on
Convective Core Dynamos in A-type Stars”, ApJ, 705, 1000

29
energy, are analyzed there. Chapter 4 is devoted to a more detailed examination of the processes
responsible for the generation and maintenance of such strong magnetic fields and closes withreflections on the significance of achieving such super-equipartition dynamo behavior.

2.1

Ap Stars and Magnetism
The unusual chemical abundances of peculiar A-type stars were first noted by Maury (1897),

who discovered a strong Si II doublet in the spectrum of α2 CVn. Subsequent observations of thisand similar stars over the last century have shown these stars to exhibit strong and variable spectral
lines (relative to solar values) in Si and certain other elements (e.g., Sr, Hg). Most, if not all, of the
Si and Sr-Cr-Eu pecularity classes possess equally variable and unusually strong magnetic fields(Babcock 1947, Mestel 1999). These magnetic fields are generally thought to be of primordial
origin, resulting from the diffusive evolution of the magnetic field threading the initial molecular
cloud. However, convection in the cores of these stars may also play a significant role. Our interests
concern the nature of the interaction between such a primordial magnetic field and a core-dynamo
harbored within the convective interior of an A-type star.

2.1.1

Observations of Magnetic Fields

Observations of magnetic fields in Ap stars are carried out primarily through measurementsof the longitudinal (line of sight) field, deduced by measuring the wavelength shifts between spectral
lines of opposite circular polarization (Mathys 2001). Typical longitudinal field strengths for themagnetic Ap stars are a few hundred Gauss, but field strengths ranging from 20,000 G down to alower threshold of ∼ 300 G have been observed (Aurière et al. 2007). Some inferences about thenature of the structure of magnetic fields on these stars has been made through complementarymeasurements of the magnetic field modulus (also termed “the surface field”). However, Doppler
broadening makes such measurements possible only for the slowest rotators (Mathys 2001). Observations of the slow rotators suggest a lower cutoff for the mean modulus of about 3,000 G, well
above the detectability limit (Mathys et al. 1997, Freyhammer et al. 2008). For stars on which

30
both the longitudinal field and the field modulus have been measured, the field modulus is typicallysimilar in magnitude to the longitudinal field. A comparable longitudinal field and field modulus
hint at the presence of a large-scale ordered field in these slow rotators. A spatially intermittent
field would have a much smaller ratio of the longitudinal field to the modulus, as is the case for
the Sun, where variations in the line-of-sight field across the solar surface tend to cancel each other
out.
Field strengths for the Ap stars, when variable, are observed to vary at the stellar rotation
rate (Deutsch 1956, Preston 1971). The Ap stars appear to preferentially rotate more slowly than
the normal A-stars, with the possible exception of the early A-stars (Abt 2000, Royer et al. 2007).
Within the class of Ap stars, however, there is little or no correlation between magnetic fieldstrength and rotation rate (e.g. Kochukhov & Bagnulo 2006). Such a lack of correlation is contrary
to what one might expect from fields generated by a core dynamo. Within the class of Ap stars,
rotational periods range from decades to a fraction of a day, with typical periods on the order of afew days (e.g. Borra et al. 1982, Mestel 1999).
The variability and apparent large-scale organization of the magnetic fields observed in Ap
stars has led to oblique rotator models (ORM) (Stibbs 1950, Mestel 1999). The ORM envisage alarge-scale magnetic dipole, inclined to the rotation axis, frozen into the atmosphere of the star,and dragged along with the star as it rotates. This model is highly idealized and modifications toa simple dipole are often required. The modulus measurements of Mathys et al. (1997) suggest
a significant deviation from a purely dipolar geometry, with only one minimum and maximum
present during a rotational cycle (as opposed to the two expected for a dipolar field). More recently
Kochukhov et al. (2004) showed that the spectrum of 53 Cam was well fit only after employing
higher order multipoles (spherical harmonic degree ℓ up to 15).

2.1.2

Possible Primordial Origin of Magnetism

The apparently static nature of the magnetic fields of Ap stars has led to the developmentof the fossil field theory. As originally pointed out by Cowling (1945) in the context of the Sun,

31
a typical Ohmic decay time for a global-scale stellar magnetic field is longer than the lifetime ofthe star. The observed magnetic fields may then be a remnant of an interstellar field threading
the cloud of gas from which the star formed. Recent observational surveys suggest that Ap stars
exhibit markedly less magnetic flux as they near the end of their main-sequence lifetimes, possibly
due to diffusive decay or reorganization of magnetic fields in the radiative zones of these stars(Kochukov & Bagnulo 2006, Landstreet et al. 2008). The fossil theory may also help explain theslow rotation of Ap stars (relative to their nonmagnetic counterparts) through magnetic brakingof the star earlier in its lifetime (Mestel 1999). Such braking has been observed in magnetic Aestars (the pre-main sequence progenitors of Ap stars) which exhibit rotation rates about six times
slower than the non-magnetic Ae stars (Alecian et al. 2008).
The structure of a typical fossil field has been a subject of speculation for some time. On
the basis of stability arguments, these fields have generally been thought to be comprised of someform of linked poloidal and toroidal fields. Purely toroidal configurations have been shown by
Tayler (1973) to be unstable to perturbations of azimuthal order m = 1, with the perturbations
growing on the time scale needed for an Alfven wave to travel along the magnetic field line. Purely
poloidal field has been shown to be unstable to non-axisymmetric perturbations (e.g. Wright 1973,
Markey & Tayler 1973, Flowers & Ruderman 1977). The analysis of Braithwaite (2007) suggeststhat rotation decreases the rate of growth of instability, but does not prevent the onset. Likewise
Brun (2007) has studied the role of rotation on the onset of the poloidal instability and foundthat the most unstable wavenumber m increases with the rotation rate. A stable magnetic field
configuration in the radiative zone must necessarily contain more twist than purely poloidal or
toroidal configurations afford. Prendergast (1956) showed that such a linked, poloidal-toroidal
configuration satisfied the criteria for equilibrium, but stopped short of proving the stability ofsuch a field. More recently, Duez & Mathis (2009) have derived a linear stability analysis of themagnetic field configuration in a barotropic star, confirming the mixed toroidal-poloidal nature of
the stable field.
The numerical simulations of Braithwaite & Spruit (2004) and Braithwaite & Nordlund (2006)

32
for a nonrotating star have shown that a linked poloidal-toroidal field configuration is indeed stableon time scales commensurate with those of the star’s evolution, and appears to be the preferred
equilibrium condition for an initially random field. They find that a randomly oriented fossil field,
present from the time of the star’s appearance on the main sequence, will slowly relax into atwisted torus shape. The diffusion of this torus through the radiative zone, into the nonconducting
atmosphere of the star (which cannot support a twisted field) results in the appearance of a largely
poloidal field at the stellar surface. They estimate this diffusion time to be around 2 × 109 years, atime somewhat longer than the main sequence lifetime of the star. The appearance and subsequentdiffusion of this fossil field may also help to explain the observed middle-aged nature of the magnetic
Ap stars (Kochukhov & Bagnulo 2006).

2.1.3

Dynamic Origins for Fields

The cores of A-type stars (roughly the inner 15% by radius), along with a thin shell at thesurface, are known to be convective, and such core convection has long been suspected of achieving
dynamo action (Krause & Oetken 1976). Brun, Browning & Toomre (2005 hereafter BBT05) haveinvestigated this possibility through 3-D numerical simulations and found that the cores of A-stars
sustain vigorous dynamo action, with typical field strengths reaching near equipartition values withthe convective kinetic energy. Core dynamos may have bearing on the surface fields observed in
Ap stars through the production of magnetically buoyant structures. This possibility has been
investigated through simulations by MacGregor & Cassinelli (2003) for more massive O and B
stars, demonstrating that flux tubes may rise to the stellar surface in a time significantly less than
the star’s main sequence lifetime. However, MacDonald & Mullan (2004) have pointed out thatthe inclusion of realistic compositional gradients into the model of MacGregor & Cassinelli would
severely hinder the rise of flux tubes, ultimately requiring a very high field strength if a magnetic
structure is to reach the stellar surface. Given that these fields are assumed to be produced bydynamo action in the core, their strengths would need to be well in excess of equipartition values,
possessing magnetic energies at least ten times greater than the convective kinetic energy.

33
Convection may also play a role in the early evolution of these stars’ magnetic fields. The
fully convective Hayashi phase that A-stars are thought to undergo on their journey to the main
sequence holds the potential for magnetic dynamo action. Fully convective dynamos have beenfound, for instance by Dobler et al. (2006) and in 3-D modelling of M-type stars (Browning &Basri 2007; Browning 2008). Surveys of Ae stars by Alecian et al. (2008) reveal the existence ofobservable magnetic fields in these stars following the Hayashi phase (but prior to the developmentof a convective core). Either these fields must have been present prior to a Hayashi phase andsurvived, or were generated by the convective motions present during that epoch. Dynamo activity
during a Hayashi phase thus seems likely and may augment or alter any primordial magnetic fields.
Separately, the field instabilities noted in §1.2 could lead to dynamo action through interaction
with prominent differential rotation, as suggested by Spruit (2002), though some complications withthis mechanism have been raised by Zahn et al. (2007). The dynamo scenario envisioned by Spruitis unlikely to be operating in main-sequence Ap stars, since the surfaces of these stars are observed
to be rotating as solid bodies. This mechanism may, however, have had a role in field generation
at earlier epochs.

2.1.4

Interaction of Core Dynamo with Fossil Field

While a primordial magnetic field seems to be the likely source for the surface magnetism
observed on Ap stars, convection may also be implicated. Whether during the Hayashi phase earlyin the star’s lifetime, or during the main sequence phases, convection and associated dynamo action
may have a role in the magnetic fields that are observed. This leads to several basic questions. Whatis the nature of the interaction between a primordial field and one generated contemporaneously
by a core dynamo? Does a strong fossil magnetic field help or hinder dynamo action, and does the
dynamo have any significant effect on the configuration of the fossil field? Kinematic modeling by
Moss (2004) suggests that a sufficiently strong exterior field could hinder the dynamo action in thecore. Although instructive, this axisymmetric model does not capture the complex 3-D processes
likely to be occurring within the core.

34
We are thus motivated to turn to detailed modeling of the coupling of a fossil field andconvective core dynamo within an A-type star. The nonlinear 3-D simulations used here allow us
to assess both the growth and equilibration of the magnetic fields, and also the feedback of suchfields on the nature of flows within the core. Our results provide some insights into the interaction
between an exterior fossil field and the stellar core dynamo, demonstrating that the presence of a
fossil field holds out possibilities for enhanced dynamo action.

2.2

Simulation Elements

2.2.1

Anelastic MHD Approach

Our numerical simulations are carried out using the anelastic spherical harmonic (ASH)code which is described in detail in Clune et al. (1999) and in Brun et al. (2004). ASH isa pseudo-spectral code designed to perform efficiently on massively parallel supercomputers andsolves the 3-D magnetohydrodynamic (MHD) equations in a rotating spherical shell. The anelastic
approximation assumes that fluid motions are subsonic and that perturbations to thermodynamic
variables are small compared to their mean, horizontally-averaged values at a given depth in the fluid(Gough 1969, Gilman & Glatzmaier 1980). Within the anelastic approximation, the thermodynamic
variables are linearized about their spherically symmetric and evolving mean state with density ρ̄,pressure P̄ , temperature T̄ , and specific entropy S̄. Fluctuations about this state are denoted as
ρ, P , T , and S. In the uniformly rotating reference frame of the star, the MHD equations are
expressed as∇ · (ρ̄v) = 0,(2.1)

∇ · B = 0,

(2.2)

35




∂v
ρ̄
+ (v · ∇)v + 2Ω0 × v = − ∇P + ρg∂t1
+ (∇ × B) × B − ∇ · D − [∇P̄ − ρ̄g],
4π

(2.3)

∂S
+ ρ̄T̄ v · ∇(S̄ + S) = ∇ · [κr ρ̄cp ∇(T̄ + T )∂t4πη+ κρ̄T̄ ∇(S̄ + S)] + 2 |j|2c 1
+ 2ρ̄ν eij eij − (∇ · v)2 + ρ̄ǫ,
3

(2.4)

∂B
= ∇ × (v × B) − ∇ × (η∇ × B).
∂t

(2.5)

ρ̄T̄

The velocity v expressed in spherical coordinates is v = (vr , vθ , vφ ) relative to a frame rotating atconstant angular velocity Ωo , g is the gravitational acceleration, B = (Br , Bθ , Bφ ) is the magnetic
field, j = c/4π(∇ × B) is the current density, cp is the specific heat at constant pressure, κr is theradiative diffusivity, η is the effective magnetic diffusivity, and D is the viscous stress tensor given
by



1
Dij = −2ρ̄ν eij − (∇ · v)δij ,
3

(2.6)

where eij is the strain rate tensor, and ν and κ are effective eddy diffusivities. We have chosen to
represent the energy generation by nuclear burning of the CNO cycle with a volumetric heating termρǫ (see §2.2.3 below for details). This set of equations is closed by assuming the thermodynamic
fluctuations satisfy the linear relations

P
T
P
S
ρ
=
− =
− ,
ρ̄
cp
P̄
T̄
γ P̄

(2.7)

P̄ = Rρ̄T̄ ,(2.8)

assuming the ideal gas law

where R is the gas constant. In this MHD context, the anelastic approximation also filters out thefast magneto-acoustic modes. Alfven modes, as well as the slow magnetosonic modes are retained.

36
The divergence-free nature of the mass flux and of the magnetic field is maintained through theuse of a toroidal-poloidal decomposition wherebyρ̄v = ∇ × ∇ × (W er ) + ∇ × (Zer ),(2.9)

B = ∇ × ∇ × (Cer ) + ∇ × (Aer ).

(2.10)The streamfunctions W and Z, as well as the magnetic potentials C and A, are functions of
all three spatial coordinates and time; er is the radial unit vector.
This set of equations in the six variables W , Z, C, A, S, and P requires 12 boundary conditions to be well posed. We have chosen to impose the following boundary conditions throughout
the simulation:(1) Impenetrable top and bottom surfaces in the deep shell: vr = 0(2) Stress-free top and bottom: (∂/∂r)(vθ /r) = (∂/∂r)(vφ /r) = 0(3) Constant entropy gradient at top and bottom: ∂ S̄/∂r = constant, ∂S/∂r = 0(4) Perfect conductor top and bottom: Br = (∂/∂r)(Bθ /r) = (∂/∂r)(Bφ /r) = 0
These conditions ensure that no mass or angular momentum is lost from the system via either
radial momentum flux or torques arising from viscous or magnetic stresses at the upper and lower
boundaries of the deep spherical shell being studied. Emergent flux through the top and the bottom
remains constant in time, but differs between the two boundaries to account for energy generationin the core. The leakage of magnetic energy from the domain is prevented by forcing the Poynting
flux to vanish at the boundaries. Other boundary conditions, especially for the magnetic fields,
could be imposed, but at this stage these appear to be the most neutral, particularly at the lower
boundary where strong magnetic fields are generated.

37
2.2.2

Implementation of ASH

Seeking to deal with full spherical global domains, ASH is a large-eddy simulation (LES) code,with subgrid-scale (SGS) descriptions for dynamics occurring on scales below the spatial resolutionof the simulations. We have chosen to focus on the larger scales of flow, believing that they are
most likely to be responsible for establishing the mean properties of core convection and dynamo
action. In the solar context, this approach has been reasonably successful in making contact withhelioseismic deductions of the solar interior differential rotation (Brun & Toomre 2003, Miesch,
Brun & Toomre 2006, Miesch et al. 2008), and has been shown in the A-star context (BBT05) to
produce sustained dynamo action. Here we treat the SGS terms most simply as enhancements to
the kinematic viscosity as well as to the thermal and magnetic diffusivities. Much as in BBT05,
we have defined our eddy viscosity and diffusivities to be a function of the mean density alone(and hence radius), independent of horizontal position or time. Specifically, we have taken thesequantities to be proportional to ρ−1/2 .
ASH’s numerical implementation involves expanding all variables in spherical harmonics
Ylm (θ, φ) in the horizontal directions and Chebyshev polynomials in the radial direction. To gain
higher resolution at the interface between the convective and radiative zones, and thus to better
resolve the penetrative convection occuring in this region, we have employed a stacked Chebyshev
scheme. In doing so, we have split the computational domain into two separate regions (radially)and performed separate Chebyshev expansions for each region. The spherical harmonic expansion
is truncated at degree ℓmax , with all azimuthal orders m retained in a triangular truncation, ensuring that we have uniform resolution over spherical surfaces. For our simulations, we have takenℓmax = 170, corresponding to 256 mesh points in the latitudinal direction (Nθ ), with the longitudinal mesh having Nφ =2Nθ , and in radius Nr = 49 + 33 = 82. The time evolution of our simulations
is computed using an implicit, second-order Crank-Nicolson scheme for the linear terms and an
explicit second-order Adams-Bashforth scheme for the advective, Coriolis and Lorentz terms. The
computational demands of these calculations are substantial, and thus the ASH code has been

38
optimized to run on massively parallel supercomputers, using the message passing interface (MPI)
to communicate between different computational nodes.

2.2.3

Modeling the A-type star

We have chosen to model the inner 30% by radius of a main-sequence A-type star of 2M⊙
rotating at four times the mean solar rate of Ω0 = 414 nHz. The mean rotation period is thus
about seven days. The emerging luminosity of this star is 19L⊙ . The inner 2% of the star by radius
has been excluded to avoid the coordinate singularity at r = 0. This model consists of a convective
core occupying the inner 15% by radius of the star with an overlying radiative zone extending tothe outer radius of our computational domain. The density contrast across the full domain is 21.5,
and that across the convective core is 2.04.
All models presented here were initialized using the statistically mature case C4m of BBT05.
The convection and associated dynamo action involved complex flows and magnetism that were
continuously evolving and rebuilding structures, but the overall system had equilibrated in terms of
its time-averaged properties such as energies and differential rotation. Case C4m was evolved froma one-dimensional stellar structure model (at an age of 500 My) computed with the CESAM stellar
evolution code (Morel 1997) using realistic microphysics and employing a classical mixing-length
formalism to describe the convection. A necessary modification to this model was the softening of
the entropy gradient in going from the convective core to the radiative zone. The stark contrast in
entropy between these two regions would otherwise result in the driving of high frequency internal
gravity waves which would severely limit the time steps. The principal input parameters are thesame for the simulations presented here as for their progenitor case. As in BBT05, the energy
8generation term ǫ was implemented as ǫ = ǫ0 T6 , with ǫ0 = 6.0 × 10−9 erg g−1 s−1 . Furthurmore, wehave taken the thermal diffusivity κ to be 4.0×1012 cm2 s−1 at the top of the domain. Similarly there
the SGS viscosity is ν = 1.0 × 1012 cm2 s−1 and the magnetic diffusivity is η = 2.0 × 1011 cm2 s−1 .
Thus our effective Prandtl number Pr=0.25, and our magnetic Prandtl number is Pm=5.0. These
values were chosen to achieve a high magnetic Reynolds number Rm at reasonable computational

39
cost and to facilitate comparison with BBT05. Using rms velocities and magnetic field strengths
at mid-core and a length based on the core depth, the Reynolds number Re for our primary case A
following equilibration is 136, and the corresponding magnetic Reynolds number Rm is 680. The
associated Rossby number for this case is Ro=6 × 10−3 .

2.2.4

Imposing a Fossil Magnetic Field

All simulations here have been initialized by superimposing external magnetic fields of varying
geometries onto the existing field structure of case C4m of BBT05 at a time corresponding to day
2,000 of that simulation. Our choices for the initial magnetic configuration in each case have been
motivated by stability considerations and the need to sort out the role that the toroidal and poloidal
components of the imposed magnetic field have in the overall response.
We have chosen a mixed poloidal-toroidal magnetic field (case A), a poloidal field (caseB), and a purely toroidal magnetic field (case C) for our initial external field configurations. To
disentangle the effects of spatial scale of the field, we have also examined external field configurations
with higher-order multipolar structure (lacking a toroidal component) in which each hemispherereceives zero net magnetic flux across the interface between the convective core and the radiativezone.
The toroidal field for cases A and C was defined as(r − r0 )2
Bφ = A exp −
2σ 2



,

(2.11)

where r0 lies along the equator at 0.225R, σ is taken to be 2r0 , and the amplitude A is taken to be
30 kG.
Our mixed field case A included the same functional form for Bφ but with an added poloidal
field consistent with a current threading through the center of the magnetic torus. The strength ofthe poloidal field was adjusted so that the ratio of energy in the poloidal field to that in the toroidal
field was 1:9, as suggested by the results of Braithwaite & Nordlund (2006). Placing such a mixed
field into C4m resulted in only a 10% increase to the total magnetic energy of that system. The

40
a

vr

b

-/+ 20Br-/+ 40c

B-/+ 59dve-/+ 20B-/+ 7

Figure 2.1: Flow and magnetic field properties of the core convection dynamo in the progenitor case
C4m. (a) Snapshot of radial velocity at mid-core (r=0.10R), showing some columnar alignment ofwavering rolls with rotation axis (upright). Regions of upflow are shown in light tones, downflows
in dark tones; ranges in m s−1 . (b) Companion snapshot of radial magnetic field at mid-core, withstronger field roughly coincident with downflows. Positive fields denoted by light tones, negativeby dark; ranges in kG. (c) Azimuthal magnetic field Bφ at mid-core with some extended ribbonlike structures running east-west. (d) Mean zonal velocity hvφ i in radius and latitude possessing acentral column of slower rotation. Prograde flow shown in red tones, retrograde in blue tones; ranges
in m s−1 . (e) Mean azimuthal magnetic field hBφ i, with complex structure evident and exhibiting
no preferred sense in either hemisphere. Positive field in red tones, negative in green/blue; ranges
in kG.

purely toroidal and purely poloidal cases were also studied for comparison. In those cases, the field
geometry and strength was identical to that component in the mixed field case. Due to stability
considerations, it is the mixed field case A that appears to be the most physically relevant. As
such, we have focused the bulk of our computational efforts and analyses on that case.

2.3Core Dynamo in Progenitor
Some key aspects of the progenitor dynamo simulation case C4m are summarized in Figure

2.1 illustrating the typical flows and fields achieved at one instant in time. These flows have beensampled at mid-core (r = 0.1R) and are displayed on spherical shells. In addition, we have sampled
the differential rotational profile and mean fields realized in this system (Figs. 1d and 1e), where
the azimuthal velocity vφ and magnetic field Bφ have been averaged over longitude and displayed
as a function of radius and latitude.
Case C4m of BBT05 involves highly time-dependent flows with complex, vigorous convectionthat penetrated substantially into the overlying radiative zone. These flows are characterized by41
convective motions spanning multiple scale heights that serve to couple widely separated parts ofthe computational domain. These convective motions give rise to a nearly adiabatically stratified
core slightly prolate in shape. The prolate shape owes to the greater effect of the Coriolis force
on radial velocities at the lower latitudes (see discussion in BBT04). Surrounding the convective
core is a further region of overshooting, roughly spherical in shape, which does not substantiallymodify the stable stratification of the radiative zone. A striking change in case C4m compared toits own hydrodynamic progenitor is the significant suppression of the differential rotation in thecore. The hydrodynamic case has a prominent column of slow rotation maintained primarily by
Reynolds stresses arising from the convection that transport angular momentum equatorward andradially outward. These counteract viscous stresses that transport angular momentum poleward
and radially inward.
The addition of a seed dipole magnetic field to that hydrodynamic case to yield case C4m,
as discussed in BBT05, resulted in a marked change to the differential rotation of the system. The
dynamo action achieved by this system produces Maxwell stresses that tend to transport angularmomentum poleward, opposing the Reynolds stresses and ultimately disrupting the differential rotation. There remains a weak column of slower (retrograde) rotation (Fig. 2.1). The persistent
dynamo action realized in this simulation yields nearly equipartition magnetic energies, with thoseabout 90% of the kinetic energy. The magnetic energy in this system is largely comprised of fluctuating (non-axisymmetric) fields, with the axisymmetric poloidal and toroidal fields contributing
little to the overall energy balance. Both the magnetic fields and flows in case C4m are highlyintermittent in time, with no evident preferred field polarity.
The addition of a twisted dipolar magnetic field, consistent with the expected geometry ofa primordial field, to case C4m leads to remarkable changes in the core dynamo. Magnetic fields
achieved in the core become super-equipartition in nature, with an associated energy roughly tenfold
that of the convective motions. We turn now to a detailed exploration of this dynamo state in thefollowing chapter.

Chapter 3

Surprising Implications of a Primordial Magnetic Field

3.1

Modified Dynamo Action
A preliminary assessment of the dynamo action achieved in our 3-D simulations is provided

by examining the evolving global energy balances achieved in each case. The energy associatedwith the kinetic energy (KE) is defined as
1
KE = (ρ + ρ)(vr2 + vθ2 + vφ2 ),
2

(3.1)

and magnetic (ME) energy asME =1(B 2 + Bθ2 + Bφ2 ).
8π r

(3.2)

The azimuthal velocity vφ is taken relative to the rotating frame. We find it useful to examine
a number of decompositions of ME. To facilitate contact between our work and previous dynamo
modeling, we are often interested in mean fields and flows and the fluctuations about those means.
We thus adopt the following notation for the fluctuating (non-axisymmetric) velocitiesv ′ = ((vr − hvr i)er + (vθ − hvθ i)eθ + (vφ − hvφ i)eφ .

(3.3)

Angular brackets denote averages in azimuth (longitude). The fluctuating magnetic field vector
B ′ is defined likewise. The magnetic energy is further decomposed into those portions associated with the mean toroidal fields (MTE), the mean poloidal fields (MPE) and the fluctuating
nonaxisymmetric fields (FME). We define these quantities as
MT E =

1hBφ i2 ,
8π

(3.4)

43

Figure 3.1: Superimposing a fossil field into the progenitor case to initiate case A. (a) Azimuthallyaveraged Bφ at the initial instant and (b) the companion poloidal field. White lines denote poloidal
field lines, and the underlay the poloidal field modulus. (c) Magnetic energy density shown volume
rendered, with the imposed field forming the outer torus of moderate field strength. Field lines
correspond to those of the fossil field only.

44

1(hBr i2 + hBθ i2 ),
8π

(3.5)

1
[(Br′ )2 + (Bθ′ )2 + (Bφ′ )2 ].
8π

(3.6)

MPE =
and
FME =

We also find it useful to distinguish between the magnetic energy in the fluctuating φ componentof the field (FTE) and that in the fluctuating r and θ components (FPE), with
1
(B ′ )2
8π φ

(3.7)

1
[(Br′ )2 + (Bθ′ )2 ].
8π

(3.8)

FTE =
and
FPE =

The superposition of our twisted fossil magnetic field into case C4m, involving the introduction of both toroidal and poloidal components, and denoted here as the mixed-field case A,
is shown in Figure 3.1 at the outset. Magnetic field has been added primarily into the radiative
zone, but there is some threading of the field near the equator into the convective core. Thereensues a notable departure of behavior in our mixed-field system from that of the progenitor case
C4m. This is readily apparent in the systematic evolution of kinetic and magnetic energies (Figs.
3.2a, 3.3). The KE, ME, and MTE shown here have been averaged over the full computational
domain. We have also studied the early response of the system to either the introduction of theexternal purely poloidal field (case B, Fig. 3.2b) or of a purely toroidal field (case C, Fig. 3.2c).
Variations in ME and KE occur on the dynamical time scale (i.e. the overturning time of roughly
100 days) in each case. It is evident that the toroidal field has little effect, whereas the purely
poloidal field leads to a more rapid early growth in ME than the mixed-field case A. This suggeststhat the presence of a poloidal component, and not the twist (helicity) of the external field, is an
important ingredient leading to significant growth in overall ME. However, since it has been argued
that neither a purely toroidal nor a purely poloidal field would survive in the radiative exterior, weshall hereafter concentrate on the stable mixed-field case A.

45

10

a

Case A(Mixed Field)

5
ME
3
KE

Energy Density (107 erg/cm3)

1
10

b

Case B(Poloidal Field)

5

ME

3
KE
1
10

c

Case C(Toroidal Field)

5
3

1
0

500

1000

1500

Days

Figure 3.2: Separating the role of flux and helicity. (a) Temporal evolution of the mixed-field case
A. Volume-averaged magnetic energy ME is shown in red and KE in black. (b) Temporal evolution
for the poloidal field case B. (c) Evolution of the purely toroidal field case C. Cases A and B exhibita tendency for transition to super-equipartition. Case C exhibits no such behavior, suggesting thatthe inclusion of a poloidal magnetic field spurs the transition.

Energy Density (erg/cm3)

Energy Density (erg/cm3)

46

a

108

ME
MTE

KE
107

0

2

4

6

8

10
12
Days (Thousand)

b

108

108

107

14

16

18

20

14.5

15.0

c

107

8.0

8.5

9.0
Days (Thousand)

9.5

10.0

13.0

13.5

14.0
Days (Thousand)Figure 3.3: Temporal evolution of overall energy densities (volume-averaged) in case A after imposition of an external mixed magnetic field. (a) Full evolution sequence over 20,000 days (or about
2,900 rotation periods), revealing that the magnetic energy ME (red) has grown in strength to become about ten-fold greater than the kinetic energy KE (black). The mean toroidal field magnetic
energy MTE (blue) has also increased prominently in strength, attaining in (b) roughly equipartition levels with KE in the close-in view starting at 8,000 days, and (c) slightly super-equipartition
behavior in the interval starting at 13,000 days.

3.1.1

Temporal Evolution of Energies in Case A

The full evolution of case A is shown in Figure 3.3. Following the imposition of our fossil
magnetic field (constituting only about a 10% increase in ME), the overall ME (volume integrated)of case A undergoes a phase of gradual growth. This growth of ME continues for approximately
one magnetic diffusion time across the core radius (∼7,000 days), at which point the system reaches
a saturation level of M E/KE ∼10. Over the roughly three diffusion times (and many hundredsof dynamical times) we have evolved this case, the ratio M E/KE is highly intermittent in time,varying on time scales of a few hundred days. It reaches a peak of about 16 near day 10,000,
followed by a rapid decline to more typical values. Such behavior is in marked contrast to thatof the progenitor case C4m where approximate equipartition of magnetic and kinetic energies wasmaintained throughout the simulation.
The energy balance here is achieved through a roughly fourfold rise in ME, along with a47
twofold decrease in KE. The growth of ME is accompanied by a similar growth of axisymmetric(mean) toroidal fields (MTE) which generally account for about 10 − 20% of the total ME, muchas in the progenitor solution. A steady rise in MTE continues until approximate equipartition withrespect to KE is reached. By contrast, the energy in the axisymmetric poloidal fields (not shown)comprises about 5% of the total ME on average. Axisymmetric toroidal fields undergo a brief phaseof super-equipartition of their own near 15,000 days. Such equipartition and super-equipartition
behavior of the mean fields is examined in expanded temporal detail in Figures 3.3b,c respectively.
Here KE and ME both vary on time scales of one to two hundred days. By contrast, MTE varies
on time scales of 500 to 1,000 days and follows the long term trends in ME. Shorter trends in ME
are due largely to non-axisymmetric (fluctuating) magnetic fields within the core which vary ontime scales commensurate with the convection. The axisymmetric fields are stored largely within
the lower radiative zone where flows are much more quiescent. Diffusive processes play a more
dominant role in the evolution of the magnetic field there, resulting in slower evolutionary trends
for MTE. The variation in KE in the two extracted closeups of the time series (Figs. 3.3b,c) ismodest despite the different MTE, in large part because the bulk of kinetic energy for this system
is convective and concentrated within the core. Convective motions and axisymmetric fields thusinteract largely near the convective core and radiative zone interface, which includes the region ofovershooting.

3.2

Dynamics of Super-Equipartition State
The spatial structuring in the velocity and magnetic fields within the super-equipartition

states realized by dynamo action in case A is relatively complex in order to avoid strong Lorentz
forces that would otherwise seek to quench the flows. We now examine the nature of such structures
and their implications for the overall dynamics.

48

vr

r = 0.10 R

ar = 0.16 R

d-/+ 20Brb-/+ 1.5e-/+ 70Bc-/+ 70gh-/+ 20fi-/+ 140

Figure 3.4: Snapshot of radial velocity vr (upper panels) and the two magnetic field components
Br and Bφ (middle and lower panels ) for the mixed-field case A at about day 15,000. The tilted
orthographic views (left column) and Mollweide projection are sampled at mid-core (r=0.10R; left
and center columns), and in the region of overshooting near the core and radiative zone interface(r=0.16R; right column). Using a common color table, velocity ranges are indicated in m s−1 andmagnetic field sense and strength in kG; yellow tones are positive, dark violet negative.

49
3.2.1

Overall Properties

Strong magnetic fields generated in our mixed-field case A are present fairly uniformly
throughout the convective core. The convection continues to exhibit broad upflows and downflows spanning much of the convection zone, with smaller scale motions present in the region ofovershooting. A snapshot of the evolving flows and magnetic fields is shown in Figure 3.4, sampling radial velocities and both radial and azimuthal magnetic field components both at mid-core(r=0.10R) and in the overshooting region (r = 0.16R) at around day 15,000. Velocities and magnetic fields at mid-core are shown in both orthographic view (with the north pole tilted into view),
and also rendered in Mollweide projection, with the full spherical surface thus visible. In the latter
rendition, lines of latitude are denoted by horizontal lines (equator at middle), and those of longitude by curved lines. The radial velocity patterns (Fig. 3.4a,d) indicate roll-like motions roughlyaligned with the rotation axis and involving about four to six rolls. These rolls extend from pole
to pole, maintaining their coherence across both hemispheres. Such behavior contrasts with thatof the progenitor where convective rolls tend to waver in latitude (Fig. 2.1a). Flows in case A exhibit somewhat smaller amplitudes than those of the progenitor case. Typical convective velocities(based on rms velocity at mid-core) for case A are ∼20 m s−1 versus ∼29 m s−1 before the external
field was imposed. The radial velocity pattern in the region of overshooting (Fig. 3.4g), compared to mid-core, involves smaller scales and less alignment, consistent with intermittent upward
penetration by plumes, which involve roughly tenfold weaker velocities.
The accompanying radial magnetic fields for case A at mid-core (Figs. 3.4b,e) and in theregion of overshooting (Fig. 3.4h) exhibit similar patterns to that seen in the convective flows. Thesemagnetic structures are less spatially intermittent than their counterparts in BBT05. Azimuthal
magnetic fields display similar behavior (Figs. 3.4c, f, i) and exhibit much less banding within thecore than that observed in the progenitor case. Typical magnetic field strengths at mid-core areroughly 80 kG, or about 20% higher than those present prior to imposing the fossil field. Regions ofprominent azimuthal and radial fields exhibit approximate antisymmetry about the equator. The

50
simplified magnetic geometry within the core in case A is suggestive of a columnar, helical magnetic
field topology unlike that seen in the progenitor case C4m.
Near the edge of the convective core, some rotational shearing still exists in the transition
from weak differential rotation in the core to nearly solid body rotation in the radiative exterior.
The strong Bφ (Fig. 3.4i) fields trace out the boundaries of the largest convective rolls. Shearing
and stretching in the overshooting region allow these structures to grow into ribbons of field thatwrap their way poleward around the cylindrical convective rolls. Magnetic field strengths within
the centers of these ribbons are about 200 kG.

3.2.2

Complex Interplay of Flows and Field Structures

The intricate configurations of the strong magnetic fields realized in case A are clarified by
Figure 3.5 showing 3-D volume renderings of the magnetic energy and tracings of the magnetic
field lines, sampling the same time instant as in Figure 3.4. Magnetic energy throughout theconvective core and the region of overshooting is presented in Figure 3.5a. The regions of strong
magnetic energy density (denoted by green/orange tones) reach from equator to pole in long helicalarcs. Toroidal structures are visible here as well, with magnetic energy tending to wrap aroundconvective columns parallel to the rotational axis rather than the entire core. These impressions
are accentuated in magnetic field line tracings (Fig. 3.5b), showing complex wrappings around theconvective cells. Other volume renditions involve a narrow slice in the equatorial plane, sampling
again the magnetic energy density (Fig. 3.5c) and the kinetic energy density (Fig. 3.5d). Herewe view the equatorial plane as if from the north pole. The strong flux ropes (about 200 kG
in strength) extend throughout much of the core, tracing out the pattern of the most prominent
convective rolls found there. Regions of blue indicate strong (∼ 300 kG) regions of magnetic field.
These strong fields are part of extended helical arcs (Fig. 3.5b). The accompanying kinetic energy
rendition (Fig. 3.5d) reveals that the strongest fields are laterally displaced from the fastest flows,
recognizing that otherwise Lorentz forces may tend to suppress the flows unless the velocity andmagnetic exhibit some measure of alignment (see §4.1).

51

Figure 3.5: Volume rendering of magnetic stuctures (with rotation axis vertical) near day 15,000
accompanying Fig 3.4. (a) Magnetic energy density realized in the convection zone and region ofovershooting; yellow/green tones indicate high values, and red tones low values. (b) Magnetic fieldlines traced throughout the convection zone. Yellow/red tones correspond to high field strengths,
violet/blue tones to low. Strong flux ropes coincide with regions of high magnetic energy. (c)Magnetic energy density rendered in the equatorial plane viewed as if from the north pole, withgreatest energy levels shown in blue/green. (d) Companion kinetic energy density rendered in theequatorial plane showing that the fastest flows are positioned adjacent to the sides of the regionsof strongest field, with blue/green again highest energy density.

These dynamic structures in magnetic energy evolve on time scales of about 100 days. A senseof these intricate changes is provided by Figure 3.6 by sampling the evolution of magnetic fields inthe equatorial plane over a 350 day interval. Here close inspection reveals that magnetic structuresare seen to rise radially, stretch, and ultimately dissipate. At times large bundles of magnetic flux
cleave into each other, forming tubular structures that encircle much of the convection zone at theequator.

3.2.3

Statistical Properties of Super-Equipartition State

Topological differences in magnetic fields and flows between case C4m and case A may befurther characterized through probability distribution functions (pdf) of v and B, as shown in
Figure 3.7 for the velocity and magnetic field at mid-core. Magnetic fields (Fig. 3.7b,d) peak athigher values and exhibit a noticably greater range of variability in case A than in case C4m. This
broadening, related to the spatial intermittency, may be quantified by measuring the kurtosis K(the fourth moment) of these curves. A large value of K indicates a broad distribution, and a low
value indicates a narrow one. We find that K for Br has transitioned from a value of 4.4 in case C4m

Figure 3.6: Sequence of views showing the evolution of magnetic energy density in narrow volume renderings of the equatorial plane,
starting from day 15,000 (as in Fig. 3.5c). Green/yellow tones denote high values, red tones low values. Magnetic features stretch,strengthen, and at times join with those created by adjoining convective rolls.

52

53

Figure 3.7: Time-averaged pdf’s for case A (black) and case C4m (red) of velocities (Vr and |V |)and magnetic fields (Br and |B|) sampled on a spherical surface at mid-core (r=0.10R). Case A is
characterized by narrower velocity distributions and broader magnetic field distributions than case
C4m.

to 6.7 in case A. The convective core of case A is thus pervaded more uniformly by strong magnetic
field than in case C4m. The converse is true for the the velocity distributions (Fig. 3.7a,c), which
are narrower and peak at lower values in case A relative to case C4m. Correspondingly the K forvr has transitioned from 4.4 in C4m to 4.1 in case A. The most vigorous downflows in case A are
somewhat weaker than those in case C4m which exhibits a more prominent negative tail in vr .

3.3

Growth of Axisymmetric Helical Fields
The super-equipartition state achieved in case A involves the generation of strong axisymmet-ric (mean) magnetic structures well beyond that achieved in the absence of an external fossil field.
Figure 3.8 shows the profiles with radius and latitude of the mean toroidal and poloidal magnetic
fields for the greater part of the overall evolution. Two prominent twisted structures of opposite
helicity have developed near the equator for our mixed case. Both torii persist on time scales of54-/+ 35t = 3,100d

t = 8,550d

t = 9,350d

t = 20,050d

0/70

Figure 3.8: Temporal evolution of the mean toroidal (upper ) and poloidal (lower ) magnetic fields
for case A. Each sample involves azimuthal averages, accompanied by a 100 day temporal average(centered about the indicated times). Helical structures of opposite polarity wax and wane overthe course of the 17,000 days sampled. Magnetic field strength is shown in color; contours indicate
poloidal field lines for Bpol , and zero field strength for Bφ . Units in kG.

55order the diffusion time, but the northern, negatively signed torus tends to vary in strength morethan its southern counterpart, disappearing completely after about day 10,000. These oppositely
signed toroidal tubes with a common, enveloping poloidal component have typical field strengthsof a few tens of kG, with peak mean field strengths of around 100 kG near their centers. For
comparison, overall peak field strengths of 350 kG are achieved within the core. Evolution of thesystem for multiple diffusion times leads to the emergence of a single positively signed toroidal flux
structure in the radiative envelope.
The prominent equatorial torii that develop mainly in the radiative exterior at low latitudes
reside in a region previously characterized by strong rotational shear. The energy initially in the
differential rotation (DRKE) is insufficient to account for the bulk of the magnetic toroidal energy(MTE) contained in these structures, since the DRKE comprised only about 10% of the total KE
initially. Values of MTE near the end of the simulation (then comparable with KE) are roughly
five times the initial value of DRKE. Rather than reflecting field structures that truly encircle
the convective core, the torii of Figure 3.8 represent the mean of a more complicated magnetic
geometry as seen in Figure 3.6. Thus while an Ω-effect is certainly at work on the poloidal fields
in this region, non-axisymmetric motions are responsible for the bulk of the mean field generation
here. We discuss these generation mechanisms further in §4.1.

3.4

Nature of Differential Rotation
The coupling of convection with rotation tends to redistribute angular momentum so thatdifferential rotation is established, yet the structure and amplitude of such rotational shearing flows
is difficult to predict. Browning, Brun & Toomre (2004) explored the differential rotation established by core convection in the purely hydrodynamic version of our A-type star. The convection
yields a slow column of differential rotation throughout the core. With the inclusion of magnetism,
and the realization of dynamo activity as in case C4m, this differential rotation is considerably
weakened, but the slowly rotating (retrograde) central column of fluid (Fig. 2.1d) is retained. The
radiative zone close to the equator exhibits a prograde mean zonal flow hvφ i.

56

v

(nHz)

a

b
45 N
A

0

C4m

-8

8

Figure 3.9: Differential rotation realized in Case A. (a) Zonal velocity vφ averaged over longitude
and over 200 days around day 15,000. (b) Radial cuts of mean angular velocity Ω at low (0◦ ) andmid latitudes (45◦ N) for case A (black) and case C4m (red). The radiative zone is in a state ofnear solid body rotation for case A. The system frame rate Ω0 is 1,655 nhz.

Imposing mixed or poloidal external magnetic fields alters the differential rotation until theradiative zone is nearly in solid body rotation (Fig. 3.9a,b). The slow central column of retrograde
flow in the progenitor case is largely absent, replaced by a region of slightly prograde rotation close
to the rotation axis. Magnetic torquing resulting from correlations between the imposed poloidal
field, and the subsequently developing toroidal fields, transports much of the angular momentum
from the radiative zone (seen as the faster vφ in Fig. 2.1d) into the convective core. Radial andlatitudinal shear throughout much of the domain is diminished at late times, particularly near thecore and radiative zone interface. A weak retrograde jet has formed coincident with the location of
the large helical magnetic structures (evident in Bφ in Figure 3.8) formed close to the interface. The
super-equipartition state arising from the imposition of a mixed magnetic field has thus significantly
diminished the differential rotation present in the progenitor case.
Weak departures from solid body rotation in the outer radiative envelope of case A are
generally columnar in nature, with lines of constant Ω cutting across field lines in the radiativeenvelope. There the angular momentum transport is dominated by magnetic torques and meridional

57
circulations. Magnetic torques realized in case A are much stronger relative to case C4m owing
to the large axisymmetric fields realized in the radiative envelope of case A. However, transportof angular momentum by meridional circulations prevents a pure state of isorotation as suggested
by Ferraro’s law (Ferraro 1937) along magnetic field lines, ultimately maintaining a differential
rotation similar in nature (but weaker in amplitude) to that achieved primarily through meridional
circulations and viscous transport in the progenitor case C4m.
Before proceeding, we note that our perfect conductor boundary conditions, which force
magnetic field lines to close within the domain may also play some role in sustaining the nearly solidbody rotation observed in the radiative zone. The progenitor case C4m employed a radial magnetic
field boundary condition which may have favored the columnar differential rotation realized in thatcase. The closing of our magnetic field lines within the domain will necessarily disrupt a tendency
toward such columnar differential rotation in the radiative zone. In the super-equipartion regime,
the core dynamo has a tendency to build strong horizontal magnetic fields at the lower boundary,
and a perfect conductor lower boundary condition is more appropriate for case A. The appropriate
boundary condition for the outer boundary is less clear. Were the strong toroidal fields generated at
the core-radiative zone interface to become magnetically buoyant, adopting a radial magnetic fieldboundary condition at the top of the domain may become problematic as the (likely toroidal) field
is forced to change its geometry over the last few grid points. However, we have not yet observed
the buoyant rise of magnetic field in case A, and we note that it may thus be interesting to examine
what effects of a radial magnetic field boundary condition at the outer boundary may have on the
differential rotation achieved in the radiative zone.

3.5

Distribution of Magnetic Energy
Magnetic energy in case A is largely contained within the convective core of the A-type star.

However, convective overshoot, meridional circulations, and diffusion serve to transport some magnetic field into the overlying radiative zone. We now examine several components of the magnetic
energy and their distribution throughout the computational domain. Various time-averaged mag-

58
109
Energy Density (erg/cm3)

FME
KE

108

107

MTE

MPE

106

105
0.05

0.10

0.15
r/R

0.20

0.25

Figure 3.10: Variation with proportional radius of horizontally-averaged kinetic and magnetic energy densities for case A. Profiles are averaged over 100 days near day 15,000. The penetrative
convective core boundary is roughly at r=0.16R.

netic and kinetic energy densities are shown as a function of radius in Figure 3.10. These quantities
have been averaged over the full sphere at each radius and then temporally averaged over the interval spanning 15,000 to 15,100 days. Magnetic fields are super-equipartition compared to the kinetic
energy throughout the core and the region of overshooting owing largely to fluctuating magnetic
fields (FME). The mean field energies MTE and MPE are largely sub-equipartition throughout thecore. However, the mean toroidal fields peak to super-equipartition values near the interface, owing
to the development of large toroidal bands of magnetic field there. Magnetic energies of all types
fall off substantially outside of r∼0.16R, the approximate outer edge of the region of overshooting.
Magnetic energies in the remainder of the radiative zone are achieved largely through diffusive
processes. The energy balance both inside and outside the core is detailed in Table 1. All energies
are quoted relative to the average KE in each region.

3.6

Growth of Large-Scale Magnetic Structure
The addition of our fossil field has led to the development of larger-scale magnetic structures

than those present in the progenitor simulation. Figure 3.11a depicts the time-averaged magnetic

59
TABLE 1
Magnetic to Kinetic Energy Ratios
RegionConvective Core Radiative Envelope
M E/KE
7.38
5.18
F M E/KE
6.36
3.75
M T E/KE
0.64
1.22
M P E/KE
0.38
0.21

Table 3.1: Time and volume averaged magnetic energies compared to kinetic energies. Temporal
averaging covers days 15,000 to 15,100 as in Fig. 3.10. Volume averaging is carried out over theconvective core (r ≤ 0.15R) and over the radiative envelope (r > 0.15R).

energy spectrum (with spherical harmonic degree ℓ) soon after the imposition of an external field incase A. Shown also is the spectrum for a more mature stage (∼ day 15,000) during the simulation
after super-equipartition magnetic fields have developed. Individual spectra were created using
spherical shells at mid-core and averaged over about 200 days in each case. Magnetic energy has
increased on scales as small as spherical harmonic degree ℓ = 60, but the larger scales (ℓ ≤ 10)
exhibit the most growth. The mixed external field has thus led to the growth of the large-scale
components of the magnetic field, while leaving the energy in the smaller structures relatively
unchanged. Moreover, these structures are largely non-axisymmetric, with axisymmetric structurescomprising only about 20% of the magnetic energy on average.
The accompanying kinetic energy spectra are shown in Figure 3.11b . As the dynamo reaches
super-equipartition levels, kinetic energy decreases by a factor of about two across much of thespectrum. This trend reflects the overall decrease of kinetic energy shown in Figure 3 and thetransition of the system to a less turbulent state.
The growth of large-scale magnetic structure is also apparent in volume renderings of themagnetic energy in Figure 3.12, showing snapshots of the magnetic energy density just prior toinserting the external field and during the super-equipartition phase (day 15,000). The renderingof the full convective core conveys a good sense of the helical convection in the core. The ribbons
of magnetic energy seen at early times (Fig. 3.12c) have become noticebly wider in the superequipartition regime (Fig. 3.12d). Similar morphology is apparent in renderings of ME in the

60aLate Times

1010

Energy Density (erg cm-3)

108

Early Times

106
Magnetic Energy
1010b
Early Times

108

Late Times

106

Kinetic Energy
104
1

10
Spherical Harmonic Degree

100

Figure 3.11: Magnetic and kinetic energy spectra sampled at mid-core in case A for early and latetimes. (a) Time-averaged initial ME spectrum (dashed line) and the spectrum realized after 15,000
days (solid line). (b) Accompanying KE spectra. The imposition of mixed fields in case A has
yielded a dynamo with more prominent global-scale magnetic fields than in progenitor case C4m.

61

Figure 3.12: Volume rendering of magnetic energy density viewing the entire convection zone (a)prior to the insertion of the fossil field, and (b) 15,000 days later. The rotation axis is vertical.
Structures are noticably larger in scale in the later super-equipartition regime. Companion imagesof magnetic energy density viewing the equatorial plane from the north, showing (c) progenitor
state and (d) mature super-equipartition state.

equatorial plane. The early time (Fig. 3.12c) exhibits much smaller-scale fields than the later time
(Fig. 3.12d). The shift in peaks of KE and ME in the energy spectra of Figure 3.11 suggests thatthe dynamo now running in case A is much more global-scale in nature than that in the progenitor
case C4m. We turn now to a detailed examination of the mechanisms by which case A sustains
such a large-scale, super-equipartition state.

Chapter 4

Maintaining the Strong Dynamo

4.1

Sustaining a Super-Equipartition State
The super-equipartition state achieved here in which magnetic energy ME is tenfold greaterthan the convective kinetic energy KE is quite remarkable, for one might anticipate that the feedback from the Lorentz forces would prevent this. In contrast, progenitor case C4m was nearly inequipartition by this measure. So how does case A achieve and sustain such strong magnetic fields?
The shift towards larger-scale convective motions and magnetic field structures suggests that amore global field/flow topology is necessary for this state. The high magnetic field strengths realized in case A might be anticipated to yield Lorentz forces that react back on the flow, ultimately
quenching the flows that generate such fields through an α-quenching process (Cattaneo & Hughes
1996). For our energy balance to survive, the flows must organize themselves in a fashion thatminimizes these Lorentz forces while still maintaining sufficient induction to balance the dissipative
processes present in the system (which also serve here to provide some field reconnection). Thismay be accomplished by separating regions of strong flow from those of strong magnetic field, orby some tendency to align j and B and thereby lessen Lorentz force feedbacks. We turn now to a
discussion of these issues and of the dominant mechanisms responsible for magnetic field generationin case A.

63

Figure 4.1: Dual rendering of kinetic energy (red) and magnetic energy (blue) near day 15,000.
Regions where both quantities are strong (relative to their respective rms values) appear in green.
(a) Mollweide view at the interface between the core and the radiative zone (r=0.15R), and (b)at mid-core, (r=0.075R). (c) Partial view of the equatorial plane from the pole. The dotted line
denotes the edge of the convective core, and the outer dashed line r=0.19R, capturing a portion of
the radiative zone (including the region of overshooting).

64

Figure 4.2: Comparing partitioning of fast flows and strong fields in a snapshot sampled in theequatorial region. (a) Instantaneous flow tangent lines (day 9,000). Coloring in blue/green indicatesthe fastest speeds (∼30 m s−1 ) and regions of violet/pink that of weaker flows (∼5 m s−1 ). Six
prominent rolls are present. (b) Patterns of stronger magnetic energy are similar to those of thefaster flows though somewhat displaced. Regions of bright green denote ∼250 kG fields, reddish
regions ∼30 kG fields. (c) Flow lines and magnetic energy density rendered together. Regions offast flow and strong fields exhibit notable alignment.

4.1.1

Surviving Lorentz Feedbacks

We can assess the nature of alignment and displacement of field and flow by examining theirlarge-scale organization as illustrated in Figure 4.1. Here instantaneous snapshots of magnetic
energy and kinetic energy near day 15,000 are overlaid at two depths and in an equatorial cut.
Regions of strong ME are shown in blue, strong KE in red, and regions where both are strong(relative to their rms values) are shown in green. Clearly the patterns in KE are similar to thoseseen in ME. However, some of the features are laterally displaced, with arcs of ME running adjacent
to arcs of KE. Some patches of overlap are visible as well, particularly in Figure 4.1c which suggestthat some strong magnetic structures merge with fast velocity structures. For such overlapping
regions to survive, Lorentz forces must not be strong enough locally to disrupt the flow. This, inturn, implies some alignment of j and B. The current density j depends, through Ohm’s law, on
both the electric field and on the alignment of flows and fields (through v × B). Minimizing theLorentz force locally thus typically implies some alignment of v and B.
The large-scale alignment of flows and fields may be explored through volume renderings

65

Figure 4.3: Probability distribution functions (pdfs) of the relative angle θ between v and B
measured throughout the convective core for cases A (black) and C4m (red), averaged over about
50 realizations spanning 100 days. Flows and magnetic fields exhibit more alignment in the superequipartition regime of case A than in case C4m. The alignment of flow and magnetic fields incase A is more pronounced in regions of strong flow and field, shown by plotting pdfs only for sites
where local values of both v and B exceed their rms values (labeled A1) or twice those values (A2).

as in Figure 4.2a where instantaneous streamlines are depicted in the equatorial region. The
corresponding magnetic energy density (Fig. 4.2b) and a combined rendering of these two quantities(Fig. 4.2c) are illustrated as well. The large, cloverleaf-type patterns in the magnetic energy
density are typically tangent to the local magnetic field vector. Rendered in the same image, thetwo patterns trace each other closely, with alignment most evident in regions of strong magnetic
field. Another measure of the alignment of flows and magnetic fields is achieved by examining
the relative angle θ between v and B throughout the core. Figure 4.3, examining probability
distribution functions (pdfs) of that angle, shows that our evolved case A typically exhibits morealignment than the progenitor case C4m. By weighting these by the local KE and ME, thus
concentrating on regions of fast flow and strong field, the alignment becomes more pronounced.
This may be seen by forming pdfs where both v and B exceed a threshold of either just their rms
values (curve A1) or twice those values (A2). For the latter, the pdf peaks at θ ∼10◦ , compared tosampling all sites which yields a broad maximum centered on 25◦ .

Figure 4.4: Columnar convection for case A visualized using instantaneous streamlines near day 15,000. Streamlines are colored by velocity
component along the rotation axis vz . Blue (yellow) tones indicate northward (southward) motion; equatorial plane indicated by light
blue. (a) View of convective core from the north. The core is typically dominated by four to six prominent rolls. (b) Individual columnar
cell from the roll complex. The columnar motion extends across the equatorial plane. Tilting of the orbits due to the spherical boundary of
the core is visible. (c) Interior view of the same roll. Many rolls in this system possess an axial flow. Such flows freely cross the equatorial
plane as seen here. (d) Rendering of that roll along with a neighboring roll. Axial flows tend to link neighboring rolls near the edge of theconvective core.

66

67
4.1.2

Global Connectivity of Structures

Figure 4.2 showing an equatorial cut has emphasized that the core is typically occupied byfour to six cylindrical convective rolls, with their axes primarily aligned with the rotation axis,
as suggested by Figure 3.4d. The geometry and connectivity of the flows in this system play asignificant role in how the strong magnetic fields are built and sustained. Turning to analysis offlow streamlines shown variously in Figure 4.4, we see that the rolls typically possess a narrow core,with a prominent upward or downward velocity (say vz ) along the rotation axis. Within each roll,surrounding the fast core is a broader region of circular flow in a plane nearly perpendicular to
the rotation axis. As these convective motions overshoot into the radiative zone, buoyancy forces
brake the radial component of the flow. The resulting tilt of the roll motions reflects the slope of
the intersection of the roll with the edge of the convection zone and is similar to that observedin simulations of rotating spherical shells with rigid boundaries (e.g., Busse 2002). These outerregions of the rolls exhibit a weaker vz than in the fast core of the roll, though of the same sense,
thus yielding helical flow lines.
It is quite striking that the fast axial core flow in one roll can often be seen to connect near theradiative interface to the core flow (of the opposite sense) in an adjacent roll. The axial core flows
extend along with the overall roll structure across the equator, though there are some examples inwhich the sense of the axial flow changes at the equator. These axial core flows can thus serve toeffectively connect the two hemispheres, in addition to linking neighboring rolls.
The coupling of the rolls is also evident in the time sampling of both streamlines and magnetic
energy density shown in equatorial cuts in Figure 4.5. Close study shows that there are frequent
events in which streaming flows extend across the full domain, thereby coupling opposite sides of theconvective core. It appears that the genesis of such streaming flows comes from two counter-rotatingrolls on one side of the convective core (yielding a zone of converging flow toward the rotating axis)
becoming correctly phased with another set of rolls on the opposite side that are driving someflow away from the axis. Such large-scale flows appear to be crucial to establishing the largest

Figure 4.5: Sampling the evolving flow streamlines (a-d) and accompanying magnetic energy density (e-h) close to the equatorial planein four time instants each separated by about 50 days (starting at day 9,000). Violet tones indicate positive motions in the y-direction,and yellow tones negative motions. Regions of strong ME are shown in yellow/green tones.

68

69
magnetic structures built by the dynamo action. As these streaming flows couple distant portionsof the core, they stretch and advect the fields into the configurations seen in the the evolving ME
patterns. Figure 4.5 shows that broad swathes of magnetic field extending across the convective
core are nearly coincident with the fast core crossing flows. Thus it appears that attaining globalscale magnetic field within a core that possesses little differential rotation, yet senses Coriolis forces
significantly, relies on the large-scale connectivity among the evolving system of rolls.
Although the convective rolls are fairly complex and variable in time, the more prominentrolls maintain their coherence for a few tens of rotation periods. The prominent axis-crossing
motions serve to connect rolls on opposite sides of the core, stretching magnetic field from one rollto another as they do so. The resulting magnetic field structures thus thread through multiple rolls,
connecting across the rotation axis and along the periphery of neighboring cells. The topology has
thus become much more global in scale, resulting in the notable spectral change seen in Figure
3.11.

4.1.3

Details of Magnetic Field Generation

As the magnetic energy in case A is manifest predominantly in the fluctuating (non-axisymmetric)
magnetic fields, we focus on the generation of FME here. The generation of energy in the fluctuating fields may be computed by taking the dot product of the induction equation with the fluctuating
magnetic field B ′. Contributions to the generation of FME arise from correlations between fluctuating flows and fields with both their fluctuating and mean counterparts. However, while mean
flows and magnetic fields contribute to the generation and destruction of FME, we find that the
primary balance is struck between terms involving shearing, diffusion and advection associated withthe fluctuating components of the flow and magnetic field so that∂(F M E)1 ′≈0≈B · [(B ′ · ∇)v ′| {z }∂t
4π
shearing−∇ × (η∇ × B ′ ) −(v ′ · ∇)B ′]
|
{z
} |
{z
}
diffusion
advection

(4.1)

70
holds approximately. We examine the contributions of these individual terms to the generation anddissipation of the constituents of FME (namely FPE and FTE) in Figure 4.6.
The local enhancement and diminishment of FPE by advection (Fig. 4.6a) reflects thecolumnar nature of the rolls. Orbital motions about the axes of these rolls efficiently carry FPE
from regions of generation (via shearing as in Fig. 4.6b) to the rest of the core. Generation of FPE
occurs primarily near the rotation axis and near the outer boundary of the core. The noticably
weak columnar region of generation at mid-core is roughly coincident with the location of the rollaxes. Orbital motions of the rolls are largely parallel to the equatorial plane throughout most of
the core. As buoyancy brakes the radial component of these motions in the region of overshooting,the resulting motions tend to conform to the spherical boundary of the core. The ensuing tiltingof the orbits (see Fig. 4.4) sets up a region of shear between the rotation axis and the outer corethat serves to generate field parallel to the core boundary. Additional generation of field in theouter core also occurs through shearing in the azimuthal direction that serves to convert Bφ into
Br and Bθ there. Generation of FPE near the rotation axis, somewhat stronger than that seen in
the outer core, is accomplished primarily through generation of Br due to shearing of Bφ . The close
proximity of inflows and outflows of the columnar rolls near the rotation axis sets up regions ofstrong shearing of radial flows (see Fig. 4.5a,c) that efficiently stretches the strong azimuthal field
generated in this region into radial field. Core-crossing flows in particular, such as that prominently
visible in Figure 4.5d, have a similar effect in this region.
The effect of advection on FTE (Fig. 4.6d) is similar to that seen with FPE in that circulations generally carry energy away from the region of strong generation near the rotation axis (Fig.
4.6e) to the outer core where less generation is occurring. Generation of FTE occurs primarily near
the rotation axis as well. The convergence of multiple rolls at the rotation axis sets up gradientsof vφ in the direction perpendicular to the rotation axis, resulting in strong shear in this region.
Moreover, core-crossing motions such as in Figure 4.5d will effectively shear radial field orthogonal
to these motions, generating azimuthal field in the process. An additional effect, not included in
the decomposition of Figure 4.6, is the action of mean shear near the rotation axis (see Fig. 3.9a)

71

FPE

a

b

c

-/+ 30-/+ 60-/+ 60

Advection

Shearing

Diffusiondef-/+ 10-/+ 20-/+ 20FTE

Figure 4.6: Generation and dissipation of fluctuating magnetic energy for case A. (top row ) Generation of FPE from (a) advection, (b) shearing and (c) diffusion. (bottom row ) Generation of
FTE likewise from (d) advection, (e) shearing and (f ) diffusion. Generation terms are shown inerg cm−3 s−1 and have been averaged in azimuth and in time over 200 days near day 15,000.

72
acting on the non-axisymmetric Br and Bφ . Very near the rotation axis, the generation of FTE
through this shearing effect is comparable to that due to the non-axisymmetric flows. We notethat no similarly strong effects of mean flows on generation of FPE is observed. An additional
region of strong generation of FTE is present at the outer boundary of the convective core. In
contrast to FPE, the geometry of this region is largely spherical in nature owing to the transition
from relatively strong azimuthal flows in the core to a region of nearly solid body rotation in theradiative zone.
The effects of diffusion on FPE and FTE (Figs. 4.6c,f ) is strongest near the rotation axis
and weakest in the outer core where advection mitigates the effects of generation more readily. The
net effects of diffusion in the regions of generation is thus to destroy FTE and FPE. Beyond theregion of overshooting, however, diffusive processes serve to carry energy from the outer edge ofthe convective core into the radiative zone where generation of magnetic field is minimal.

4.2

Multipolar External Fields
We have thus far investigated the effects of an external field involving only a dipolar (ℓ = 1)

magnetic field. However, measurements of magnetic fields in Ap stars suggest that the surface
fields of these stars involve multipolar components with spherical harmonic degrees greater thanℓ = 1. It is appropriate to examine the effects of more complex magnetic topologies on the system.
Does the scale of the imposed magnetic flux matter, and if so, how? We have thus examined two
additional external fields configurations, with ℓ = 2, m = 0 and ℓ = 4, m = 0. In each case, thestrength of the fossil field was adjusted so that the integrated unsigned magnetic flux across theconvective core boundary was identical to that in the mixed-field case A. Unlike the mixed-field
case, neither of these two cases possessed a net flux of magnetic field through either hemisphere.
Moreover, here there was no magnetic flux linking the northern and southern hemisphere of thestar across the equator.
We show the evolution of the energies for these two cases in Figure 4.7. In each a transition
to super-equipartition is evident. However, the growth of the ℓ = 2 case (Fig. 4.7a) is faster than

73

10

a

=2

8
ME

Energy Density (107 erg/cm3)

6

4
KE

2
6

b

=4

5
4
3
2
1
0

1000

2000

3000

4000

5000

Days

Energy Density (erg/cm3)

109

c

=2
=4

108

progenitor

107

106
1

10
Spherical Harmonic Degree

100

Figure 4.7: Multipolar cases. (a) Temporal evolution of the quadrupolar ℓ = 2 case. Volumeaveraged ME is shown in red and KE in black. (b) Evolution for the ℓ=4 multipolar case. (c)Spectra at mid-core for both cases at day 4000. ℓ=2 is shown in black, and ℓ=4 in red. The
progenitor case is shown in blue (dotted) for reference.

74
that of the ℓ = 4 case (Fig. 4.7b). Each of these exhibit growth that is in turn slower than the
dipolar cases. The mechanisms that cause the growth clearly have a dependence on length-scale.
Magnetic energy spectra for these cases (Fig. 4.7c) are qualitatively similar to those of the mature
mixed-field case A. In both instances, larger-scale fields have grown more than the smaller-scale
fields.
Evolution of the axisymmetric fields for the ℓ = 2 and ℓ = 4 cases are shown in Figures 4.8
and 4.9 respectively. The effects of weak differential rotation (Ω-effect) initially operating in theradiative zone are evident in snapshots of the toroidal fields. Such banded toroidal structures give
way to those produced by processes in the core as the initial poloidal field component diffuses away
and the differential rotation subsides. In each case, the imposed magnetic flux is reprocessed by
the core until a roughly dipolar field configuration is realized within the core. Continued pumpingof this field into the radiative zone, combined with diffusion of the imposed fossil field placed there,
leads ultimately to a magnetic field geometry which is largely dipolar in the radiative zone in bothcases. This behavior suggests that linkage of the initial fossil field to that inside the core is not
crucial to achieving a super-equipartition state. Rather, it seems that the scale over which magnetic
field lines are able to connect distant regions of the core may be key.

4.3

Conclusions and Perspectives
Our work has revealed that a remnant primordial field linking the convective core and theradiative envelope of an A-type star may significantly impact the dynamics of convection within
the core. By placing a modest fossil field into the radiative zone involving a 10% increase in the
total magnetic energy, the dynamo responds by a five-fold increase in the overall magnetic energy
in the system. With an accompanying two-fold decrease in convective kinetic energy, the dynamo is
now in a super-equipartition state with the ratio of ME/KE of order ten. The bulk of the magnetic
energy (about 85%) achieved through the resulting dynamo action is contained in the fluctuating(non-axisymmetric) magnetic fields, with about 10% in mean (axisymmetric) toroidal fields andabout 5% of the energy contained in the mean poloidal fields.

75

B-/+ 20 kG

0/30 kG|Bpol|
Figure 4.8: Mean field evolution for the quadrupolar ℓ = 2 case. (upper row ) Azimuthally-averaged
Bφ at three instants in time spanning 4,000 days. (lower row ) Corresponding mean poloidal fields.
White lines denote poloidal field lines. Solid lines are counterclockwise, and dashed lines areclockwise. The colored underlay indicates the strength of the mean poloidal field modulus.

76

B-/+ 20 kG

0/30 kG|Bpol|
Figure 4.9: As in Fig. 4.8 showing the ℓ = 4 multipolar case. (upper row ) Azimuthally-averaged
Bφ at three instants in time spanning 4,000 days. (lower row ) Corresponding mean poloidal fields.
Magnetic fields in the ℓ = 4 case also evolve toward a dipole in the radiative zone.

77
The mean magnetic fields built in this regime are substantial, with peak mean toroidal fieldstrengths of around 105 G (about one-third of the peak fluctuating fields found in the core). The
mean toroidal fields typically manifest as large, oppositely signed torii encircling the convective core
just above and below the equatorial plane. These bands wax and wane in strength over time, butmaintain their approximate anti-symmetry about the equatorial plane. The mean poloidal fields
in case A exhibit a prominent dipolar component with strengths of around 2×104 G near the edgeof the core. While the strength of the poloidal field is variable in time, we have observed no flip ofthe dipole moment in the roughly three magnetic diffusion times spanned by the simulation.
We have examined the sensitivity of the super-equipartition transition to the structure of theimposed mean field as well. Toroidal external fields with some threading through the core haveno effect. Rather, the transition requires an external field with a poloidal component. Systems
with smaller-scale poloidal fields (ℓ = 2 and ℓ = 4) take longer for a transition to occur but reachend states similar to case A and case B (ℓ = 1). All of these systems are characterized by superequipartition between the magnetic and kinetic energy, diminished differential rotation throughout
the core, and a tendency to develop a mean poloidal field configuration that is dipolar in nature.
The sign of the dipole moment for the ℓ = 2 case is the same as that of case A, but that for theℓ = 4 case is oppositely directed, suggesting that the external field may bias the system in somefashion.

4.3.1

Sustaining the Strong Field Dynamo

The maintenance of such strong magnetic fields is an intricate matter, particularly in lightof the strong Lorentz forces that potentially accompany them. Organization of flows in a fashionthat minimizes Lorentz feedbacks is thus crucial to this balance. The velocity configuration in case
A is characterized by four to six prominent convective rolls aligned with and extending along therotation axis. Convective motions in these rolls are largely perpendicular to the rotation axis except
in the roll center. There pressure gradients induced by tilting of the roll orbits at the sphericalboundary of the convective core induce flows along the roll axis. Some of these axial flows originate

78
and terminate at the equatorial plane, while others flow freely between the two hemispheres (Fig.
4.4). The substantial reduction in differential rotation in this system has allowed for axis crossing
streaming flows. When rolls phase properly, inflow (i.e. toward the rotation axis) from one roll
may cross the rotation axis, merging with the outflow of another roll.
The convective state achieved in case A thus exhibits motions that are much more global innature (Fig. 3.11) than case C4m, and the core is in general more topologically connected as a result,particularly through the intermittent core-crossing motions. The magnetic fields of case A follow
a similar trend, exhibiting a much more global-scale field topology than in case C4m. Magnetic
structures typically wrap around the boundaries of the convective columns in large arching helicalbands of field. Moreover, core crossing flows display a clear magnetic signature (Fig. 4.5), withlarge bundles of magnetic flux stretching across the core, linking rolls on opposite sides of the core.
We find that the coexistence of our convection with these strong magnetic fields relies on agreater coalignment of the flows and magnetic fields than that achieved in case C4m. Magnetic
field structures tend to mimick patterns visible in streamline renditions of the flow-field, suggestiveof a guiding effect of the magnetic field on the flow. However, regions of strong field are typicallyseparate from (albeit often adjacent to) regions of strong flow. This is particularly true of the fast
axial flows in the center of the convective columns. We find these regions to be largely devoid ofmagnetic energy.
Our primary case A bears some similarities to phenomena observed in geodynamo modeling.
In the absence of magnetic fields, convection in those models also exhibits a strong tendency for
Taylor columns aligned with the rotation axis (e.g. Roberts & Glatzmaier 2000). As with ourcase A, the convection within their cells is helical, with axial flows setting up circulations linking
neighboring columns at high latitudes and near the equatorial plane. In the presence of magnetism,
this convection maintains its columnar nature, and the associated helical motion serves to create anefficient α-effect as magnetic field is stretched and twisted around the convective rolls (e.g. Olsenet al. 1999). Super-equipartition is found in these models as well, with typical values of ME/KE
of 10-20, but with some reaching as high as 103 (Kuang & Bloxham 1999, Glatzmaier et al. 1999).

79
This regime, referred to as the strong field regime, involves comparable Lorentz and Coriolis forces.
As the inner 35% of the Earth’s core is solid, geodynamo models have an inner sphere that flows
may not cross. Dynamics closer to the rotation axis than the tangent cylinder associated with thisinner sphere can be considerably different than that occurring within the outer convective rolls.
Olson et al. (1999) for example note the presence of thermally driven upwellings along the rotation
axis in the tangent cylinder that tend to modify the mean poloidal field, effectively lowering thepoloidal flux through the core boundary at the poles. We observe no such effects in our A-star
simulations as the tangent cylinder is much less pronounced, with flows crossing the rotation axis
freely above and below our small inner cutout.
The force balance achieved in geodynamo models typically involves some balance between
the Lorentz, Coriolis, and buoyancy forces (Christensen & Aubert 2006). Such a force balance
was not present in the progenitor C4m which was largely geostrophic in nature. However, thesuper-equipartition state realized in case A does exhibit such a triple force balance in the direction
parallel to the rotation axis. Pressure plays a stronger role in the horizontal, where the predominantbalance is between pressure, Coriolis, and Lorentz forces.
Inertial forces are weak in geodynamo models. Due to the rapid rotation of the Earth andthe relatively weak convective velocities achieved, geodynamo models are much more rotationallyconstrained than our A-star models. A typical Rossby number, which measures the inertial force
relative to the Coriolis force, is about 10−5 in the geodynamo versus about 10−3 for our A-star.
We find the inertial term in case A, while diminished relative to the progenitor case C4m, to be
non-negligible in our simulations, typically contributing to the force balance at the 10% level.
Saturation in geodynamo models can involve the advection (driven by the axial circulations
within the rolls) of poloidal field into the centers of anticyclonic rolls. Equilibration between
Lorentz and Coriolis forces is achieved in these rolls through an outward directed Lorentz force thatcounteracts the inward Coriolis force (Sakuraba & Kono 1999, Olson et al. 1999), a phenomenon
that we do not observe in our super-equipartition cases. Instead for our case A, the magnetic fields
exhibit a tendency to encircle the rolls, leaving the roll centers relatively devoid of magnetic field.

80
The encircling field has the net effect of squeezing individual rolls. For a given intersection of tworolls pervaded by magnetic flux, Lorentz forces may work in tandem with the Coriolis force on
one side of the intersection and against it on the other. The similarities of our super-equipartition
dynamos to the geodynamo are nevertheless striking. It seems that by imposing a modest external
field we have enabled our system to transition toward the strong-field dynamo regime.

4.3.2

Possible Field Emergence

We have demonstrated that the presence of a fossil field may induce the core dynamo of an Astar to transition to a more laminar but stronger dynamo state. From an observational standpoint,
however, it is prudent to ask what implications such a state may have for magnetic fields at thesurface. While we have not observed the rise of buoyant magnetic structures in this simulation,
we do find that more efficient generation of mean fields at the edge of the convection zone has
led to the development of 105 G fields there. Even stronger fields may be required for structures
to become buoyant and eventually reach the surface, as has been considered when studying suchpossibilities for the more massive O and B stars (MacGregor & Cassinelli 2003; MacDonald &Mullan 2004). They suggest that fields with super-equipartition strengths would be needed toachieve field emergence. It is thus of some interest that we have found similar strong-field dynamo
states for our less massive A-type stars.
We have yet to ascertain if our core dynamos have attained the upper limit for magnetic fieldstrengths. Indeed, periods of brief super-equipartition are apparent even with the mean fields in
Figure 3.3. What sets the strength and temporal extent of these phases is unclear, and may rely
on fortuitous phasing of the convective rolls and magnetic field near the core boundary. Moreover,the development of sufficiently strong flux tubes may rely on adopting much lower diffusivities in
the radiative zone, such as through hyperdiffusivity. Our eddy diffusivity parameterization for the
SGS motions in this region may simply be too diffusive to allow the buildup of 106 G or greater
fields.
Such considerations encourage detailed studies of what is required for magnetic structures

81
to become unstable and reach the surface in these A-stars. Our simulations have revealed thata distinctly super-equipartition state can be achieved by the core dynamo, yet judging from thecomplexities realized with geodynamos, there may well exist a wide range of strong field states thathave yet to be explored.

Chapter 5

Helioseismic Ring Analyses with 3-D Inversions

5.1

Introduction
The local helioseismic technique of ring analysis (Hill 1988) has yielded many valuable insights

into the nature of flows in the upper convection zone and their interaction with magnetism there.
Ring analyses assess subsurface velocities below a localized region of the Sun by measuring the
Doppler shifts in the resonant acoustic modes induced by the presence of a flow. Different modes
average the subsurface flow differently in depth, and frequency splittings measured using different
modes can be inverted in depth to determine the depth dependence of the flow. By tiling the solar
disk with a number of overlapping analysis regions, a solar flow map in longitude, latitude, anddepth may be constructed. In their current formulation, ring analyses employ only one tiling sizein the generation of such flow maps. We have developed a novel method for generating these flow
maps by self-consistently incorporating measurements made using tiles with multiple horizontal
resolutions.
Mapping the surface of the Sun using ring analysis measurements has already begun to reveal
the intricate nature of flows in the upper convection zone. One of the most notable examples is thepresence of evolving large-scale meandering flows, now known as solar sub-surface weather (SSW),
that exhibit a clear interaction with magnetic features through their inflows and outflows (e.g.,
Haber et al. 2002). The vortical nature of these flows in the presence of active regions has beena subject of particular interest due to their likely role in the twisting and shearing of magnetic
field, thought to lead to eruptive events such as flaring and coronal mass ejections (Haber et al.

83
2004; Komm et al. 2005, 2009, Hindman et al. 2006, 2009). When applied to the measurementof meridional circulations, ring analysis has shown that these global-scale motions also possess aclear linkage to magnetic activity, varying in strength and location over the course of the solar cycle(Haber et al. 2002; González Hernández et al. 2008; Basu et al. 2010).
The rich potential for ring analysis to illuminate the variation of flows in the upper convection
zone has been only partially exploited, however, owing to the single resolution and 1-dimensionalityof the applications thus far. In a ring analysis study, a single analysis region size, or tile size, ischosen, and the disk is tiled with a number of these analysis regions, overlapping one another byroughly half their horizontal extent. The resulting flow maps are constructed with a horizontal
resolution determined by that tiling size, and at each position on the disk, a given measurement
represents a fairly uniform horizontal average over this tile. The choice of tiling size is determined by
the scientific objective. Small analysis regions permit the measurements of the shallowest modes andafford finer horizontal resolution. Measurements from larger tile sizes employ longer-wavelength,
more deeply penetrating modes, and sample depths of the convection zone inaccessible to thesmaller tilesizes, albeit with lower horizontal resolution. Thus, for the study of deep meridional
circulations (e.g., González, Hernández et al.(2000) a tile size of 30◦ might be selected, while forstudies of the detailed interactions between flows and magnetic fields in the upper convection zone,
a tile size of 2◦ might be used (e.g., Hindman et al. 2006). Thus when high resolution is chosen, it
is done so at the expense of depth sensitivity and vice versa.
Many applications of ring-analysis stand to benefit if deep mapping of the subsurface flow canbe accomplished alongside high resolution mapping of the near-surface layers. Studies that examine
the small-scale nature of convection along with its large-scale means stand to benefit most. As wediscuss in detail in §5.2, 2◦ measurements have larger errors than the lower resolution 16◦ measurements, primarily owing to the smaller amount of spatial averaging involved in computing the2◦ power spectra. Power spectra are a nonlinear data product, and it is unlikely that the noisy
2◦ measurements will yield the same large-scale mean as the 16◦ measurements when spatiallyaveraged to the larger scales. Moreover, the incorporation of measurements from multiple tile

84
sizes can lead to a natural progression of horizontal resolution with depth. One might then assess
how the surface inflows observed near active regions transition into outflows at depth and how thesmall-scale supergranular patterns around active regions establish the larger-scale surface inflow(e.g., Haber et al. 2004; Hindman et al. 2009). Searches for the solar giant cells of convection may
also benefit through a better understanding of the depth dependence of the observed north-southalignment of supergranules, thought to arise from the organizing effect of deep-seated, larger-scale
motions (Lisle et al. 2004; Featherstone et al. 2006). Moreover, a detailed, unified description of
subsurface vorticity would undoubtedly benefit attempts to predict active region flaring using ring
analysis which is currently accomplished with low-resolution 16◦ measurements (e.g., Komm et al.
2010; Maurya & Ambastha 2010)
The self-consistent incorporation of multiple resolutions of ring-analysis measurements into a
unified solution describing the solar flow can thus be highly advantageous, but it is nontrivial. For agiven analysis-region size, complementary measurements are produced using different radial modeorders with different spatial wavenumbers. These measurements represent unique spatial averages
over the solar subsurface flows, and velocity maps obtained from fits to the different modes will
necessarily disagree with one another. Moreover, the same modes from different tile sizes will
have different horizontal averaging functions. A simple averaging together of results from different
analysis regions sizes is thus impossible. To that end, we have developed a 3-D inversion techniquewhich we have termed Adaptable Resolution Ring Diagram Inversions (or ARRDI), that allows
us to simultaneously invert data from a number of different modes and region sizes to obtain aself-consistent map of the horizontal flow over the full disk of the Sun and for a range of depths.
We proceed now with a detailed discussion of the measurement technique in §5.2, followed by adiscussion of the calculation of the sensitivity kernels §6.1. The inversion algorithm is discussed in§7.1, and its averaging properties in §7.2, followed by a comparison of the results realized using our
new method and previous methods in section 8.1.

85

5.2

Ring-Analysis Technique

5.2.1

Obtaining Doppler Measurements

Modern Doppler measurements of the entire solar surface are generally obtained throughone of two world-class observatories. For space-based observations, local helioseismic techniques
rely primarily on the Michelson Doppler Imager (MDI; Scherrer et al. 1995) aboard the Solar andHeliospheric Observatory (SOHO). MDI measures Doppler shifts and line-of-sight magnetic fieldusing five narrow filters spanning the wings and core of the Ni 6768 Å line. Using a CCD arraywith 1024×1024 pixels, MDI generates Dopplergrams and magnetograms with a 4” resolution (1.45
Mm per pixel at disk center in the Sun) at a 1 min cadence. The primary drawback to observations
made aboard SOHO, which is positioned at the Earth-Sun L1 Lagrange point, concerns telemetry
limitations to the data transfer. Full-resolution Dopplergrams are returned for only two to threemonths of every year, and lower resolution Doppler data is transferred during the remainder of theyear.
The data transfer limitations of MDI can be partly overcome through measurements madeusing Doppler data from the Global Oscillations Network Group (GONG; Harvey et al. 1996),
comprised of six ground-based telescopes sited around the world, with a combined duty cycle ofaround 90%. GONG employs 1024x1024 CCD arrays to generate full-disk Dopplergrams of the
Sun using the same Ni line as MDI. Dopplergrams are generated with a 2.5” resolution and 1 min
cadence, and line-of-sight magnetograms are captured every 20 min.
Helioseismic investigations in the near future will rely heavily on spaced-based observation
with the Helioseismic and Magnetic Imager (HMI; Kosovichev et al. 2007). HMI is aboard the Solar
Dynamics Observatory (SDO) and will return full-disk Dopplergrams and vector magnetogramswith 1” resolution at a 45 s cadence using a 4096×4096 ccd array. Measurements of Doppler
shifts and magnetic fields obtained using this instrument will use the FeI line at 6173 Å, which ismore magnetically sensitive (i.e. higher Lande factor) and thus better suited for the generation ofvector magnetograms than the Ni line used in MDI and GONG. The ring-analysis inversion package

86
described in Chapter 5 has been developed for integration into the HMI data pipeline.

5.2.2

Measuring Flows with Ring Analysis

Our ring-analysis measurements are made using Doppler-velocity data obtained with theMichelson Doppler Imager (MDI) aboard the Solar and Heliospheric Observatory (SOHO). Eachanalysis operates upon a sequence of images constructed by tracking a square patch of the Sun,16◦ in size, for 1664 min (27.7 hr) with a 60 s cadence. The tracking slides the analysis region acrossthe solar surface at the Snodgrass (1984) differential rotation rate corresponding to the central
latitude of the patch. This tracking largely removes the prominent signal associated with thesurface differential rotation of the Sun, allowing us to sample the smaller-scale convective patternsof interest. In total, 189 such patches, overlapping their nearest neighbors by 7.5◦ are used to tile
the solar disk, covering the central portion of the disk within 60◦ of disk center. As it is tracked, thedata is remapped to a Cartesian grid using Postel’s projection. This projection has the property
that great circles passing through the center of the analysis region are remapped as straight lines,
and distances along those lines are preserved, thus minimizing distortions in the power spectrum
and enabling the measurements to be made at higher latitudes than would otherwise be possible(Haber et al. 1995; Bogart et al. 1995).
The tracked mosaic of 16◦ regions is then further subdivided into similarly overlapping mosaics of 4◦ and 2◦ tiles. The pixel resolution of the MDI instrument is 0.125◦ at the solar equator.
The horizontal dimension for these different tiles is thus a power of 2 (in pixels) in each case,enabling their efficient computation through fast Fourier transforms. Alternative tiling schemes
could also be employed involving intermediate sized tiles of 8◦ (for example), but we have foundthat the higher resolution 4◦ measurements provide good depth coverage down to about 10 Mm.
Beyond this depth our primary interests have concerned the large-scale flows around active regions,which we find to be well sampled using 16◦ measurements. Higher resolution tiles, such as 1◦ and0.5◦ could also be employed, but fitting frequency shifts in these smaller power spectra is difficult,owing in part broadening of the ridges in frequency as we discuss further in §6.1.3, and may require

6

2

a

87

b

ky (Mm-1)

1
0
-1

4

-2
-2

-1

0
1
kx (Mm-1)

2

-1

0
1
kx (Mm-1)

2

3
2
2

c

1
ky (Mm-1)

Frequency (mHz)

5

1

0
-1

0.5

1.0 1.5k (Mm-1)

2.0

-2
-2

Figure 5.1: Structure of a typical solar power spectrum. (a) Example of a 16◦ power spectrum
as a function of frequency and horizontal wavenumber. Power has been averaged over angle inwavenumber space, with regions of high power indicated by dark tones, and low power by lighttones. Power is concentrated along curved ridges. (b) When cut at a constant frequency, theseridges appear as concentric rings, shown here at 3 mHz. (c) Subsurface flows produce a shift ofthese rings about their centers, as demonstrated here using a 1.5 km s−1 eastward flow.

the development of alternative fitting algorithms. For the time being, we thus proceed with a setof measurements made at the 16◦ , 4◦ and 2◦ resolutions.
Each tracked data cube is apodized spatially and temporally. The temporal window functionsmoothly tapers the duration to 1612 min, while the spatial window function is circularly symmetricand reduces the analysis region to disks that are 15◦ , 3.75◦ , and 1.875◦ for the 16◦ , 4◦ , and 2◦ tiles
respectively. Power spectra (in space and time) are then computed from the apodized images. Suchpower spectra possess curved ridges of power in wavenumber and frequency, which when cut at aconstant frequency appear as concentric rings of power (Fig. 5.1). Flows passing beneath one ofthese patches of the Sun produces a Doppler shift in the wavefield, causing these rings to shift
about their centers. As described in Schou et al. (1998), at each total wavenumber k in the power
spectrum, we fit the observed power P about a single ridge of power to a model that accounts for

88
the frequency splittings produced by an underlying flow,
P =

b
AΓ
+ 3.
2
2
k
(ω − ωn + u · k) + Γ(5.1)

Here A is the amplitude of the mode, Γ is the half width, ωn is the resonant frequency of the mode,and b is a smoothly varying background. The wavenumber vector k is given byk = kx êx + ky êy ,(5.2)

where êx and êy are unit vectors in the east-west and north-south directions respectively. Similarly,
the velocity vector u is given byu = ux êx + uy êy ,(5.3)

where ux and uy represent a spatial average over the underlying horizontal flow in the east-west
and north-south directions. A single velocity measurement u is thus characterized by 5 parameters:
it’s location on the solar disk (i.e., the central longitude φ and latitude λ of the analysis region),
the size T of the analysis region, the mode order n, and the wavenumber k. For the remainder ofthis chapter, we adopt the notation ui as shorthand forui ≡ u(ki , ni , Ti , φi , λi ),(5.4)

with the subscript i denoting a unique combination of these measurement parameters. By generating a set of measurements ui , we can then construct a flow map over the solar disk for each modeand tile-size combination.

5.2.3

Interpretation of the Fitted Velocities

Combining the different measurements ui into a unified description of the solar subsurface
flow, which we call v, is nontrivial in that different modes average the horizontal flow differently
with depth. The details of this averaging are described for each ui by a sensitivity kernel tensor
Ki such thatux,i =XZj

Kxj,i (r)vj (r)dr

(5.5)

89
and
uy,i =XZ

Kyj,i (r)vj (r)dr,

(5.6)

j

where r is the position vector, and where we have adopted the shorthandZ
Z ∞Z ∞Z ∞dr ≡dx dy dz
−∞

−∞

(5.7)

−∞

to indicate an integral taken over all space. In practice, the off-diagonal elements of the kerneltensor K are small in relation to the diagonal elements (i.e. Kxy ≤ 0.05 Kxx ) and are insensitive toa horizontally uniform flow (see Birch et al. 2007). We thus proceed by neglecting these off-diagonal
elements, thereby yielding the more succinct form of Equations (5.5) and (5.6),Zux,i = Kxx,i (r)vx (r)dr,

(5.8)

and
uy,i =

Z

Kyy,i (r)vy (r)dr.

(5.9)

A set of measurements ui thus defines a system of integral equations that may be inverted
to solve for the underlying flow field v. Such inversions have previously been carried out using the16◦ data but not for the 2◦ and 4◦ data. Rather, for these data sets, surface flow maps have been
generated using High-Resolution Ring Analysis (HRRA). As noted in Hindman et al. (2006), thesmall range of wavelengths accessible to the small-tile measurements suggests that the different
measurements using only the f -mode ridge sample the upper convection zone in a similar manner.
An average in wavenumber over f -mode measurements is thus a sensible representation of theupper convection zone. Flow maps with 2◦ and 4◦ resolution generated using HRRA make no useof p-mode data, however, and thus do not provide a similarly high-resolution description of deeper
flows.
Traditionally, the inversion of 16◦ measurements has been accomplished through the use of
linear inversion methods wherein the velocity solution at the jth point in the inversion grid, vj ,may be written as a linear combination of the data. Specifically,vx,j =

Xiax, ji ux,i ,(5.10)

90
where the weights ax,ji are determined through the inversion procedure. Combining Equations(5.8) and (5.10), we see thatvx, j =Z

Cx, j (r)vx (r)dr,(5.11)

withCx, j (r) =Xax, ji Kxx,i (r).

(5.12)

i

The function Cx, j is known as the averaging kernel and describes how the solution vj , whenconstructed according to Equation (5.10), samples the underlying flowfield.
Two different approaches to building the Cx, j ’s have been applied to the 16◦ measurements.
The first of these, Regularized Least Squares (RLS; e.g., Craig & Brown 1969; Haber et al. 2000)
finds a solution vj that fits the measurements well, while sacrificing the best possible fit to maintainsome desired degree of smoothness of the solution. Inversions utilizing Optimally Localized Averages (OLA) have also been applied extensively to the analysis of 16◦ data (e.g. Haber et al. 2004).
In an OLA scheme, a set of ax,ji ’s is generated that minimizes the difference between the averaging
kernel Cx,j and some desired target function (typically a Gaussian in depth) while keeping the error
associated with vj below some desired value. The solution vj found through an OLA inversion may
not best fit the data, but it can yield results that are easier to interpret. For example, Haber et al.
(2004) have shown that the use of an OLA inversion method can result in averaging kernels thatlack the negative sidelobes in depth that often occur in RLS averaging kernels.
As applied to ring analysis, both RLS and OLA inversions have traditionally been 1-D innature, yielding a solution as a function of depth beneath each tile. These inversions have employedhorizontally invariant sensitivity kernels, whose structure in depth has been taken to be the kinetic
energy density of the associated mode (e.g., Haber et al. 2002). A particular limitation of thesemethods is that a 1-D inversion performed for a particular tile makes no use of information fromneighboring tiles that sample a common portion of the solar disk. This can be overcome through thedevelopment of 3-D sensitivity kernels for all tile sizes that describe both the horizontal and depth
sensitivity of a measurement to a subsurface flow. We will discuss the calculation of such 3-D kernels

91
shortly, but turn first to a characterization of the ring-analysis measurements themselves.

5.2.4

Characterization of the Data and its Errors

Subsurface flow maps realized through linear inversion of ring-analysis measurements are
intimately linked to both uncertainties in the individual measurements as well as to variations in
the density of measurements across the disk. A characterization of the measurements and theirerrors will thus facilitate a better understanding of the inversion results. For ring-analyses, thefitted velocities are those which minimize the difference between the data and the fitting functionof Equation 5.1. These fitted parameters thus correspond to a minimum in parameter space about
which a Taylor expansion is performed to determine the local curvature of this minimum. This
curvature (in the ux or uy direction) is then taken as a representative error for our measurement.
Using ring-analysis measurements generated for the 72 days spanning 1 March 2001 through 26
May 2001, we have calculated average errors and data filling factors both mode by mode and byposition on the disk for each of the three analysis-region sizes.
The presence of magnetic fields, which may suppress convection and produce asymmetries
in the wavefield, can also make the fitting of power spectra difficult by introducing perturbationsin the power spectrum not accounted for in Equation (5.1). For instance, active regions tend toabsorb acoustic power (e.g., Braun 1995), leading to a suppression of power in portions of the power
spectrum corresponding to waves travelling away from a region of strong magnetic field. Effects suchas this can significantly alter power spectra for the 4◦ and 2◦ regions, which are generally noisierthan their 16◦ counterparts owing to the smaller amount of Dopplergram data that is averaged in
the power spectrum computation. Moreover, even in the absence of a magnetic field, the frequency
linewidths increase with decreasing region size, leading to broader, less well defined power ridges
for the smaller analysis regions. These enhanced uncertainties in the fidelity of the power spectrum
translate into greater uncertainties in the measurements as the analysis-region size is decreased.
Figure 5.2 illustrates the sensitivity of 2◦ measurements to a magnetic field using data from
27 March 2001. The unsigned magnetic flux density, degraded to the resolution of the 2◦ data, is

92

a

b

c

Latitude (Degrees)

40
20
0
-20
-40

-40

-20

0

20

0

40
400

Flux Density (G)-40 -2002040Carrington Longitude (Degrees)
011
Data Count

-40

-20

0

20

50

40
160

Error (m s-1)

Figure 5.2: The effect of magnetic field strength on ring-analysis errors for 2◦ measurements from
27 March 2001. (a) Magnetic flux density degraded to the resolution of the 2◦ measurements. Dark
tones denote high field strength, and light tones low field strength. (b) Data measurement count,averaged over all attempted modes. High and low values indicated as in (a). Positions where no
successful measurement was made are indicated in red. (c) Ring-analysis errors averaged over all
attempted modes. Coloring as in (b). Regions of high field strength correlate with higher error anda lower measurement success rate.

shown in Figure 5.2a. Regions of strong flux density tend to have a lower number of successful
measurements (Fig. 5.2b). As indicated by the red coloring in Figure 5.2, regions of particularly
strong magnetic field are often devoid of velocity measurements. The average measurement errors(Fig. 5.2c) are often high in regions of strong magnetic field, with the imprint of the active region
belts clearly visible in the associated error mapping.
This general trend persists for the 4◦ measurements as well. In Figure 5.3a we show the diskand mode-integrated error, averaged in time over 72 days, relative to the magnetic field strength.
Magnetograms for each day have been degraded to the resolution of the different analysis-region
sizes, as reflected in the fact that the 16◦ measurements do not possess corresponding fields strengths
as high as those associated with the 4◦ and 2◦ measurements. Measurement errors are relatively
insensitive up to field strengths of ∼ 50 G. The measurement success rate drops accordingly as
power spectra become more difficult to fit in the presence of strong magnetic fields (Fig. 5.3b).

Average Error (m s-1)

150

Measurement Success Rate

93
a

100

50

0
10

100
Flux Density (G)

1.2

b

1.0
0.8
0.6
0.4
0.2
10

100
Flux Density (G)Figure 5.3: Average effects of magnetic field strength on measurements made using different
analysis-region sizes. (a) Measurement error (averaged over all modes and positions on the disk)
versus magnetic flux density for the 16◦ (red), 4◦ (blue), and 2◦ (black) measurements. 16◦ measurements are relatively insensitive to the magnetic field strength. 4◦ and 2◦ measurements exhibit
higher uncertainties in regions of strong field and lower measurement success rates.

We note that even in the presence of weak magnetic field, only about 40% of 2◦ measurements are
successful.
Measurement errors exhibit considerable variation between modes as well. We plot the timeand disk-averaged measurement errors for each mode in Figure 5.4 where each symbol-color combination denotes an average error for a particular tile size and mode order. Errors are plotted relative
to the lower turning point of their associated modes. In this case, the lower turning point for eachmode (where the phase speed and the sound speed are equivalent) has been determined using thesound speed profile from solar model S (Christensen-Dalsgaard et al. 1996). The evanescent f
modes have no lower turning point, and so their errors have been plotted relative to the e-folding
depth of each mode. Plotting in this manner is advantageous in that it provides a sense of theeffective penetration depth of the different measurements as well as the typical uncertainties at
those depths.
This information can be particularly useful when deciding on how the horizontal resolution of
our inversion grid should change with depth. For example, the deepest penetrating 2◦ measurement(triangles) is for a p1 mode (red) with a lower turning point of 6 Mm. However, the error is quite
high for 2◦ measurements beyond a depth of about 3 Mm. We might then expect an inversion

94

300

Average Error (m s-1)

250
200
150
100
50
0
1

10
Lower Turning Point (Mm)

16 degree
4 degree2 degree

f
p1
p2

p3
p4
p5

p6

Figure 5.4: Average ring-analysis errors as realized with MDI data. A representative error for eachtile size, frequency, and horizontal wavenumber combination was generated by averaging errors
from 72 days of ring measurements made using MDI dopplergrams from 1 March 2001 through 26
May 2001 over the full solar disk. Tile size and mode order are indicated using the coloring andsymboling scheme indicated. Errors have been plotted as a function of the lower turning point of
the associated mode. Errors are lowest for the 16◦ tiles. Within a given tile size, errors are lowest
at mid-ridge and increase with increasing radial mode-order.

95
algorithm that is aware of the errors to emphasize the 2◦ measurements when building the solution
near the surface and down to depths of 3-4 Mm. Below this depth, the 4◦ data will likely be favored,
and a natural degradation of the inversion’s horizontal resolution can occur. Morever, our deepest
measurements occur in the 16◦ dataset and only sample the upper 35 Mm of the solar convection
zone, thus setting a lower limit on the depths which we may hope to probe.
Figure 5.4 shows that within any combination of analysis-region size and mode order, theerrors exhibit a systematic variation with lower turning point. For a give mode order, the lower
turning point is inversely proportional to wavenumber. Errors thus become high at the lower
and upper ends of a particular ridge, with a minimum somewhere in between. At mid-ridge, thefrequency is closer to 3 mHz, where the solar power spectrum peaks. At the ends of the ridges,
the power is then lower than the mid-ridge and background noise can make the power spectrum
more difficult to fit. This effect is compounded at the high-wavenumber end of the ridge where
frequency linewidths become larger, and the power becomes more spread out in frequency space,thus reducing the ratio of signal to noise.
The time- and disk-averaged data-filling factor for the different modes is plotted in Figure
5.5 using the same symbols as Figure 5.4. The 16◦ measurements demonstrate a near 100% successrate irrespective of their associated errors. The 4◦ and 2◦ measurements, however tend to possesshigh filling-factors where errors are low and vice versa. Thus, in addition to having high errors
beyond about 3 Mm, few 2◦ measurements are typically available for an inversion algorithm tomake use below this depth. We see from Figures 5.4 and 5.5 that, when choosing a horizontal
resolution for our inversion grid, we might degrade it twice relative to the surface: once below 3
Mm and again below about 11 Mm where the 4◦ measurement success begins to decline. Before we
discuss the specifics of our inversion algorithm any further, however, we turn to an examination of
the development of the sensitivity kernels and their variation with the measurement parameters.

96

Data Filling Factor

1.0
0.8
0.6
0.4
0.2
0.0
1

10
Lower Turning Point (Mm)

16 degree
4 degree2 degree

f
p1
p2

p3
p4
p5

p6

Figure 5.5: Data filling factor over the solar disk for ring-analysis measurements for different modes
and different tile sizes. Tile size and mode order are indicated as in Figure 7.4. Measurements
attempted for the 16◦ data are generally successful, with some data dropouts occuring at high n.
Success rates for the 4◦ and 2◦ measurements behave similarly to the errors, with mid-ridge andlow n measurements faring the best.

Chapter 6

Properties of the Sensitivity Kernels

6.1

Development of the Sensitivity Kernels
We now outline the procedure by which a ring-analysis measurement ui can be related to ahorizontal flow v beneath that corresponding patch of the sun through the sensitivity kernel Ki by
ux,i =

Z

Kxx,i (r)vx (r)dr(6.1)

uy,i =

Z

Kyy,i (r)vy (r)dr.

(6.2)

and

The general framework for the calculation of sensitivity kernels has been established in Gizon &Birch (2002) and Birch et al. (2004), and its application to ring-analysis kernels is described in
Birch et al. (2007). Our aim here is to present a synthesis of these three papers so that a complete
description of the model and calculation is readily available when we discuss the properties of theresulting sensitivity kernels in §6.1.3.
Given that we already have two velocity variables defined (u and v), we choose to work in thefluid displacement vector ξ(r, t) when describing waves in the solar interior. The surface velocity
projected along the line of sight, our primary observable, and which we choose to call φ, is thengiven byhiφ(x, t) = A(x, t) ℓ̂ · ξ̇(x, zobs , t) .

(6.3)

Here zobs is the height at which the measurements are made, ℓ̂ is a unit vector directed along theline of sight, x is the horizontal position vector, and the overdot on ξ denotes a time derivative.

98
The function A describes the spatial and temporal windowing function applied to some patch of thesun before a power spectrum is constructed and a ring-analysis velocity measurement performed.
The observed power spectrum P is then given by
P (k, ω) = φ(k, ω)φ∗ (k, ω),(6.4)

where we have adopted the convention that a function written explicitly in terms of ω and k denotes
the temporal and horizontal Fourier transform, namelyφ(k, ω) ≡Z Z

φ(x, t)e−iωt eik·x dk dω.

(6.5)

In the presence of a flow, a perturbation to the power spectrum, δP , is introduced such that
P (k, ω) = P0 (k, ω) + δP (k, ω),(6.6)

where the subscript “0” to indicates quantities in the absence of a flow. The presence of a flow
similarly produces a perturbation δφ to the surface velocity field so that, to first order in δφ, δP isgiven by
δP (k, ω) = δφ(k, ω)φ∗0 (k, ω) + δφ∗ (k, ω)φ0 (k, ω) = 2Re [δφ(k, ω)φ∗0 (k, ω)](6.7)

Ring-analysis techniques derive a velocity from δP by fitting the power spectrum at a specifick and n as described in Equation (5.1). The fitting described in §5.2 is nonlinear, however, and if
we wish to perform a linear inversion, a linear approximation to our fitting technique is needed tocompute the sensitivity kernels. Birch et al. (2007) have demonstrated that a linear fitting method
given byui =Z ω0 +∆ωi Z 2π
ω0 −∆ωi

W i (ki , ω)δP (k i , ω)ki dθ dω,

(6.8)

0

with weights W i defined byWi ∼∂P0 (ki , ω)êi,
2P0 (ki , ω)
∂ω

(6.9)

yields results that agree well with our nonlinear fitting method for measured velocities up to 400 m
s−1 . Here êi is a unit vector in the directon of ui , ∆ωi indicates the frequency range that is used

99
when fitting a ridge of radial order ni at a particular wavenumber ki , and θ indicates the angle thewavenumber vector k makes with respect to the kx -axis so thatki = ki cosθ eˆx + ki sinθ eˆy(6.10)The only quantity in Equation (6.8) that depends on the presence of a flow is δP . Computing thelinear sensitivity of a measurement ui to the underlying flow v requires finding an expression forδP of the formδP (k, ω) =Zv(r) · f (r, k, ω)dr,(6.11)

where f (r, k, ω) is some as yet to be determined function. Combining Equations (6.8) and (6.11)
then yields (for a measurement in the x-direction),ux,i =Zvx (r)Kxx,i (r)dr,(6.12)

W x,i (k , ω)fx (r, ki , ω)ki dθ dω.

(6.13)

with the sensitivity kernel Kxx,i given by
Kxx,i (r) =Z ω0 +∆ω Z 2π
ω0 −∆ω

0Computation of the sensitivity kernel Kxx,i thus requires a model for the wave field in the absenceof a flow, φ0 , and for the perturbations to that wave field introduced by a steady flow. We turn
now to the development of such a prescription of the solar wave field which may then be used to
calculate f (r, k, ω).

6.1.1

Modeling the Solar Wavefield

Our ring-analysis measurements are made using modes with a number of different horizontal
wavelengths, the largest of which (obtained with 16◦ data) is ∼ 0.03 of the solar radius. As eventhe largest wavelengths considered are quite small compared to the solar radius, we adopt a planeparallel treatment of the solar interior for our model, neglecting any effects of sphericity. We begin
by considering perturbations about an adiabatic background state in hydrostatic balance and lookfor wave solutions in the displacement ξ. We denote thermodynamic background state (horizontally

100
invariant) variables by an overbar and perturbations about this state by a prime symbol, such thatthe density ρ is given by¯ + ρ′ (x, y, z),ρ(x, y, z) = ρ(z)(6.14)and so too for the pressure P . Here z denotes the depth, and x and y should be associated withthe east-west and north-south directions respectively. Upon carrying out this linearization andneglecting the presence of any large-scale flows, the momentum equations becomeρ̄ξ̈ = −∇P ′ + ρ′ g + ρ̄∂t (Γξ) + S(6.15)

ρ̄g = ∇P̄ .

(6.16)

ρ′ = −ξ · ∇ρ̄ − ρ̄∇ · ξ,

(6.17)

P ′ = −ξ · ∇P̄ + c̄2 (ρ′ + ξ · ∇ρ̄).

(6.18)

andThe continuity equation is given by

and energy equation by

Here g denotes gravity, γ is the ratio of specific heats, and c the local sound speed. We use solarmodel S (Christensen-Dalsgaard et al. 1996) to describe the background stratification and neglect
the variation of gravity with depth since the modes we consider sample the upper 30 Mm of theSun (see Figure 5.4). The function S describes the excitation of waves by near-surface turbulent
convection, and the operator Γ describes their damping.
We adopt phenomenological models of the damping and source functions that have been
tuned to reproduce the observed power spectra. For the damping, we adopt the model of Birchet al. (2004), which has been tuned mode-by-mode to reproduce the frequency linewidths of highdegree (kR⊙ ≥ 300) modes as measured by Korzennik et al. (2004) using MDI data. As noted in
Gizon & Birch (2002), a prescription for the source function itself is never needed, but rather one
for its covariance. Following Birch et al. (2004), we assume a vertical momentum source whose

101
covariance M has the form2

M (k, z ′ , z) = e−(ωTsrc ) δ(z − zsrc )δ(z ′ − zsrc ).

(6.19)

This form of the source covariance, along with a source depth zsrc of 100 km and source correlation
time Tsrc of 68 sec were found by Birch et al. (2004) to yield model power spectra that matchedwell with observations. We have assumed that the sources are uncorrelated spatially, resulting in
the independence of M on the horizontal wavenumber k. Our ring-analysis measurements employ
modes with wavelengths as small as 4.5 Mm which is still much larger than the measured correlation
length for solar granulation of 450 km (Title et al. 1989), making the spatially uncorrelated
approximation reasonable.
These four equations can be combined into a single operator L describing the response of ξto the source S,

Lξ = ρ̄ξ̈ − ∇[γ P̄ ∇ · ξ + ξ · ∇P̄ ] + (∇ · ξ)∇P̄ + ξ · ∇(∇P̄ ) + ρ̄∂t (Γξ) = S.

(6.20)

When a large-scale flow v is present, the effects of advection must be included in Equation (6.15),
and a perturbation to the wave operator δL becomesδLξ = 2ρ̄v · ∇ξ̇,

(6.21)

which generates a corresponding perturbation to the wave field δξ. Denoting the absence of a flowby the subscript “0”, Equation (6.20) becomes(L0 + δL)(ξ 0 + δξ) = S 0 ,(6.22)

where, for simplicity, we have neglected the effect of advection on S. Expanding into zero- andfirst-order contributions, Equation (6.22) becomes to lowest order,L0 ξ 0 = S 0 ,(6.23)

L0 δξ = −δLξ 0 ≡ δS.

(6.24)and to first order,

102
We see from Equation (6.24) that the advection of the zero-order wavefield serves as the source
function for the perturbed wavefield, and have denoted this quantity by δS accordingly. The
solutions ξ 0 and δξ can be constructed from the same set of vector Greens function Gj for theoperator L0 , which we define as
L0 Gj (x − s, t − ts ) = êj δD (x − s)δD (t − ts ).

(6.25)

Here δD is the Dirac delta function. The Green’s function Gj represents the velocity resulting froman implusive force at location s and time ts in the eˆj direction. We construct the Gj ’s through asummation over the normal modes of the system described by Equations (6.15)-(6.18) as detailed
in Birch et al. (2004). The zero-order solution ξ 0 (and so too the zero-order power spectrum) can
be constructed from the Greens function of the operator L0 ,XZ ∞ Z ∞ξ 0 (r, t) =Gj (r − r ′ , t − t′ )Sj (r ′, t′ )dr ′ dt′ .
−∞j(6.26)−∞

From Equations (6.24) we see that advection acting on the zero-order wavefield serves as the source
for the perturbed wavefield, so thatδξ(r, t) =XZ ∞ Z ∞j−∞

Gj (r − r ′, t − t′ )δSj dr ′ dt′ .

(6.27)−∞The perturbed wavefield δξ is linearly related to the steady flow v through δS, and so we now
have the pieces in place to construct the function f (r, k, ω) of Equation (6.11) and the sensitivity
kernel Ki,x .
6.1.2

Assembly of the Sensitivity Kernels

We now proceed with construction of the sensitivity kernels by rewriting δP in terms of ξ bycombining Equations (6.3) and (6.7) to yield∗δP (k, ω) = −2Re[ω δξ˜z (k, ω, zobs )ξ̃0,z(k, ω, zobs )] ≡2Zv(r) · f (r, k, ω)dr,(6.28)

where ℓ̂ is assumed to be in the êz direction for simplicity. We have further adopted the convention
that an overhat on a quantity indicates its convolution with the apodization function A so thatZ Z
˜
ξz (k, ω) =A(k − k′ , ω − ω ′ ) ξz (k′ , ω ′ ) dk ′ dω ′ .
(6.29)

103
Taking the Fourier transform of δξz , we findXZδξ z (k, ω, zobs ) =Gz,j (k, ω, zobs − z)δSj (r, ω)eik·xdr,(6.30)

jand combing Equation (6.24) with Equation (6.30), we haveZ˜δξz (k, ω, zobs ) = v(r) · F̃ (r, zobs , k, ω)dr,(6.31)

whereXρ̄(z) eik·x Ĝz,j (k, ω, zobs − z) ∇ξj,0 (r, ω).

(6.32)hi∗f (r, zobs , k, ω, zobs ) = −2Re ω 2 ξ˜z,0(k, ω, zobs )F̃ (r, zobs , k, ω) ,(6.33)F (r, zobs , k, ω) = −2j

The function f is then simply

with the sensitivity kernel for a measurement in the êx direction to a horizontal flow vx given byZ ω0 +∆ω Z 2π
Kxx,i (r) =Wx,i (ki , ω)fx (r, ki , ω)ki dθ dω.
(6.34)
ω0 −∆ω

0

Many properties of this complex formulation for K are better understood by examining
its 3-D structure, as we next do in §6.1.3. However, there are two important points to note
concerning Equation (6.33). First, the effect of advection manifests in this formulation throughthe presence of the gradient operator which operates on the zero-order wavefield. Second, thesensitivity kernel depends on the apodizing function A through the presence of ξ˜z in both f andin the weights Wi (which themselves depend on the zero-order power spectrum). Thus the spatial
and temporal windowing needed to construct our power spectra can lead to substantial averaging
of modal power in spectral space. A measurement made at a particular wavenumber ki will thus
sample the underlying flowfield differently as the analysis region becomes smaller. We turn now
to a more complete discussion of these effects and of the general 3-D structure of the sensitivity
kernels.

6.1.3

Structure of the Sensitivity Kernels

The 3-D structure of sensitivity kernels corresponding to 16◦ measurements of the 3 mHz
f and p1 modes has been examined by Birch et al. (2007). These kernels were found to be

104

a

b

Latitude (degrees)10

1.1
Sensitivity

5
0
-5
-10
-10

-5
0
5
10
Longitude (degrees)

0.5
0.0
-0.5
-10.0
-3.3
3.3
10.0
Longitude (degrees)

Sensitivity

0

0.7

c

d

3

3
Depth (Mm)0Depth (Mm)

0

6

9
-10

6

-5
0
5
10
Longitude (degrees)

9
-0.5

0.0
0.5
Sensitivity

1.1

Figure 6.1: Sample sensitivity kernel for a 16◦ measurement with radial order n = 2 and horizontal
wavenumber k = 1.1 Mm−1 (ℓ = 765). (a) Horizontal cross section of the kernel at a depth of 4.5
Mm. Bright orange tones indicate high sensitivity, and dark tones low. (b) Horizontal cut throughthe kernel center at 4.5 Mm. (c) Vertical cross-section cut through the kernel center. (d) Verticalcut through the center of the kernel (black) along with kinetic energy density of the associatedmode (red). Depth dependence of the sensitivity kernel largely mirrors that of the mode whose
measurement it represents.

105
horizontally smooth, with a vertical structure reflecting the kinetic energy density of the measured
mode. Our inversion technique requires a kernel set corresponding to a wide range of wavenumbers(200 ≤ ℓ ≤ 1100), radial orders up through n = 6, and including kernels for three tile sizes: 2◦ ,4◦ and 16◦ . We have found that the inversion kernels have properties that vary significantly between
tile sizes and along a given radial order ridge. As an awareness of some of these features is usefulin interpreting the inversion results, we turn now to a detailed discussion of the sensitivity kernel
structure.
We first note an important limitation in the calculation of our kernels. As indicated in Equation (6.3), our ring-analysis sensitivity kernels should vary with the central latitude and longitude of
the observed patch of the sun. The sensitivity kernel calculation is highly time-consuming, however,requiring approximately n + 1 days to compute a single sensitivity kernel for a measurement madefor radial mode order n using a modern desktop computer. The computation of a full kernel set(using 3 observation tile sizes and modes through p6 ) for measurements made at just one positionon the solar disk requires about one month of computing time using 8 cores. The kernel calculation
code is written in MATLAB and is not optimized to run on parallel architectures. Computing afull kernel set for each location on the disk currently require many months running on a handful ofcores. We have thus taken ℓ̂ to be directed in the ẑ direction, as would be appropriate for a measurement made at disk center (and assuming a plane-parallel geometry), and treat our sensitivity
kernels as though they have no variation across the solar disk.
An example of a typical 16◦ kernel is shown in Figure 6.1. This particular kernel is for ameasurement with radial order n = 2 and horizontal wavenumber k = 1.1 Mm−1 . The verticalstructure of the kernel is well approximated by the kinetic energy density of the associated mode
(as computed using our model). Horizontally, the structure and extent of the sensitivity kernel islargely determined by the analysis-window size. The horizontal variation of sensitivity across thekernel is minimal with the exception of a noticable rim of enhanced sensitivity at the outer edge ofthe kernel. These edge effects originate from Gibbs ringing at the edges of the apodization window.
As noted in Birch et al. (2007), this rim becomes more enhanced as the frequency range over which

106

Figure 6.2: Variation of sensitivity kernel structure with horizontal wavenumber. (upper row )Horizontal cuts of the sensitivity kernels for several wavenumbers along the p3 ridge along with(lower row ) vertical cuts through the center of the sensitivity kernel (black). At high wavenumbers,
the kernels develop enhanced rings of sensitivity along the edges. Vertical sensitivity also begins todepart from the kinetic energy density (red) at high wavenumbers.

107
the fit is carried out is increased. As the frequency range available to the kernel calculation isincreased, power from a wider range of modes is available for the reconstruction of the apodization
function profile, and the edges become sharper.
Our standard ring measurement technique encompasses a range of frequencies extending frommode center to 40% of the way to the next ridge on either side of the line center. Most sensitivity
kernels that we have computed thus exhibit an enhanced sensitivity near the edges as seen in thehorizontal structure of Figures 6.1 a,b. This basic horizontal structure exhibits little variation withdepth (Figure 6.1c), and the 16◦ kernels are thus nearly separable. The frequency spacing between
neighboring ridges in the power spectrum increases with horizontal wavenumber, and so too does
our fitting range. At higher wavenumbers, the sensitivity kernel is constructed from an increasingly
wider range of smaller length-scale modes, leading to more prominent rims at high wavenumber.
This effect is visible in Figure 6.2 (top row) where we show cross-sections of sensitivity kernels forp3 measurements.
The vertical structure of the 16◦ sensitivity kernels also exhibits some variation with wavenumber, developing noticable departures from the kinetic energy density of the mode at the higher
wavenumbers (lower row, Fig. 6.2). The reason for this behavior lies with the finite nature of
the observing window (i.e. the spatial apodization) whose effect is to blend power in horizontal
wavenumber by an amount inversely proportional to the size of the analysis region. This effect is
illustrated in Figure 6.3 where we have plotted our model power spectrum at 3 mHz along with
the Fourier transform of the 2◦ and 16◦ apodization windows. The averaging in wavenumber resulting from a 16◦ measurement tends to smear power across a distance only slightly larger than
the width of a ridge. For the 2◦ measurements, this effect is much larger and can blend powerfrom several distinct ridges. This effect is exacerbated at higher wavenumbers where the slope of1

a power spectrum ridge goes as ∼ k− 2 and the ridge flattens out, and an averaging over somefixed range in wavenumber will average more power. The broadening of the sensitivity peaks at
depth with increasing wavenumber in Figure 6.2 is thus a reflection of the fact that a ring-analysis
measurement taken at a specific wavenumber is in fact sensitive to both shallower and more deeply

108
penetrating modes.
1.0

1.0

a

0.8

b

0.8

Power

Power

0.6
0.4

0.6

0.4
0.2
0.2

0.0
-0.2
0.4

0.6
0.8
1.0
1.2
Horizontal Wavenumber (Mm-1)

1.4

0.0
2

3

4
Frequency (mHz)

5

6

Figure 6.3: Effects of a finite observing region on the observed power spectrum. (a) Model power
spectrum at 3 mHz (black) along with the Fourier transform of the window functions for 2◦ measurements (blue) and 16◦ measurements (red). Substantial averaging of power in wavenumber
occurs for small analysis regions. (b) Model power spectrum (black) cut at horizontal wavenumber
k = 1.39 Mm−1 (ℓ = 967) after convolving with a circular apodization window for the 2◦ (blue)and 16◦ (red) measurements. Averaging in wavenumber leads to an effective smearing of power infrequency.

As indicated in Figure 6.3, for 2◦ measurements the apodization leads to substantial averaging
of power in wavenumber. A representative sensitivity kernel for the 2◦ measurements is shown in
Figure 6.4. This kernel is for the same radial order and wavenumber as that shown in Figure 6.1 (p2 ,k = 1.1 Mm−1 ) but exhibits a vertical structure very different from the kinetic energy density andpossesses lobes of negative sensitivity near the surface. Negative sidelobes near the surface such as
this can result from a measurement being sensitive to power from shallower modes. For example,
the observed p2 ridge is actually constructed by some average that includes (among other modes)
both p2 and p1 modes. A near surface flow in the x-direction might shift the p1 ridge towards p2 on
the positive kx -direction, but away from it negative kx -direction. After averaging in wavenumber,the observed p2 ridge could peak at a lower frequency in the positive kx -direction, but remain
unaffected in the negative kx -direction. This asymmetric shift in the central frequency of the ridge
between the two sides of the power spectrum, induced by a shallow flow in the positive x-direction,

109

a

b
1.1

1.5

Sensitivity

Latitude (degrees)

3.0

0.0
-1.5
-3.0
-3.0 -1.5 0.0 1.5 3.0
Longitude (degrees)

0.5
0.0
-0.5
-3

-1
1
3
Longitude (degrees)

Sensitivity

-0.2

0.6

c

d

3

3
Depth (Mm)0Depth (Mm)

0

6

6

9
-3 -2 -1 0 1 2 3Longitude (degrees)

9
-0.5

0.0
0.5
Sensitivity

1.1

Figure 6.4: Sample sensitivity kernel for a 2◦ measurement with radial order n = 2 and horizontal
wavenumber k = 1.1 Mm−1 (ℓ = 765). (a) Horizontal cross section of the kernel at a depth of 4.5
Mm. Bright orange tones indicate high sensitivity, and dark tones low sensitivity. (b) Horizontal
cut through the kernel center at 4.5 Mm. (c) Vertical cross-section cut through the kernel center.
(d) Vertical cut through the center of the kernel (black) along with kinetic energy density of theassociated mode (red). Owing to apodization effects, the kernel structure can deviate considerablyfrom the kinetic energy density of the mode for the 2◦ analysis regions.

110
0

0.0

a

0.0

b

0.0

c

d

-2
-2.5

-2.5

-2.5

-5.0

-5.0

-5.0

-7.5

-7.5

-7.5

Depth

-4

-6

-8

-10
-0.5

0.0 0.5
Sensitivity

1.0

-10.0
-0.5

0.0 0.5
Sensitivity

1.0

-10.0
-0.5

0.0 0.5
Sensitivity

1.0

-10.0
-0.5

0.0 0.5
Sensitivity

1.0

Figure 6.5: Variation of the vertical structure of a 2◦ sensitivity kernel for the p2 mode withhorizontal wavenumber of 0.26 Mm−1 (ℓ = 181). Horizontally integrated sensitivity is shown inblack, and the accompanying kinetic energy density (normalized to 1) is shown in red. The kernel
has been computed using a zero-order wavefield with (a) all modes present, (b) only p2 modes, (c)p2 and p1 modes, and (d) p1 , p2 , and p3 modes present. Averaging in wavenumber, arising from thefinite spatial extent of the observations, can mix power from different mode orders, substantiallyaltering the vertical structure of the sensitivity kernels for the 2◦ measurements.

would be interpreted as a flow in the negative x-direction by our measurement technique.
We can illustrate this effect on the vertical structure of our 2◦ kernels by removing various
mode orders from the zero-order solution. In Figure 6.5, we show a 2◦ p2 kernel calculated by (a)including all modes, (b) by removing all modes but the n = 2 mode, by (c) including only n = 2
and n = 1 modes, and by (d) including n = 1, 2 and 3 modes. The vertical structure of the kernel
clearly reflects the kinetic energy density of the target mode when no other modes are present,
and shows varying degrees of departure from this structure as additional modes are included. It
is unavoidable that any measurement technique, ring analysis or otherwise, that employs a power
spectrum constructed from a finite data set will inherently incorporate such mode-averaging effects
into the measurements. These effects, while potentially undesired, are not necessarily a weakness

111
of the technique, provided that the sensitivity kernels account for their existence, as we have donehere.

Chapter 7

Attaining 3-D Inversions of Ring Analysis Data

7.1

Elements of the Inversion Procedure
With sensitivity kernels in hand, we now turn to the inversion of the measurements. Ourring analysis technique involves measurements made at numerous wavenumbers and overlapping
regions of three sizes on the solar disk. In total, one 26 hour period of observation with MDI can
produce data sufficient to make approximately 3 × 105 ring measurements. We have found that
inverting about one-ninth of the solar disk at once is computationally tractable. A single inversion
can in principle thus involve roughly 3 × 104 measurements. As discussed in §5.2, this collection of
measurements defines a system of integral equations describing the flow field,ui =Z

Ki (r)v(r)dr,(7.1)

that must be solved simultaneously. For brevity, we have here dropped the use of the subscripts x
and y under the assumption that the x-velocity measurement ux,i is assumed to be independent of
a flow in the y-direction, vy , and vice versa. Equation (7.1) may then apply to either measurement
direction without loss of generality. Discretizing Equation (7.1) leads to a set of linear equationsui =X

Kij vj drj .

(7.2)

j

In order to simplify the representation of the linear system, we employ a 1-D decomposition of our
3-D inversion grid, using the index j as shorthand for a single location in the solution grid so that
Kij ≡ Ki (xj , yj , zj ).

(7.3)

113
Here drj denotes the corresponding volume weights used to evaluate the integral. Adopting theconvention that matrices are indicated with Roman font, 1-D matrices by lower case letters, and
2-D matrices by upper case letters, we can writeu = Av,(7.4)

Aij = Kij drj .

(7.5)

with

7.1.1

Regularized Least Squares Approach

Given N measurements, we can then in principle invert A to solve for the solar flow fieldat N points. However, an additional complication lies in the inherent error in the measurements
caused by uncertainties in the measurement procedure so that Equation (7.4) becomes(u + δu) = A(v + δv).

(7.6)

Spurious errors δu in u can lead to departures δv from the true solution v. Dependingon the kernel structure, an erroneous measurement ui can “pollute” the entire solution v, rather
than just a single value vj . Moreover, the presence of errors can lead to inconsistencies between
measurements with overlapping sensivity kernels, effectively rendering the matrix singular. We
thus require additional information to stabilize the solution. Regularized Least Squares (RLS; c.f.
Craig & Brown 1969) is one such approach that involves both our estimate of the errors on u andour intuition concerning the expected structure of the solution, namely that it should be smooth
to some degree. Rather than solving Equation (7.4) directly, the RLS approach seeks a trade-off
between a good fit of the solution to the data and the smoothness of the solution by minimizing
the functional
Ψ=

X 1 
i

σi2

ui −

Z

Ki (r)v(r)dr

2

+λ

Z

|Dv(r)|2 dr.

(7.7)

Here, the σi ’s are the estimated errors on the measurements. The operator D can be any functionaldescribing the smoothness of the solution. In our case, we take D to be the gradient operator
Dv ≡ ∇v.

(7.8)

114
Other smoothing operators are possible as well (the second derivative is often employed), but
we have chosen the first derivative for simplicity and owing in part to the complexity of the 1dimensional data representation of the three-dimensional grid. Here λ is a scalar parameter used
to tradeoff between the goodness of fit and the smoothness of the solution. The discretized formof Equation (7.7) becomes

Ψ=

X 1
i

σi2



ui −

X
j

2

Aij vj  + λ

X X
j

k

Djk vk

!2

,

(7.9)

with D indicating the discretized form of the gradient operator D. We can find the solution v thatminimizes Ψ by differentiating Ψ with respect to vj . Adopting the shorthand that
Bij ≡

Aij
ui
, di ≡ ,σi
σi(7.10)

we find an equation for v of the form
(BT B + λDT D)v ≡ Mv = BT d.

(7.11)

v = M−1 BT d,(7.12)The solution v is then given by

and the solution may be obtained by inverting the matrix M. As the dimension of the kernel matrix
B is nsolution × nmeasurements , the dimension of M is nsolution × nsolution , where nsolution is thenumber of points in the solution grid and nmeasurements likewise indicates number of observations.
We obtain the inverse of the matrix M using the parallel algorithms available through the Scalapack
Fortran package. These inversion routines employ an LU-decomposition whose computation timescales as n3solution , and it becomes important to carefully select the horizontal and vertical resolutionof the inversion grid when solving for the flow beneath large patches of the sun. The selection ofan appropriate gridding scheme is discussed further in section 7.2.
The matrix M contains the associated kernels, weighted by their errors, for all measurements
included in the inversion. We select which measurements to include in M by defining an inversion
region of arbitrary spatial extent and then assemble M using all measurements whose associatedLatitude (Degrees)

115
15

A

B

10
5
0
0

10
20
Longitude (Degrees)

30

40

Figure 7.1: Illustrating the need for a modified RLS inversion to bring adjacent inversion regions
into agreement. Inversion regions “A” and “B” are indicated by the red and blue dashed-dot lines
respectively. Solid circles denote the bounds of 16◦ measurements. A standard RLS inversion
carried out in region “A” cannot account for the overlapping information from region “B” (shadedin gray).

sensitivity kernels lie entirely within that region. This approach can be problematic in light of theoverlapping nature of our measurements.
For example, suppose we wish to carry out an inversion of two adjacent regions comprisedof two overlapping 16◦ measurements each. This situation is illustrated in Figure 7.1 where tworegions are indicated in red (region “A”) and blue (region “B”), and the bounds of the two 16◦ measurements (which are apodized to 15◦ ) for each region are indicated by the circles. The 7.5◦ overlap
between adjacent measurements leads to a similar overlap between the two inversion regions. Information from the measurement centered at 15◦ longitude in “A” cannot be included directly in
the inversion for region “B” because the associated sensitivity kernel would extend beyond thebounds of the inversion region “B”. Depending on the noise characteristics of the measurements,
the solution for “A” may not agree with that from “B” in the region of overlap.
This disagreement can be problematic when analyzing flows that span multiple subregions.
One solution to this problem is to simply invert the measurements for the full disk at once. However,
our experience has shown that the inversion of a square subregion greater than 45◦ in size can require
memory in excess 32 GB. Moreover, inversions for such large regions can take in excess of 1 day tocompute. Such computations can be difficult to accomplish using a modern desktop with a typical116
core count (4-8) in some time period less than the observing time. To this end, we have devised an
iterative method for bringing inversion results from overlapping subregions of the disk into close
accord with one another. As the resulting flowmap results from a quilting-together of its individual
subregions, we refer to this technique as a stitched inversion and turn now to a description of its
implementation.

7.1.2

Stitched Inversions

The stitched inversion method approximates the solution to the standard RLS problem (Eq.
7.7) by subdividing the inversion domain into several smaller subdomains and inverting for theirunderlying flows as would normally be done. A second set of inversions is then carried out iteratively
where the mismatch between adjacent regions in their overlapping interval is penalized through an
additional regularization term.
Let us assume there are N subregions of the solar disk whose underlying flows we wish
to self-consistently obtain by inversion to yield a full-disk solution. Let the subscript α denote
an individual subregion with α = 1, 2...N and the index n denotes the iteration from which thesolution results so that vαn is the solution for subregion α at the nth iteration. An initial solution,vα0 for each subregion is first obtained by solving our standard RLS equation
Mα vα0 = BT α dα .

(7.13)

Subsequent solutions vαn are generated by introducing an additional regularization term
Rαn to Equation (7.13). This term penalizes the difference between the new solution vαn and theprevious solution vαn−1 at points where the different subregions overlap. We take Rαn to be
Rαn = ΛXX
(vj,αn − vj,β n−1 )2 ,(7.14)

β6=α jwhere Λ is an adjustable regularization parameter, and the sum in j is taken over all points common
to regions α and β. Adding Rαn to the least squares Equation (7.13) and carrying out the usual

117
200
vx (m s-1)

a

10-1

100

b

100
0

0

-100

-50

-200

10-2

-100

e

100

50

vx (m s-1)

vx (m s-1)

d

10-3

0

0

-100
10

-4

0

-50

-200
2

4
6
Iteration

8

10

c

50

vx (m s-1)

100

-100
180

200
220
240
Longitude (Degrees)

260

180

200
220
240
Longitude (Degrees)

260

Figure 7.2: Convergence properties from a full-disk inversion using the iterative stitching method.
In this case, the flow field for most of the solar disk was reconstructed using 9 square sub-regions
37.5◦ in horizontal extent overlapping, by 7.5◦ with their nearest neighbors. (a) Stitching regularization R-value integrated over the full disk (relative to its value after the zero-order inversion).
Convergence is achieved within three to five iterations. (b, c) Non-stitched horizontal velocity at
one latitude for three overlapping sub-regions (central region colored in red) at 0.5 Mm and 9 Mm
in depth respectively. (d, e) Horizontal flow for these same subregions resulting from the converged
iterative stitching procedure. Good agreement in the region of overlap is reached at both depths
following convergence.

differentiation with respect to vj leads to a new set of equations
Qα vαn ≡ (Mα + ΛXNαβ )vαn = BT α dα + Λ

β6=α

X

Nβα vβ n−1 .

(7.15)

β6=αHere the Nαβ ’s are diagonal matrices whose elements are defined as

Njj, αβ =



 1


 0

if coordinate j in region α also occurs in region βotherwise.

Other prescriptions for Rαn are possible. In particular, one might weight Rαn by the errors
on vα n−1 from the previous iteration. However, the formulation of Equation (7.14) yields matrices
Qα that are independent of the iteration n. We are thus required to perform only two matrix
inversions for each region: the initial inversion of the matrix M, Equation (7.13), and the inversionof the matrix Q, Equation (7.15). Subsequent solutions are generated by updating the right-handside of Equation (7.15) and multiplying by Q−1
α .

118

-10

-20

vx (m s-1)

0

10

0

-20

140

150

160
170
Longitude

180

50

140
40.0 m/s

150

160
170
Longitude

180

c

0
-50
100

-10

vx (m s-1)b10LatitudeLatitude

a

100

50

d

0
-50
-100
140

150

160
170
Longitude

180

Figure 7.3: Comparision between an iterative stitched inversion (a) and its companion full inversion(b) at a depth of 9Mm (with horizontal resolution interpolated by a factor of two). Blue lines denoteboundaries of the four overlapping subregions used in the iterative inversion. Each square subregion
was 22.5◦ in horizontal extent. (c) Sampling of results along the latitude line indicated in red
in panel (a) for the full inversion (black symbols) and zero-order results for two overlapping subregions of the iterative inversion (black and red lines). (d) Same as in c, but displaying the results
for the two sub-regions after the stiched solution has converged. Following the iterative inversion,flows in both regions show good agreement with the full inversion.

In practice the solutions converge quickly, typically requiring three to five iterations. In
Figure 7.2 we illustrate the convergence properties for a full disk inversion carried out using thistechnique. The disk was tiled with nine 37.5◦ subregions with nearest neighbors overlapping by 7.5◦ .
An integrated Rαn for all nine regions was calculated over 10 iterations, with convergence evident
around the fourth iteration. Figure 7.2b, d illustrates the pre- and post-convergence properties ofthe solution. Solutions from neighboring regions typically agree well at the surface following thezero-order inversion. Solutions in the upper layers are determined largely by the 2◦ measurements.
Unlike the 16◦ measurements, very few of these smaller measurements are omitted near the edgesof the inversion and neighboring inversions contain redundant information near the surface. At the
9 Mm depth, however, the missing information from the 4◦ and 16◦ measurements is importantand neighboring inversions agree well only after the iterative procedure has converged.
We are unable to invert an entire disk to check the validity of this technique. However,
we have tested our iterative results on a smaller scale by decomposing a 37.5◦ region into four
22.5◦ regions overlapping by 7.5◦ . A comparison between the results for the iterative stitched119
inversion relative to the full inversion at 9 Mm in depth is shown in Figure 7.3. The two flow fields
show good agreement with one another. A plot of vx at constant latitude (see red line Fig. 7.3a andFigs. 7.3c, d) shows that good agreement both between neighboring regions and the full inversion
is reached once the procedure has converged.

7.2

Averaging Properties of the Inversion Algorithm
The solution values vj are linearly related to the measurement values ui so that each solution

point vj may be written as
vj =X

aji ui ,(7.16)iwhere the coefficients aji are found through Equation (7.12) withaji =1 X −1
Mjℓ Aiℓ .
σi2(7.17)ℓ

The solution is thus a weighted sum of the data whose details depend on the measurements errors
and the inverse matrix M−1 . The error on the solution, ǫj , can be written in similar fashion interms of the measurement errors σi as
ǫ2j =

X

a2ji σi 2 .

(7.18)

i

Note that the aji ’s do not depend on the measurements themselves, only upon the measurement
errors and the sensitivity kernels. Inversions made using different measurement values, but withidentical errors and identical mode sets, will yield a result by averaging the raw data in exactly thesame manner. Interpretation of the results generated using the 3-D inversion procedure requires
an understanding of the aji ’s and the localized averaging of measurements they represent. We thus
devote the next section to describing in detail the averaging properties of our inversion algorithm.

7.2.1

Data Combination Properties of the Inversion Algorithm

We can gain a sense of how the inversion scheme combines data by using the set of average
MDI errors described in §5.2.4 in an inversion. We have carried out an inversion of a 15◦ × 15◦ region

120
with a horizontal resolution of 0.9375◦ and a grid extending to 25 Mm in depth, assembling theinversion matrix M by assuming that all modes are fit with an uncertainty given by the average
error for that mode. We find that a convenient parameterization of the relative contribution of thedifferent data points to the solution is given by the quantity Aj , which we define as
Aj (n, k, T ) ≡P

2i δD (T − Ti , k − ki , n − ni ) aji
P 2
.
i aji

(7.19)The quantity Aj is thus a measure of the disk-integrated contribution, mode-by-mode, of thedifferent measurements to the solution point vj .
The mode-by-mode contribution, Aj , for the solution vj at the horizontal center of theinversion domain and at four different target depths is illustrated in Figure 7.4. For a given target
depth in the inversion grid, the principle contribution to the result generally comes from modes
whose lower turning points are within one Mm or so of that depth. The associated analysis region
size (indicated by the different symbols) of the highest Aj ’s are typically those that offer the finest
available horizontal resolution at that depth. Near the surface, the solution is comprised primarilyof 2◦ f -mode measurements. The influence of the 2◦ measurements extends to approximately 3 Mm,with the 4◦ measurements beginning to contribute appreciably beyond this depth. Such behavior is
to be anticipated in light of Figure 5.4, where we see that the 2◦ errors begin to rise sharply beyond
about 3 Mm, thus making the 4◦ data at that depth more favorable to the RLS scheme. The modes
employed in the 4◦ measurements in turn penetrate only as deeply as 13 Mm (not shown). At thedepth of 15 Mm then, the 16◦ measurements are primarily responsible for forming the solution.
The changing contributions from the different tile sizes means that the horizontal resolution of
the inversion solution naturally degrades with depth (even in the absence of varying horizontal
resolution of the inversion points with depth). This natural change in resolution is evident uponinspection of the averaging kernels resulting from the inversion.

121

Target Depth : 0.50 Mm

Target Depth : 3.00 Mm

0.300

0.075

0.225

0.056

0.150

0.038

0.075

0.019

0.000
0

1

2
3
4
Lower Turning Point (Mm)

5

0.000
0

Target Depth : 7.00 Mm

4
6
8
Lower Turning Point (Mm)10Target Depth : 15.00 Mm

0.050

0.100

0.038

0.075

0.025

0.050

0.013

0.025

0.000
0

2

5
10
15
Lower Turning Point (Mm)

20

0.000
0

16 degreefp3

4 degreep1

p4

2 degree

p2

p5

5

10
15
20
Lower Turning Point (Mm)25

p6

Figure 7.4: Depth-averaging properties of the inversion scheme. The Aj (n, k, T )’s for a solution
point at the horizontal center of the inversion grid, at four target depths, are plotted as a functionof the lower turning point of the contributing modes. Tile sizes are indicated by distinct symbols,and different radial orders by the indicated colors. At target depths close to the surface, 2◦ measurements with a shallow lower turning point are largely used to construct the solution. At mid
depths, 4◦ measurements dominate with contributions from the 2◦ and 16◦ data. Deeper in theconvection zone, 16◦ measurements comprise the bulk of Aj .

122

1.5 Mm
5.0 Mm

0.5
9.0 Mm

13.0 Mm

0.0

0.8

b

0.0 Mmc

0.6
Integrated Sensitivity

1.0aIntegrated Sensitivity

Integrated Sensitivity

1.5

0.4

0.2

0.6
0.4
0.2
0.0

0.0
-0.5

-0.2
-6

-4

-2 0 2Longitude

4

6

-6

-4

-2 0 2Longitude

4

6

-6

-4

-2 0 2Longitude

4

6

Figure 7.5: Horizontal structure of the averaging kernels. (a) Variation of the horizontal structurewith target depth. (b) Horizontal averaging kernel structure for a target depth of 5 Mm using an
inversion grid with horizontal resolution at that depth of 0.9375◦ (black line), 1.875◦ (red line), and3.750◦ (blue line). The averaging kernels achieved are substantially degraded when the horizontal
resolution of the inversion is below that of the measurements. (c) Kernel structure for an inversion
point lying between two measurements (red). The kernel is much broader than one that is achieved
when the inversion point coincides horizontally with a measurement (black). The vertical blackline indicates the target longitude.

7.2.2

Averaging Kernels

The averaging kernels for a given solution point are generated from the aji ’s as described by
Equation (5.12). Thus, as the Aj ’s begin to favor measurements from one tile size over another, theeffective width of the averaging kernels changes accordingly. We examine the horizontal structureof the averaging kernels in Figure 7.5a. Near the surface, these averaging kernels are roughly 2◦ inwidth, with negative sidelobes on either side of the maximum. These sidelobes arise from theinversion subtracting measurements on either side of the target longitude and latitude, and adding
those measurements that fall onto the target location. This results in a narrower central peak, but
one with small negative sidelobes on either side. At mid-depths, the averaging kernels are roughly
4◦ in width, with smaller negative sidelobes. Beyond 13 Mm, the horizontal structure is reflectiveof the 16◦ measurements.
The horizontal structure of the averaging kernels also depends on the horizontal resolution of
the inversion grid. As the resolution is degraded, sensitivity kernels are binned accordingly in thehorizontal direction, effectively resulting in a spatial uncertainty in their associated measurements.

123
We examine the effect of varying the horizontal resolution in Figure 7.5b, where we plot horizontal
averaging kernels for the solution point at the domain center and a depth of 5 Mm. The averaging
kernel realized when adopting a horizontal resolution of 1.875◦ at this depth is indicated in blackand that for a grid with 0.9375◦ resolution in red. At this depth, the 4◦ measurements are thepredominant contributors to the averaging kernel. These observations are spaced apart by 1.875◦ .
The coarser grid thus has the same horizontal resolution as the measurments, and the finer gridhas one interpolating point between measurements. The averaging kernel resulting from the finer
grid is less sharply peaked about the solution point than that resulting from the coarser grid, butpossesses less pronounced sidelobes on either side of the central peak. However, the differences
between the two are small.
The finer horizontal mesh also allows for the reconstruction of the flow at a point halfway
between two measurements. The averaging kernel for such an interior solution point at a depth of
5 Mm is shown in Figure 7.5c (red). The adjacent averaging kernel for a solution point coinciding
with a measurement is shown in black. With no 4◦ measurement centered at the target longitude,
the averaging kernel formed is essentially an average of those for the two measurements made on
either side of the solution point. The resulting kernel is thus much broader than one targeted at thelocation of a 4◦ measurement, and yields information that is more easily obtained through simple
interpolation of the lower resolution inversion results. When using our current tiling scheme, thebenefits of adopting a higher horizontal resolution than the measurements themselves can provide
at a given depth are thus minimal and unlikely to outweigh the computational costs associatedwith increasing the number of solution points.
We explore the depth structure of various averaging kernels in Figure 7.6a. The averaging
kernels tend to peak near their target depth, decaying below that point. Side lobes of positive
sensitivity near the surface are apparent for many of the kernels sampled, increasing in strength as
the target depth approaches the surface. These kernels also become wider at deeper target depths.

124
2.0

a

0.0 MmSensitivity

1.5
1.0

1.0 Mm
2.0 Mm

0.5

4.0 Mm

6.0 Mm

9.0 Mm

0.0
-0.5
0

5

10

15

Depth (Mm)10bKernel Width (Mm)Center of Gravity (Mm)

20
15

10

5
0
0

5

10
15
Target Depth (Mm)

20

25

c

8
6
4
2
0
0

5

10
15
Target Depth (Mm)20

25

Figure 7.6: Depth properties of the averaging kernels. (a) Sampling of the vertical structure of
the averaging kernels. (b) Center of gravity (solid line) achieved at different target depths. For
reference, a slope of one is indicated by the dashed line. (c) Kernel width at the different target
depths. Averaging kernels become less localized with increasing depth. It is difficult to create well
localized averaging kernels below a depth of 15 Mm when using the observations described in §5.2.

We quantify the vertical location of the kernel by the center of gravity, zcg , which we take to beR|C(r)| z dr,
zcg = R|C(r)| dr(7.20)

with z the depth below the surface. The center of gravity for averaging kernels at the different
target depths is shown in Figure 7.6b. We see that the center of gravity is generally somewhat lessthan the target depth (except near the surface). Moreover, beyond 15 Mm the inversion is unable
to build kernels centered more deeply. The vertical width of the kernels (which we define as the
distance between the 0.25 and 0.75 percentile points) increases steadily with depth as well. Kernels
are about 1 Mm in width near the surface, and about 8 Mm in width beyond the 15 Mm depth.

125
1.00

a

b

c

d

e

f

0.75
0.50
0.25
0.00
1.00
0.75
0.50
0.25
0.00
0

5

10
15
20
Target Depth (Mm)

0

5

10
15
20
Target Depth (Mm)

0

5

10
15
20
Target Depth (Mm)25

Figure 7.7: Variation of the tile-size contribution parameter B(T, z), with changing regularizationin (a-c) the horizontal direction and (d-f ) the vertical direction. 2◦ measurements are indicated inblack, 4◦ data in blue, and 16◦ data in red. Changing the regularization alters the range of depths
over which measurements from the different tiling sizes are emphasized by the inversion algorithm.

7.2.3

Effects of Regularization

The inversion results are dependent on the value assumed for the regularization tradeoff
parameter λ. It can be convenient to apply one value of λ, hereafter λh , to the horizontal portionof the regularization operator D, and a separate value to the vertical part, hereafter λv . As theregularization is increased, the solution becomes smoother, and the averaging kernels broader.
Consequently, data is combined differently as the degree of regularization is varied. Some insight
into effects of varying the regularization can be gained by examining the parameter B, given by
!
P
2
X
δ
(T
−
T
)
δ
(z
−
z
)
a
D
i
D
j
ji
i
P
B(T, z) ≡,(7.21)
2i δD (z − zj ) aji
jwhich measures relative contribution of measurements obtained from the three different analysisregion sizes to the solution at each depth.
The variation of B with changing horizontal regularization is demonstrated in Figure 7.7a-c,
and with changing vertical regularization in Figure 7.7d-e. As vertical regularization is increased,
the solution becomes more averaged in depth, and the regions of contribution for the different
tile sizes become less distinct. For sufficiently large vertical regularization, the influence of the

126
10

a
rms Error (m s-1)

rms Error (m s-1)

10
8
6

10-1

4
2
0
0

1
10
5

10
15
20
Target Depth (Mm)

25

b

8
6
10-2

4
2
0
0

10-1
1

5

10
15
20
Target Depth (Mm)25

Figure 7.8: Variation of 3-D inversion errors with changing horizontal and vertical regularization(a) Variation of rms error with depth using an λv = 0.1 and different horizontal regularizationsof λh = 10 (dashed), λh = 1 (solid) and λh = 0.1 (dotted). (b) Inversion errors realized using a
λh = 1.0 and different vertical regularizations of λv = 1 (dashed), λv = 0.1 (solid), and λv = 0.01
(dotted). As regularization is increased, the inversion scheme averages over more data and theresulting errors tend to decrease. Errors at all depths are typically less than 10 m s−1 .

2◦ measurements is felt down to 15 Mm, even though the deepest lower turning point associatedwith any of the 2◦ measurements is 5 Mm. The most noticable effect of increasing the horizontal
regularization is to increase the usage of the 4◦ and 16◦ measurements at increasingly shallower
depths. These effects can be better understood through examination of the associated averaging
kernels.
The effect of these different regularization parameters on the inversion errors is illustrated in
Figure 7.8, and their effect on the averaging kernels is shown in Figure 7.9. Increasing regularizationin either direction generally decreases the inversion errors. This is accomplished by averaging overmore data, which in turn leads to broader averaging kernels as can be seen in Figure 7.9. There we
see that increasing the horizontal regularization leads to increased horizontal width in the averaging
kernels. As the horizontal regularization is increased, the inversion scheme seeks to horizontallysmooth the solution. This is accomplished near the surface, for instance, by forming the solution at
one location by averaging 2◦ measurements over a broader range in longitude and latitude. Further
smoothing is accomplished by incorporating more 4◦ and 16◦ measurements into the near-surface
solution, as is apparent in the top row Figure 7.7. A similar pattern is followed at the mid-depths
where fewer 4◦ measurements and more 16◦ measurements are incorporated into the solution as
the horizontal regularization is increased.

127
Vertical regularization works similarly. As the vertical regularization increases, the averaging
kernels broaden in depth. The horizontal structure of the averaging kernels shows little change,
however, as λv is varied by a factor of one hundred. Averaging of the data is thus performed over a
wider range of modes as vertical regularization is increased, but not over a wider spatial extent on
the disk. This averaging over modes is why the 2◦ data is used at 15 Mm for λv = 1 as shown in
Figure 7.7. There, the solution at depth is smoothed vertically in part by incorporating the shallow
modes accessible only to the 2◦ measurements.
With no known flow to compare our results against, the choice of regularization parameters is
unfortunately somewhat arbitrary, even after consideration of the data-combination and averagingkernel implications. While the errors do vary with the regularization, they are not incredibly
sensitive to it, remaining largely within the 5-10 m s−1 for all parameter combinations explored
here. We choose to adopt values of λh = 0.1 and λv = 0.1 (Fig. 7.8a, dotted line), basing thischoice on the reasonably focused kernels they provide in the horizontal and in the vertical. This
choice also produces flowmaps that compare well with 1-D inversion of the 16◦ measurements, as
discussed in the next section. Morever, for this parameter combination, the deepest effects of the2◦ and 4◦ measurements roughly correspond with the lower turning points of their most deeply
penetrating modes.

128

1.00

1.00

a

0.75
Sensitivity

Sensitivity

0.75
0.50
0.25
0.00

-6

-4

-2

0
2
Longitude

4

0.25

-0.25
0

6

1.00

c

5

10
Depth (Mm)

15

20

5

10
Depth (Mm)

15

20

d

0.75
Sensitivity

0.75
Sensitivity

0.50

0.00

-0.25

1.00

b

0.50
0.25
0.00

0.50
0.25
0.00

-0.25
-6

-4

-2

0
2
Longitude

4

6

-0.25
0

Figure 7.9: Sensitivity of the averaging kernels to variations in the regularization for an averaging
kernel targeted at the 5 Mm depth. (a,b) Variation of the (a) horizontal and (b) vertical structureof the averaging kernel with changing horizontal regularization. λv has been taken to be 0.1. λhof 0.1, 1.0, and 10.0 are indicated by dashed, solid and dotted lines, respectively. (c,d) Variationof the (c) horizontal and (d) vertical structure of the averaging kernel with changing verticalregularization and λh of 1.0. λh of 0.01, 0.1, and 1.0 are indicated by dashed, solid and dottedlines, respectively. Variations in λv change the vertical structure of the averaging kernels whileleaving the horizontal structure relatively unchanged. Similarly, λh primarily affects the horizontal
structure of the averaging kernels.

Chapter 8

Exploring Convection and Magnetism with ARRDI Inversions

8.1

Comparison of ARRDI Results with Prior Assessments
Prior to the development of our 3-D inversion algorithm, flow maps generated from ring-

analysis measurements have relied on 1-D (depth) RLS or OLA inversion of the 16◦ data (e.g.,
Haber et al. 2004), and on f -mode fits from 2◦ and 4◦ measurements (e.g. Hindman et al. 2006).
As ARRDI combines all three levels of data in a self-consistent fashion, we should see points ofsimilarity between solutions generated using 1-D inversions and those generated using ARRDI. Thisis particularly true at key depths where the 3-D inversion tends to emphasize some small subsetof the data when constructing the solution. A comparison between solutions generated using 3D and 1-D techniques will thus serve to highlight its advantages over and differences relative tothe previous schemes. We devote this section to a comparison of the ARRDI results with thosegenerated using f -mode fits and 1-D OLA inversions. We employ MDI data from 24 April 2001 forthese comparisons. The presence of a particularly large and long lived active complex, NOAA AR
9433, in this dataset allows us to compare flows in regions of both quiet and active sun.

As discussed in §7.2.3, we adopt regularization parameters λh = 0.1 and λv = 0.1 for thiscomparison. In depth, our inversion grid extends from the surface to a depth of 25 Mm, with 0.5
Mm resolution in depth near the surface transitioning to 3 Mm resolution near the bottom of thegrid. Horizontal resolution is 1.25◦ near the surface and transitions to 3.75◦ at depth. A complete
list of target depths and their horizontal resolution is given in Table 8.1. Horizontally, we tile thefull disk with a 3×3 matrix of square 45◦ regions overlapping by 7.5◦ in longitude and latitude and130

Figure 8.1: Comparison of the ARRDI results to 2◦ f -mode fits using data from 24 April 2001.
(a) 3-D inversion solution at 1 Mm and (b) solution resulting from 2◦ f -mode fits averaged inwavenumber. Cellular structures are similar in the two flow maps, though the inversion results
exhibit higher velocities (possibly due to negative sidelobes in the sensitivity kernels). White
lettering indicates locations of the smaller target regions shown in Figs. 8.2 (and at other depths
in Fig. 8.4).

131
TABLE 8.1
Inversion Grid
Target Depth (Mm) Horizontal Resolution (degrees)
0.0
1.25
0.5
1.25
1.0
1.25
1.5
2.50
2.0
2.50
3.0
2.50
4.0
2.50
5.0
2.50
6.0
2.50
7.0
2.50
9.0
2.50
11.0
3.75
13.0
3.75
17.0
3.75
19.0
3.75
22.0
3.75
25.0
3.75

Table 8.1: Target depths and horizontal resolution of the inversion grid for the 3-D inversions of
the April 24, 2001 data.

combined with the stitching algorithm.

8.1.1

Surface Comparisons with the 2◦ data

In the upper 2 Mm of the convection zone, the f -mode measurements dominate the data (seeFigure 5.4). At a depth of 1 Mm, where our horizontal resolution is 1.25◦ , we thus expect our 3-D
inversion to emphasize the 2◦ f -mode data in generating the solution. The ARRDI -generated flow
map at 1 Mm is shown in Figure 8.1a. An average magnetogram from that period is underlayed,with regions of strong positive magnetic field indicated in red, regions of strong negative field ingreen, and weak field in yellow. Velocity vectors are indicated by blue arrows. Convective patterns
associated with the largest supergranules (∼ 3◦ in size) are clearly visible at this depth and tend
to possess boundaries that trace out regions of strong magnetic field. For comparison the flow
field generated using 2◦ f -mode fits is shown in Figure 8.1b, where we have averaged together themeasured velocities over all wavenumbers fit along the f-mode ridge.

132
Convective patterns in the 2◦ fits are similar to those in the ARRDI flow map, but amplitudes
are noticably lower by about 33% for the fitted velocities. Averaging over the seven wavenumbers
fit along the f -mode ridge to generate the 2◦ map is primarily responsible for this effect. Moreover,each individual 2◦ data measurement has a sensitivity kernel that is 2◦ in width horizontally.
Surface flows realized through the 3-D inversion can in principle build a narrower horizontal averagethrough the appropriate combination of neighboring measurements as discussed in §7.2, leading toan amplification of the flows recovered through inversion relative to those obtained by averaging
f -mode fits.
A more detailed comparison of the two flowmaps is given by Figure 8.2 where we focus on
one quiet and two active sun regions labeled A-C in Figure 8.1. In Region A, there are several
locations indicated by small blue dots (i.e. zero length arrows) where no f -mode fits were obtained.
The ARRDI flowmap is able to reconstruct the flow at these locations using other mode orders,tile sizes, and where necessary, the horizontal regularization. This leads to an interesting effect inFigure 8.1b, where the supergranules bordering the active complex appear to possess stronger flows
moving away from the active region than they do into it. In many instances, this is due to theabsence of velocity fits rather than an actual decrease in velocity amplitude. Such an effect is notapparent in the ARRDI results. Regions B and C highlight the similarities of the two results forboth active and quiet Sun in detail. When data dropouts are not a factor, areas of convergence anddivergence correlate well between the two flowmaps, with only small differences between the two.
The correlation between these two data sets can be quantified using the Spearman rank correlation
coefficient ρ. A value of 1.0 for this coefficient indicates perfectly correlated data, -1.0 indicates
perfectly anti-correlated data, and zero indicates no correlation. We find the correlation coefficientfor flows in the east-west direction, ρew , to be 0.77 and that for flows in the north-south direction,
ρns , to be 0.74.

133

Latitude (degrees)

Latitude (degrees)

Region A

Region B
35

40

20

30

35

15

25

30

10

20

25

5

15

20

25

35

40

20

30

35

15

25

30

10

20

25

5
150

155

200 m/s

160

165

170

15
130

135
140
145
150
Carrington Longitude (degrees)aRegion C

25

20
160

165

170

175

180

b

Figure 8.2: Closeups of the regions of interest indicated by white lettering in Figure 8.1 for the(a, upper row ) ARRDI results and the (b, lower row ) 2◦ f -mode fits. Detailed examination of theflow structures in quiet sun (Region C) reveals good agreement between the two datasets. In more
active regions (such as A,B), data dropouts are present in the 2◦ f -mode data set (visible as small
blue dots). The 3-D inversion uses 4◦ and 16◦ degree data to reconstruct the solution where there
is an absence of 2◦ data.

134

Figure 8.3: Comparison of ARRDI solution to fitted data a depth of 4 Mm using ring-analysis
measurements taken on 24 April 2001. (a) Flows as realized through 3-D inversion.(b) Flows asrealized by averaging in wavenumber over velocity fits from 4◦ p1 measurements. As in Figure 8.1a,White lettering indicates locations of the smaller target regions shown in Figure 8.4. Unlike thecomparisons at 1 Mm, flow structures exhibit comparable magnitudes between the two data sets,
but differing morphologies.

135
8.1.2

Comparisons at Depth: Flowmaps and Their Relation to the 4◦ and 16◦ Data.

As discussed in §7.2, at the mid-depths of the inversion, the solution is heavily influenced
by the 4◦ measurements. In particular, Figures 7.4b and 5.4 suggest that at the 4 Mm depth, the4◦ p1 measurements are primarily responsible for the construction of the solution. In Figure 8.3, wecompare the ARRDI flowmap at 4 Mm with one generated by averaging the 4◦ p1 measurements
over all horizontal wavenumbers. Such a comparison can be more deceptive than those using thef -modes as the solution at this depth cannot be attributed to any one radial mode order. For
instance, the 2◦ p1 and 4◦ p2 measurements may contribute non-negligibly to the solution at 4 Mm
as well. Nevertheless, the overall patterning of the convection agrees well between the two flow
fields, albeit with some notable exceptions.
One of the more prominent differences between the two flow maps occurs in Region C, shownin detail in Figure 8.4. The presence of the strong northward flow visible in the data is not replicatedby the inversion procedure. In this instance, the ARRDI results are much weaker in amplitude andsomewhat smaller scale in nature. Agreement between the fits and the inversion results is better inRegions A and B, where the differences are largely associated with the manner in which adjacentregions of converging flow join to one another. However, we note that the 3-D inversion does not see
the prominent outflow from the sunspot in Region A (upper right). The Spearman correlation rank
coefficients for these two datasets are found to be ρew = 0.54, and ρns = 0.55. This is somewhat
lower than the f -mode comparisons, but some differences such as this are to be expected as ARRDI
incorporates measurements made from modes besides p1 when generating a solution.
At greater depths we see that the ARRDI results agree well with the 16◦ OLA inversions.
The resulting flow maps are compared at a depth of 14 Mm in Figure 8.5. For the ARRDI results,
we have plotted the flow at the target depth of 17 Mm where the center of gravity in depth forthe averaging kernels obtained is 14 Mm (see Figure 7.6). The outflow around the active region is
replicated well by the 3-D inversion, as are the large-scale westward flows that pass by the northernand southern boundaries of the complex. Velocity amplitudes are somewhat larger in the 3-D results

136

Latitude (degrees)

Latitude (degrees)

Region A

Region B
35

40

20

30

35

15

25

30

10

20

25

5

15

20

25

35

40

20

30

35

15

25

30

10

20

25

5
150

155

100 m/s

160

165

170

15
130

135
140
145
150
Carrington Longitude (degrees)aRegion C

25

20
160

165

170

175

180

b

Figure 8.4: Closeups of the regions of interest indicated by white lettering in Figure 8.3 for the (a,upper row ) 3-D inversion results and the (b, lower row ) 4◦ p1 -mode fits. Noticable differences exist
between the two analyses, such as the absence in the 3-D inversion of the strong northward outflow
from the sunspot in Region A.

137
a

b

Latitude (degrees)

40

30

20

10

0
120
130
30 m/s

140

150
160
Longitude

170

180

120

130

140

150
160
Longitude

170

180

Figure 8.5: Flowmaps at a depth of 14 Mm resulting from (a) the ARRDI algorithm and (b) 1-D
OLA inversion of 16◦ MDI ring data from 24 April 2001. Both solutions exhibit similar outflows
around the active complex, suggesting that the 3-D solution at depth is largely determined by the16◦ data. The 3-D inversion solution exhibits somewhat stronger flows at this depth, much as in
the near-surface solutions of Figure 8.1.

than in the 1-D inversions, and the ARRDI flowmap exhibits smaller scale horizontal structure thanthat achieved through the 1-D inversions. Here, as with the 2◦ comparisons, the 3-D inversion has
used the horizontal overlap of the measurements to build narrower horizontal averaging kernels,
yielding a flowmap with effectively higher resolution and larger resulting velocities. The correlation
coefficients here indicate good correlation between the 16◦ 1-D inversions and the ARRDI results,
with ρew = 0.71, and ρns = 0.80.

8.1.3

Establishment of Large-Scale Means

The discussion thus far has focused mainly on how the flow fields compare on the smaller
and more vigorous scales of convection. What about the large-scale means of these flows? In
particular, how do the flow fields obtained through 3-D inversion compare with the 16◦ scale atthat depth when averaged appropriately? In Figure 8.6a, we examine the large-scale means of the
ARRDI results at 1 Mm by smoothing the resulting flow map to the 16◦ scale. We have similarlysmoothed the flow map resulting from averaging in wavenumber over the 2◦ f -mode fits. While

138

a

16 Degree OLA (1.0 Mm)

3D Inversion (1.0 Mm)

2 Degree f-mode Fits

120 130 140 150 160 170 180

120 130 140 150 160 170 180
Carrington Longitude (degrees)120 130 140 150 160 170 180

Latitude (degrees)

40
30
20
10
0

20 m/s16 Degree OLA (4.0 Mm)

3D Inversion (4.0 Mm)

120 130 140 150 160 170 180

120 130 140 150 160 170 180
Carrington Longitude (degrees)

4 Degree p1-mode Fits b

Latitude (degrees)

40
30
20
10
0

20 m/s120 130 140 150 160 170 180

Figure 8.6: Establishment of the large-scale mean flows. (upper row ) Flowmaps at 1 Mm resulting
from 1-D OLA inversions of 16◦ data, ARRDI results smoothed to the 16◦ resolution, and smoothed2◦ f -mode fits, all using data from 24 April 2001. Inflows into the active complex are reproduced
when the ARRDI results are smoothed to the scale of the 16◦ data. Similar smoothing of the2◦ f -mode fits produces an outflow. (lower row ) Flows obtained at a depth of 4 Mm from 1-D
OLA inversion of 16◦ data, ARRDI results smoothed to 16◦ resolution, and smoothed p1 mode fitsobtained with 4◦ data. Similar inflow structures are present in the 1-D and 3-D inversion results.
The 4◦ fits exhibit similar outflows to the 2◦ fits. Despite the similarities in small-scale structure
between the 2◦ and 4◦ fits and the ARRDI results, the large-scale means seem to be established by
the 16◦ data.

139
the non-smoothed ARRDI results agreed relatively well with the 2◦ f-mode fits at this depth, thelarge-scale means do not. The 2◦ fits display a strong outflow at this scale whereas the inversion
results clearly show an inflow. This is no doubt due in part to the lack of strong supergranularinflows into the active region and the absence of 2◦ data at some locations in the active complex.
The smoothed 3-D inversion results agree well with the 1-D 16◦ OLA inversions at this depth, andreproduces the general sense of the inflow into the active region.
Similar behavior in the results is apparent in Figure 8.6b where we examine the solutions at
4 Mm and compare with the fits to the 4◦ p1 data (which have not been inverted). Once more thefitted data possesses a strong outflow feature when smoothed to the 16◦ scale that is not sharedby the 3-D inversion results. The inflow into the active region found in the ARRDI

solution issomewhat weaker than in the 16◦ results, but is present, as is the northeasterly outflow emanatingfrom the trailing edge of the active complex.
The large-scale means resulting from the ARRDI inversion are thus determined primarily by
the 16◦ measurements, while the small-scale structure generally reflects that of the higher resolution
measurements. From examination of Figure 5.4, we see that the typical errors for the 2◦ and4◦ measurements in the upper 4 Mm of the convection zone are on the order of 60 m s−1 and 40
m s−1 respectively. A mean flow with an amplitude of 20 m s−1 such as that of Figure 8.6 is thus
well within the errors of the high resolution measurements. The errors on the 16◦ measurements
are substantially smaller, however, with a typicaly value of 10 m s−1 at these depths. A reversalof the situation would not work. If the large-scale flows of the solution were set by the 2◦ or
4◦ measurements, their 20 m s −1 outflow would not agree with the 16◦ measurements to within
their errors. The flowmaps resulting from ARRDI thus represent a multi-scale solution to theinversion problem (Eq. 7.1) that is self-consistent to within the errors of the measurements.

140

8.2

Sunspot Outflows
Comparisons between 3-D inversions with ARRDI and previous ring-inversion procedures

suggest that the ARRDI algorithms yield sensible results. As such we devote the remainder of thischapter to an exploration of ARRDI’s potential through its application to flows around sunspots.
We focus on sunspots viewed with MDI during January 2002. As with the 2001 data used forcomparisons in the previous section, the Sun was highly active during this period, and a numberof active regions and sunspots were visible on the disk. One particularly circular sunspot in active
region NOAA AR 9787 was observed during this time, and its axisymmetry makes it well suited forthe study of average outflows. Outflows around the sunspot in NOAA AR 9787 have been mapped
at depth using time-distance by Gizon et al. (2009), and so analyses of this dataset enables the
direct comparison between inversion results obtained with ARRDI to those obtained through timedistance.
Our ring-analysis measurements were obtained by tiling the solar disk each day with a mosaic
of 16◦ analysis regions overlapping by 7.5◦ as described in §5.2. Each of these 16◦ analysis regions
was then further splintered into a similarly structured mosaic of 2◦ and 4◦ regions overlapping by
0.9375◦ and 1.875◦ respectively. The tiling resolution of these measurements is thus somewhat
higher than those used for our initial assessment of ARRDI. We have adjusted the horizontal
resolution of our inversion grid accordingly, adopting a horizontal resolution of 0.9375◦ at thesurface. In order to ensure that we gain the full resolving power of the 2◦ measurements, we maintainthis high-resolution down to a depth of 3 Mm before degrading to the resolution appropriate forthe 4◦ measurements. Details of the inversion grid are provided in Table 8.2.

8.2.1

Ubiquitous Sunspot Outflows at Depth

A sampling of the flows and magnetism present during the period spanning 11 Jan 2002
through 25 Jan 2002 (Carrington Rotation 1985) is shown through synoptic mapping of the flows
at a depth of 5 Mm in Figure 8.7. This map has been generated by combining full-disk stitchedDepth = 5.0 Mm
40

Latitude

20

0

-20

-40
100

150

200
Longitude

250

300

100 m/s

Figure 8.7: Synoptic map of Carrington rotation 1985 spanning the period 11 Jan 2002 through 25 Jan 2002 at a depth of 5 Mm as
measured using ARRDI. Daily flowmaps obtained using full-disk stitched inversions were averaged in time over the this period. Each flow
vector represents a seven-day average of the flow at that latitude and Carrington longitude. Sunspots, visible as dark green and red spots
on the map, are home to strong 100 m s−1 outflows. The torsional oscillations are visible at this depth as strong longitudinal flows in theactivity belts.

141

142

-8

-8

-10

-10

-12

-12

-14

e

2.0 Mm

0

-5

-10

-14

-16
200

202

c

Latitudeb-6

Depth (Mm)

0.5 Mm204
206
Longitude

5.0 Mm-16
208
300 m/s
-6

-8

-8

-10

-10

-12

-12

-14

-14

200

202

204
206
Longitude

202

d

-6

-16

200

-16
208
100 m/s

204
206
Longitude-15

208
300 m/s

0

100
200
300
Radial Velocity (m s-1)

f

13.0 Mm

400

0
220 m/s
Depth (Mm)

Latitude

a
-6

5

10
-50 m/s

200

202

204
206
Longitude

208
40 m/s

15
0

10

20
30
40
Radius (Mm)

Figure 8.8: The depth variation of a sunspot outflow. (a-d) Horizontal flow beneath NOAA AR
9783 on 17 Jan 2002 at four representative depths as determined with ARRDI. Concentric circles
indicate radii of 18, 30, and 42 Mm relative to the sunspot center. (e) Depth variation of theradial outflow (shown positive), averaged in angle around the sunspot center, is indicated at radiiof 18 Mm (black), 30 Mm (blue), and 42 Mm (red). Solid horizontal lines correspond to the depths
indicated in a-d. (f ) 2-D view of the sunspot outflow with radius and depth. 150 m s−1 outflows
are visible at the surface and diminish in strength around 3 Mm. Outflows increase to 350 m
s−1 around 5 Mm before decaying steadily with depth. The sunspot outflows decay steadily with
distance at all depths. A near-surface inflow is visible at the perimeter of the surface moat flow.
Flow fields have been interpolated by a factor of two in all cases for ease of viewing.

143
TABLE 8.2
Inversion Grid For Sunspot Inversions
Target Depth (Mm) Horizontal Resolution (degrees)
0.0
0.93750.5
0.93751.0
0.93751.5
0.93752.0
0.93753.0
0.9375
4.0
1.875
5.0
1.875
6.0
1.875
7.0
1.875
9.0
1.875
11.0
1.875
13.0
1.875
17.0
3.75
19.0
3.75
22.0
3.75
25.0
3.75

Table 8.2: Target depths and horizontal resolution of the inversion grid for the 3-D inversions ofsunspots from CR 1985.

inversions for each day and averaging the flowmap from multiple days at each Carrington longitude.
The Carrington grid rotates 15◦ per day, and when viewed in this manner, each vector represents
a 7-day average of the flow. Flows with a horizontal scale of 4◦ are most readily visible in this flow
map, and a distinctive feature of sunspots at this depth (visible as dark concentrations of red or
green in the magnetogram), are their 100 m s−1 outflows.
A detailed view of the subsurface outflow beneath a typical sunspot (in NOAA AR 9783 )from this rotation is illustrated in Figure 8.8. Flow maps at four representative depths ranging
from the near-surface down to 13 Mm are shown in Figures 8.8a-d. Reference circles are plottedat radii of 18 Mm, 30 Mm, and 42 Mm with respect to the sunspot center, and the radial outflows(averaged in angle about the sunspot) are plotted in Figure 8.8e. Near the center of the sunspot,
100-200 m s−1 outflows are present at all depths, diminishing in strength with distance from thesunspot center. In these near-surface layers, the most distinctive convective patterns are roughly3◦ (36 Mm) in horizontal extent, and so it is likely that we are sampling the largest supergranular144
cells of convection. These cells extend across the outer reference circle at 42 Mm, contributing to
the net inflow observed in the upper 4 Mm at that radius (Fig. 8.8e, red line).
An interesting feature of the NOAA AR 9783 sunspot is the apparent presence of two components to the subsurface outflow. This is most visible in Figure 8.8f where the average outflow
is shown with depth and radius. The outflow is characterized by a near-surface component which
decays with radius much more quickly than a deeper component which peaks around 5 Mm. The
transition between these two flow components occurs around 2 Mm. We see from Figure 8.8b thatflows at this depth are similar in morphology to those near the surface, but that they are weakeroverall than the near-surface flows. The presence of a deep outflow, more expansive in radius than
the surface moat flow, seems to be a robust feature of the sunspots in this time interval. Further
examples of this behavior are provided in Figure 8.9. The deep outflow is not always stronger than
the surface moat flow at the 18 Mm radius in these examples, but it is the dominant component of
the flow at a distance of 30 Mm in most instances.
The weakening of the surface outflow occurs at a depth where reliable fitting of the 2◦ measurements becomes difficult. It is thus somewhat unclear if this diminishment is real or is in fact
an artifact of the effective transition in the horizontal resolution of the measurements from 2◦ to
4◦ at this depth. We also note that the transition region between the two components lies close
to the depth where the horizontal resolution of our inversion is degraded by a factor of 2. However, we have run several inversions varying the location of this transition and find similar results.
Thus while we cannot rule out the existence of an averaging effect resulting from the measurement
procedure, it seems unlikely that this two component structure is an artifact of the inversion grid
resolution.

8.2.2

Temporal Evolution of Sunspot Outflows

These subsurface outflows exhibit time variability as their associated sunspot evolves as well.
The time-evolution of the sunspot in NOAA AR 9783 is illustrated in Figure 8.10. Flow vectors at
a depth of 7 Mm overlay a magnetogram in the upper row of Figure 8.10 and the radial outflow is

145
1.5 Mm
-4

20

32
-6

Latitude

18

30
-8

16

28
-10

14

26
-12

12

24

150

152

154

156

158

160

128

130
132
134
136
Longitude (Degrees)

266

268

270

272

274

266

268

270

272

274

300

500 m/s

7.0 Mm
-4

20

32
-6

Latitude

18

30
-8

16

28
-10

14

26
-12

12

24

150

152

154

156

158

160

128

130
132
134
136
Longitude (Degrees)

200 m/s
0

0

5

5

5

10

10

10

15
-100 -50 0 50 100 150 200
Radial Outflow (m s-1)

15-200 -100 0
100 200
Radial Outflow (m s-1)Depth (Mm)

0

15
-50

0 50 100 150 200 250
Radial Outflow (m s-1)

Figure 8.9: Sampling sunspots outflows in CR 1985 for three sunspots as they crossed the central
meridian of the disk. Outflows at radii of 18, 30, and 42 Mm are indicated as in Fig. 8.8. Zero
outflow is indicated by the dashed line. (left column) Flows beneath NOAA AR 9788 on 22 Jan 2002
possess surface and deep outflows similar in nature to that of NOAA AR 9783 in Fig. 8.8, with adiminishment near the 3 Mm depth. While many of the spots we have examined possess these dual
outflows, not all do so. (central column) The sunspot NOAA AR 9787 (23 Jan 2002) exhibits a lesssimple outflow structure, and (right column) that of NOAA AR 9779 (13 Jan 2002) decays more orless steadily with depth.

146
plotted for each day in the lower row. This particular sunspot exhibits a fairly persistent outflow
as it rotates across the disk. The near-surface and deep components of this flow wax and wanewith time, as does the near-surface inflow which is always visible at the 42 Mm radius. With thepossible exception of the first day in the series, the deep outflow persists at greater distances than
the near-surface outflow for all days.
We have also used ARRDI to follow the evolving flow of a much shorter lived sunspot, thisone in NOAA AR 9788 , as it emerges near the western limb and disrupts just before rotating off theeastern limb. A seven-day time series of the emergence and disruption is provided in Figure 8.11
where vector maps illustrate the flows at 5 Mm in depth. Some hint of an outflow is present on 19
Jan just as the spot emerges, with a much more prominent outflow developing on the second day(20 Jan). The development of the outflow appears to emanate from the surface and the interior of
the spot before developing into the more familar two-component structure on the third day. This
two-component structure persists for a total of three days before the sunspot begins to decay, withthe decay of the surface moat flow apparently preceeding that of the deeper outflow (24 Jan).

8.2.3

Sunspot Outflows as Revealed through ARRDI and Time-Distance Inversions

The sunspot in NOAA AR 9787 was the focus of time-distance investigations reported on in
Gizon et al. (2009). Outflows measured with ARRDI beneath NOAA AR 9787 will thus allow
a comparison between these two techniques. We have averaged the radial outflows beneath thesunspot in NOAA AR 9787 over the six days spanning 20 Jan 2002 through 25 Jan 2002. A similar
temporal averaging of the outflows beneath this spot was performed in Gizon et al. (2009). The
time-distance analyses were performed using the f -mode ridge, and the p1 through p4 ridges.
Center-to-quadrant travel times were then measured and inverted using an OLA scheme to recover
the east-west and north-south flows. These were in turn averaged in azimuth about the sunspot toyield the radial outflow. Further details can be found in Gizon et al. (2009).
The results of the time-distance inversions at four depths are shown in Figure 8.12a, alongwith those measured using moving magnetic features (MMF, dashed line). Results obtained from

14 Jan 02

16 Jan 02

17 Jan 02

19 Jan 02

20 Jan 02

21 Jan 02

Latitude

-5

-10

-15
200

205

210

200

205

210

200

Depth (Mm)

150 m/s

205

210

200

205

210

200

205

210

200

210

Longitude (Degrees)

0

0

0

0

0

0

5

5

5

5

5

5

10

10

10

10

10

10

15

15
-50

205

100

250

15
-50

100

250

15
-50

100

250

15
-50

100

250

15
-50

100

250

-50

100

250

-1

Radial Outflow (m s )Figure 8.10: The slow evolution of a persistent sunspot. (upper row ) Horizontal flow beneath NOAA AR 9783 at a depth of 7 Mm as
determined using ARRDI. Concentric circles indicate radii of 18, 30, and 42 Mm relative to the sunspot center. (lower rower ) Radial
outflow averaged in angle about the sunspot center for radii of 18 Mm (black), 30 Mm (blue), and 42 Mm (red). A zero outflow is indicatedby the vertical dashed line. This particular sunspot exhibits a surface moat flow and a deeper outflow during all six days over which it
was observed. The deep outflow at 6 Mm attains an amplitude of 250 m s−1 before declining steadily with time. The surface flow waxes
and wanes in strength over the course of six days, diminishing in strength during the last day in the series.
147

19 Jan 02

20 Jan 02

21 Jan 02

22 Jan 02

23 Jan 02

24 Jan 02

25 Jan 02

Latitude

20

15

10
150

155

160 150155

160 150155

160 150

Depth (Mm)

175 m/s

155

160 150155

160 150155

160 150

160

Longitude (Degrees)

0

0

0

0

0

0

0

5

5

5

5

5

5

5

10

10

10

10

10

10

10

15

15
-50

155

100

250

15
-50

100

250

15
-50

100

250

15
-50

100

250

15
-50

100

250

15
-50

100

250

-50

100

250

-1

Radial Outflow (m s )Figure 8.11: Emergence and disruption of a sunspot. (upper row ) Flows beneath NOAA AR 9788 at a depth of 5 Mm as determined with
ARRDI. As with Fig. 8.10, concentric circles are plotted at 18, 30, and 42 Mm from the sunspot center. (lower rower ) Average radial
outflows at each of the radii, with 18 Mm indicated in black, 30 Mm in blue, and 42 Mm in red. No significant outflows are apparent
initially. Strong outflows develop on 20 Jan that decrease in strength with depth. A two-component structure is attained on 21 Jan,and the sunspot persists for three days before decaying. Outflow profiles during the final three days suggest that the surface moat flow
dissipates prior to the disappearance of the deeper outflow at 5 Mm.
148

149

Figure 8.12: Comparison of outflows beneath NOAA AR 9787 as revealed through time-distance
helioseismology and ring analysis. (a) Subsurface outflows measured with time-distance at various
depths (indicated) in the upper 5 Mm of the convection zone. (Adapted from Gizon et al. 2009)(b) Outflows as measured with 3-D inversion of ring-analysis measurements using ARRDI. Coloring
for the different depths is the same as in (a). Time-distance results indicate a steady increase in
the strength of the outflow down to a depth of about 5 Mm, whereas results from ring-analysis
indicate a diminishment of the outflow at 2.6 Mm before increasing again at 4.5 Mm.

ring analysis using ARRDI are displayed in Figure 8.12b. We see good agreement between thesurface moat outflow (1 Mm, red line) as measured by both techniques. The moat flow in thesurface layers peaks around a radius of 20 Mm and extends to about 40 Mm for both inversions.
The inflows beyond 40 Mm are similar in strength at this depth, but the peak values of the moat flow
are somewhat higher for the ring-analysis measurements than those achieved using time-distance.
Deep below the surface, there are considerable differences in the variation of the outflow withdepth between these two approaches. At a depth of 2 Mm (dark blue), the time-distance results
suggest that the moat outflow increases in strength relative to the surface, and the inflow at 40
Mm in radius disappears. Here, the ARRDI results differ and indicate that the outflow weakens instrength at this depth, while the inflow remains visible. An outflow persists to 4.5 Mm in depth
in both sets of results, though the outflow as measured with ARRDI is still weak relative to its
surface values. The time-distance results seem to indicate that the strength of the outflow increases
with depth and do not hint at the presence of a two component outflow such as in the ARRDI
inversions of NOAA AR 9787 and other sunspots from CR 1985. The ARRDI results thus agree

150with those obtained from time-distance in that outflows are realized at all depths, but there are
distinct differences in the depth dependence and amplitude of these outflows.

8.3

Further Reflections on 3-D Inversions with ARRDI
ARRDI holds distinct promise as a versatile tool for probing flows in the upper-convection

zone. Flow maps obtained with ARRDI exhibit good agreement with those from previous ring
analysis studies at depths where it is sensible to compare them. Particularly encouraging is theability of ARRDI to compensate for missing data from one tile-size (usually the 2◦ tiles) throughthe use of lower resolution measurements when available. This is useful for establishing large-scalemean flows in the near surface layers where noisy high-resolution measurements and missing data(such as in regions of magnetic activity) can otherwise lead to artifacts in the inferred large-scale
flows. As a “first light” application of ARRDI , we have examined the subsurface outflows aroundsunspots and found that many spots from CR 1985 possess subsurface outflows with a distinctive
two-component structure.
The discovery of this structure is surprising in part for its disagreements with the timedistance results of Gizon et al. (2009), which suggest that the outflow should monotonically increase
with depth. However, we are encouraged overall, rather than discouraged, by the comparison
between ARRDI and the time-distance results. Both approaches show that some form of outflow
persists to the 5 Mm depth, and we are not troubled by the overall amplitude differences between the
two approaches (ring-analyses show weaker flows overall). Our RLS approach necessarily trades
off between a good fit to the data and a smooth solution. By smoothing more or less (i.e. by
varying the regularization parameters), the amplitude of our recovered solutions tends to vary as
well. For instance, we have found qualitatively similar outflow structures, but with a velocity
amplitudes more in accord with the time-distance results, by decreasing the vertical regularizationfrom λv = 0.1 to λv = 0.01.
A detailed comparison between the time-distance and ARRDI inversions is difficult without
knowledge of the averaging kernels realized in the time-distance inversions. These kernels were151

Figure 8.13: Sunspot outflow as realized through 3-D numerical simulations. Red tones indicate
outflow from the sunspot, and blue tones indicate inflow. Black lines denote azimuthally averaged
flux surfaces of the sunspot’s magnetic field. The flow is seen to possess three components, theEvershed flow (red), a superficial moat flow (yellow-far upper right), and a deeper more extended
outflow (orange) underlying the penumbra. The deeper outflow is reminiscent of those measured
using ARRDI (see Fig. 8.8). Adapted from Rempel (2011).

not reported in Gizon et al. (2009), but averaging kernels obtained through time-distance with
OLA are undoubtedly different from those obtained through ring analysis with RLS. Differencesin the vertical widths alone could easily lead to differing outflow profiles. For instance, sufficiently
broad (in the vertical) averaging kernels around the 2-3 Mm depth for the time-distance inversions
could wipe out the two-component structure we have observed. Moreover, the presence of magnetic
fields can confound the interpretation of results obtained using either ring-analysis or time-distance
helioseismology through a tendency to absorb acoustic power and introduce phase shifts into the
wavefield (e.g. Braun et al. 1987; Cally et al. 2003, Lindsey & Braun 2005a). The sensitivity kernels
employed in ring analysis and time-distance methods do not account for these effects. Bearing these
considerations in mind, we suggest that ARRDI yields results that are largely consistent with thosededuced using time-distance helioseismology.
We conclude by noting that the two-component structure we observe in outflows beneathsunspots is remarkably reminiscent of sunspot outflows found through 3-D numerical MHD simu-

152
lations of a bundle of magnetic flux introduced into the near-surface convection and atmosphereof the Sun (Rempel 2011). These simulated outflows, which arise naturally from the blockage ofheat transport by the magnetic field, are detailed in Figure 8.13 and appear to possess a multicomponent flow structure of their own. The surface moat flow (light yellow region near 30 Mm)
decays with depth before joining to a second outflow around 4 Mm. The depth variation of thesesimulated outflows is similar to that observed with ARRDI, but the deep outflow seen by Rempel(2011) seems to decay more quickly in radius from the spot than does the surface moat outflow.
This contrasts with our observations of a deep outflow that extends further in radius than thesurface inflow. We caution that we have only analyzed a handful of sunspots using ARRDI andthat a much larger sample is needed to establish whether these two-component outflows truly area robust feature of sunspots.

Chapter 9

Reflections and Future Directions

9.1

Convective Cores and Primordial Magnetic Fields in A stars
Our explorations of dynamos in A-type stars have revealed that the deep interior cores of thesestars can attain two remarkably different dynamo states, characterized by equipartition magnetic
field strengths on the one hand, and super-equipartition magnetic fields on the other. We havefound that one route to transitioning from the first to the second state involves the presence of alarge-scale ordered magnetic field in the radiative zone with some linkage to the core. The possible
existence of a primordial magnetic field in the radiative zones of the Ap stars was our primary
motivation for the inclusion of such a magnetic structure. Our initial field was thus twisted andlargely toroidal in nature, but we found that the transition to a super-equipartition state did not
require such a specific geometry. It appears instead that the transition between equipartition andsuper-equipartition behavior is a robust response to the presence of an organized poloidal fieldstructure linking the radiative zone to the core. The existence of such a strong dynamo state in
the cores of these stars raises a number of interesting questions regarding the limits of magnetic
field strengths achievable, as well as the effects of such strong dynamo states on other stages of
evolution of these stars.

9.1.1

Possible Magnetic Buoyancy Instabilities

The discovery of such strong magnetic fields in our simulations raises the question of whether
magnetically buoyant fields might be achieved in the cores of the A stars. The super-equipartition

154
dynamo realized following the imposition of our fossil magnetic field has displayed a tendency forbuilding strong toroidal bands of magnetic field which are roughly 300 kG in strength in someregions. We have not observed the buoyant rise of these magnetic fields, but there are diffusive
limitations within our model which bear on the formation of strong magnetic fields in the radiativezone and their ability to rise.
As our models are large-eddy simulations, the magnetic diffusivities that we have adopted arequite high with respect to what might be expected in the deep interior of an A star. In particular,
the diffusion times across the radiative envelope of the star in our simulations are several thousanddays. This time is several orders of magnitude shorter than the realistic values. What might be theeffect of substantially lower diffusivities in our system? The ramifications for the convective core
dynamo are unclear.
Decreasing the magnetic diffusivity can certainly limit the diffusive destruction of magnetic
energy in the core and permit the development of much finer structuring in the magnetic field
topology than what is realized in our current simulations. However, the strength of the magnetic
fields in our super-equipartition state appears to be set by a magneto-geostrophic balance within
the core, in which the pressure, Coriolis, and Lorentz forces largely balance one another. Can
magnetic fields strong enough to break this balance arise in the core? This may not be necessary
to attain buoyant magnetic fields in the radiative zone.
Magnetic field strengths in the radiative zone are largely determined through a balance between the deposition of magnetic fields by overshooting motions and the slow decay of these fields
through diffusive processes. The potential for our super-equipartition dynamo to generate magnetically buoyant fields in regimes with lower magnetic diffusivities may thus be interesting as a further
avenue of study. Recent modifications to ASH have allowed us to begin simulating convection instars while employing a slope-limited diffusion scheme. Such schemes permit much lower diffusivities than those achieved currently using ASH, but complicate things somewhat in that the diffusion
is spatially variable, determined by the local gradients of the flow and magnetic field. Assessing a
representative magnetic Prandlt number for such a system can thus be difficult. However, it would

155
be fascinating to observe magnetic buoyancy in the radiative envelopes of these stars, and we planto pursue further simulations using lower diffusivities and more extensive radiative envelopes (tocapture any buoyant rise of magnetic field) in the near future.

9.1.2

Convective Overshooting With and Without Primordial Fields

An interesting aspect of core convection is that the convective motions occur in the presenceof nuclear burning. Convective overshooting into the overlying radiative envelope may then mix the
products of nuclear burning from the core into the lower radiative zone, returning fresh hydrogen to
the core in the process. Such a replenishment of the nuclear fuel supply holds potential for extending
the main-sequence lifetimes of massive stars. The efficiency of the large-scale motions accessible
to ASH in assisting such refueling might be assessed by adding a second chemical species to the
ASH equations. This would constitute a straightforward modification to the ASH code, thoughlarge disparities between the diffusion rates of H and He may severely limit the time stepping.
Nevertheless, an interesting extension of the mixing problem would be to compare the mixing
efficiencies obtained by convective motions from an equipartition dynamo to those realized in thepresence of the super-equipartition magnetic fields. If the presence of a strong primordial magnetic
field (likely the distinguishing characteristic between the Ap and normal A stars) does lead to alarge-scale dynamo in the convective core, then differences in the mixing effiencies might imply
disparities in the main-sequence lifetimes of the Ap stars and the normal A stars.
In principle, such differences could be probed observationally by examining the relative agesof Ap and normal A stars in open clusters. Pinning down the ages of Ap stars has been difficult until
recently since many of the known Ap stars are field stars, and the known cluster stars are typicallyvery faint. However, serious progress has been made in recent years to identify and catalogue
cluster Ap stars (Bagnulo et al. 2006; Landstreet et al. 2007; 2008). By constraining the ages ofcluster Ap stars, these surveys have been able to demonstrate (for example) that Ap stars exhibit
less surface magnetic flux as they near the end of their main-sequence lifetimes (Landstreet et al.
2007;2008). However, no comparison between the apparent ages of the normal A and the Ap stars

156
within a cluster has yet been conducted.
An alternative approach might be to assess the relative number of Ap-descended giant stars
within a cluster. Distinguishing between the two can be difficult, but some red giants have been
observed to possess chromospheric activity which is more intense than what might be expected
for their mass and rotation rate (e.g., Strassmeier et al. 1999; Aurière et al. 2008), and it has
been speculated that these inferred magnetic fields are primordial, rather than dynamo in origin.
Abundance peculiarities, though of a different sort than in the main-sequence Ap stars, may provide
a more reliable means for distinguishing between giants of peculiar and normal descent. As discussed
in Charbonnel & Zhan (2007) the presence of strong magnetic fields in the 104 to 105 kG range
might suppress thermohaline mixing in the descendants of Ap stars. This mixing mechanism, driven
by a µ-gradient instability, is thought to link the deep hydrogen burning shell of red giants to theouter convective envelopes. Its absence in the giant descendants of the Ap stars (owing to theirpresumably strong magnetic fields) may lead to a perceived overabundance of C and Li on thesurfaces of these stars.

9.1.3

Primordial Magnetism and Dynamos in Pre-Main Sequence Stars

One of the untouched questions of primordial magnetic fields inside massive stars concerns
the effects of a fully convective Hayashi phase that these stars are thought to undergo on theirjourney to the main sequence. Does a dynamo operate during the fully convective phase, and if so,
how does a pre-existing magnetic field alter its characteristics? Some sense of this is gained throughthe simulations presented in Chapter 2-4, but a fully convective sphere lacks the radiative envelopepresent in a main-sequence A star. The importance of the radiative envelope for attaining strong
and global-scale magnetic fields in our simulations is difficult to discern, but it undoubtedly plays
some role in the magnetic “memory” of the system. Overshooting motions deposit broad swaths ofmagnetic field into the lower radiative zone before diving back towards the deep core, revisiting theirpast works with each new overturning motion. The effects that a large-scale primordial magnetic
field in the 10 kG range might have on a fully convective dynamo with no such reservoir for the

157
magnetic field are thus unclear.
The discovery that the Ae stars, the presumed pre-main sequence (PMS) progenitors of the
Ap stars, exhibit a preference for slow rotation even before their arrival on the main sequence(e.g., Alecian et al. 2008; Wade et al. 2009) hints at the importance of the fully convective
phase in understanding the angular momentum loss of these stars. The nature of the differential
rotation and the convectively driven dynamo in this phase will undoubtedly impact the pre-main
sequence spindown through their interaction with a stellar wind. Moreover, the magnetic fieldremaining after the convective envelope has receded will ultimately serve as the seed magnetic field
for the convective-core dynamo. Does the primordial field maintain its twisted, toroidal shapethroughout the fully convective phase? Is it amplified or augmented in some fashion, or is the roleof turbulent convection in this regime instead to expedite the decay of the field? The destruction
or amplification of a primordial field in this stage may bear on why only 10% of the A stars aremagnetic (and peculiar).
Primordial magnetic fields may play a role in the early evolution of less massive stars as well,notably the the T Tauri stars. The majority of these stars are thought to be fully convective andexhibit surface fields of a few kG. Such strong magnetic fields are almost certainly generated bya dynamo operating in the midst of this convection as strong flaring observed on many of thesestars suggests that magnetic fields are undergoing constant reordering at the stellar surface (e.g.,
Gregory et al. 2010). And yet, the observations of magnetic fields in the fully convective T Tauristars suggests that they are predominantly large-scale and dipolar in nature. How are such strong
coherent magnetic fields generated in the presence of vigorous convection?
Any attempt at modeling these stars with ASH would necessarily have to ignore the circumstellar accretion disk, and a possibly important contributor to the angular momentum balance
within the star in the process. Nevertheless, it is interesting to ask what effects a primordial magnetic field may have on dynamos in other stars in the early stages of their evolution. In a broader
sense, what sort of dynamos might operate in pre-main sequence stars? Are these fully convective
objects characterized by the sweeping global-scale motions found in simulations of the A stars, or158
are they more like the convective M-dwarfs, where vigorous small-scale convection occurs in theouter layers, with much weaker but larger scale motions comprising the deep interior (e.g., Browning et al. 2008). Models of pre-main sequence stars, both of the massive and less massive varieties,
may be interesting to pursue with ASH. ASH provides the capability of pursuing dynamo calculations within the full spherical geometry and may be a useful tool for assessing how the global-scale
magnetic fields observed on some PMS stars are generated. However, we acknowledge that certain
complicating factors may arise when simulating convection in the early stages of stellar evolution.
The presence of an accretion disk or the slow gravitational contraction of the star, for example,
would certainly require some careful treatment if we are to assess the flavor of the dynamos in thesestars using ASH.

9.2

Possible Further Improvements to 3-D Inversions with ARRDI
ARRDI was designed with implementation into the HMI data analysis pipeline in mind,

and adapting the ARRDI algorithms to HMI data is straightforward in many respects. However,
exploiting the full potential of ARRDI requires that similar assessments of the measurement errors
and measurement success rates are performed using HMI data as were done with MDI. These
assessments are crucial for deciding upon a sensible inversion grid. HMI has a fourfold increase in
resolution over that currently achievable using MDI and will return Dopplergrams with a cadenceof 45 s vs the 60 s used with MDI. This increased spatial and temporal resolution will undoubtedly
affect the noise characteristics of the power spectra and the errors realized when performing ring
analyses.
The increased horizontal resolution attainable with HMI may also allow ring analyses to be
performed on smaller tile sizes such as 1◦ or possibly even 0.5◦ in horizontal extent. The radial
orders accessible to ring-analysis fits of these smaller tiles (and hence the depths to which theinversions might reach) are currently unknown. Based on the success rates associated with fitting
the 2◦ data, it seems safe to assume that at least a few f -modes might be fit, thus allowing us toprobe the surface layers with a higher resolution. It remains unclear if HMI will permit the fitting

159
of modes with higher radial order (and thus more depth penetration) than is currently possible
when using MDI data. However, efforts are now underway to develop a fitting procedure thatmeasures the Doppler shifts of multiple ridges simultaneously (D. Haber, private communication).
In particular, fitting in this manner may enable the measurement of modes beyond the p1 ridge forthe 2◦ tiles. Finer horizontal resolution, both near the surface and at depth, thus seems possible in
the near future. There are a few adjustments that might be made to the ARRDI algorithms andto the sensitivity kernels in order to fully exploit this capability.

9.2.1

Obtaining Optimal Horizontal and Vertical Resolution

It is interesting to ask if the horizontal resolution of the ARRDI inversion grids (using 2◦ , 4◦ ,and 16◦ measurements) might be improved. The effective horizontal resolution of the inversion grid,
which is manifest in the narrowness of the averaging kernels obtained, is intimately related to themeasurement tile size and the nature of the overlap between the tiles. What is the optimal tilingscheme? We have seen that for 4◦ measurements, ARRDI tends to produce averaging kernels with awidth of about 2◦ and with small negative sidelobes. These sidelobes stem from the manner in which
the inversion algorithm combines measurements. The solution at some target location is typicallygenerated by subtracting measurements made at neighboring points from measurements centeredon that target point. The resulting averaging kernel is constructed by adding and subtracting
the associated sensitivity kernels in a similar fashion, yielding a central peak with small negativesidelobes in the horizontal.
A narrower central peak can in principle be produced by increasing the tiling density of themeasurements. However, as the overlap of neighboring measurements is increased, the wavefields
they sample are largely the same. These RLS inversions operate on the implicit assumption that themeasurements and their errors are independent of one another. This overlap can lead to correlationsin the errors which should be taken into account for the more densely packed measurements, andpossibly for the current measurements which overlap by half of the analysis region size. In a similarvein, one might assess the optimal combination of tiling schemes and measurement resolutions.

160
Currently the 16◦ measurements are always centered on a 4◦ measurement, which is in turn centeredon a 2◦ measurement. Are there advantages to be gained in the averaging kernels by staggering
the different measurement grids with respect to one another? A study of these tiling scheme issues
could be useful for optimizing both the time required to generate the measurements and the timerequired to perform the inversion.
Finally, we have yet to perform a systematic study of the optimal vertical gridding used in
ARRDI. Such a task is accomplished much more easily in the horizontal direction, where eachring measurement has one of three distinct horizontal structures (corresponding to the 2◦ , 4◦ and16◦ tiling sizes), and falls upon a regular grid. However, the situation is complicated in the verticaldirection, where a continuum of eigenfunctions (associated with the different fitted modes) is available to construct a vertical averaging kernel. The details of how the vertical averaging function isconstructed thus depend on the modes available and the vertical grid. We have performed some
general testing of these issues and find that finer resolution than what was adopted in Chapter
8 yields no obvious benefits. However, we acknowledge that further exploration of the interplay
between these two effects may enable us to establish a minimally acceptable vertical resolution,
thereby increasing the efficiency of the 3-D inversion.

9.2.2

Further Possibilities for Improvements to the Sensitivity Kernels

Some work still lies ahead for the sensitivity kernels as well. Most notably, the effect of line-ofsight changes on measurements made near the limb needs to be accounted for in the kernels. As thecenters of the analysis regions are extended toward the limb, measurements become increasingly
sensitive to the horizontal components of the solar wavefield. In their current formulation, thesensitivity kernels assume that measurements are taken at the disk center. As a result, the kernel
calculation is concerned with the effects of advection on only the vertical component of the wavefield.
The inclusion of these effects will be vital, for instance, for the proper interpretation of outflows
around sunspots close to the limb.
We note that it may be interesting to develop structure inversion kernels, such as sound161
speed kernels. Conceptually there is no difficulty in using the existing kernel machinery to build
sensitivity kernels for sound-speed perturbations. The ring-analysis measurement technique looks
for anisotropic shifts in the mode frequencies and would need to be modified to account for theisotropic shifts induced by sound speed variations. The kernels would then need to be readjusted
to reflect the new measurement procedure. Knowing both the temperature profile and outflow
profiles beneath sunspots will help to further characterize their morphology. However, the effectsof sunspot magnetic fields on the properties of acoustic waves can complicate the interpretation of
the results.
Another limitation of the sensitivity kernels is the lack of any treatment of magnetism intheir calculation. What effect the presence of magnetism may have on ring-analysis measurements
is unclear at this time. Perhaps the most noticeable effect of a magnetic field is to absorb acoustic
power. This absorption is largely frequency dependent, steadily increasing up to 3 mHz (e.g.,
Braun 1995). A power spectrum constructed for a region adjacent to a sunspot will exhibit adiminishment in power for waves travelling away from the sunspot relative to those travelling
toward it. The sensitivity kernels assume that the power spectrum is isotropic in the absence offlows, and other anisotropies such as those introduced by a magnetic field remain unaccountedfor. It is also difficult to say how a frequency dependent diminishment in power might alter themeasurement procedure. It could modify the slope of the power spectrum in frequency, which issomething the measurement procedure is sensitive to. We thus acknowledge that a magnetic field
can certainly affect flow measurements, and that these effects are not accounted for in the kernels
or by the measurement procedure.

9.3

Future Prospects for Ring Analysis with 3-D Inversions
Though some improvements are desirable, the ARRDI algorithms are now mature enoughto warrant their application to a wide variety of phenomena in the upper convection zone. We
expect ARRDI to find wide application in the analysis of HMI data, but note that the ARRDI
algorithms can now be applied systematically to the vast existing archives of MDI Doppler data.

162
These archives span both quiet and active Sun periods and constitute an invaluable resource forlong-term and statistical studies of various elements of dynamics in the upper convection zone thatwill complement ongoing observations with HMI.

9.3.1

Detailed Probing of Sunspot Outflows

The studies presented in Chapter 8 suggest that many sunspots possess a two-component
structure to their outflows. This discovery is surprising in that it disagrees with earlier timedistance studies which suggest that sunspot outflows increase in strength monotonically down toa depth of 5 Mm (Gizon et al. 2009). We have only analyzed a handful of sunspots to date,and so the robustness of this two-component structure remains unclear. The ability to probe thesubsurface flows in such detail allows us to ask other questions concerning the sunspot morphology
as well. Are subsurface flows similar for all sunspots, or do they differ based on the magnetic fieldstrength and the age of the spot? We have already seen some hints that the surface outflow decays
before the deeper outflow (Fig. 8.11). Moreover, the deep flows around the sunspot in NOAA
AR 9783 (Fig. 8.8) exhibit a strong vorticity at a depth of 13 Mm that is not seen in the othersunspots. Is such subsurface vorticity a common feature in many sunspots, and how does the depth
dependence of the subsurface vorticity correlate with the outflows we observe? These questions are
difficult to answer based on an examination of only a handful of sunspots. However, with ARRDI
now operational, we are in a position to address such questions through the systematic study ofsunspots in the MDI archives.

9.3.2

Probing Supergranulation in the Near-Surface Shear Layer

The subsurface mapping of supergranular motions is another interesting avenue which may
be explored through the 3-D inversion of ring-analysis measurements. Supergranular convection
occurs on scales of about 15-30 Mm, and are thought to be playing a role in the angular momentum
transport responsible for the maintenance of the near-surface shear layer (e.g., DeRosa, Gilman& Toomre 2002). The depth to which supergranular motions persist, and how their properties

163
vary with depth, must determine their ability to influence dynamics in the near-surface shear
layer. Time-distance analyses (e.g., Duval 1998; Zhao & Kosovichev 2003b) place the depth ofsupergranulation in the 8 to 15 Mm range. Frustratingly, current ring analysis inversions can
extend studies of the supergranulation down to only depths of about 5 Mm, where the deepest2◦ modes penetrate.
When using these 2◦ ring-analysis measurements, which overlap by about 1◦ , we are only justbeginning to sample the horizontal structure of the largest supergranules. The increased horizontal
resolutions accessible through HMI may allow measurements to be made at 1◦ or 0.5◦ , yielding a
relatively coarse but still reasonable resolution across the profile of a supergranule. Alternatively,
the possibility of fitting higher radial order modes in the 2◦ measurements using multi-ridge fittingwith HMI data may allow for the deeper probing of supergranulation than that currently achievableusing the f and p1 modes. The application of ARRDI to supergranulation using HMI data may
thus yield additional insights into the variation of solar supergranulation with depth.

9.3.3

Discerning the Surface Signal of the Solar Giant Cells

An interesting footnote to the subject of supergranulation is that motions at the supergranular
scale may provide insights into the nature of the deep giant cells of convective motion thought topervade the solar convection zone. These large-scale motions arise routinely in 3-D numericalsimulations of the solar convection zone and are thought to play a key role in the establishment
and maintenence of the solar differential rotation (e.g., Brun, Miesch & Toomre 2004; Miesch,
Brun & Toomre 2006; Miesch et al. 2008; Bessolaz & Brun 2011). The surface signal of the giant
cells is likely masked by the vigorous near-surface convection, but it has been suggested that theselong-lived flows impose an organizing effect on the solar supergranulation.
The correlation tracking studies of Lisle, Rast & Toomre (2004) demonstrated that supergranular cells at low latitudes tend to align with one another in the north-south direction. Morerecent studies using time-distance helioseismology suggest that the north-south alignment persistsin the polar regions as well (Nagashima et al. 2011). One theory concerning these supergranular164
patterns is that they arise along the boundaries of giant cells, which are themselves thought toalign in the north-south direction at low latitudes (though this may not explain the high-latitude
component to the supergranular alignment). A more thorough probing of supergranular cells along
these north-south lanes using ARRDI may provide further clues into the cause of this alignment.
In particular, if improvements to the ring-fitting algorithms allow us to fit 2◦ modes that penetrate
as deepy 10 Mm (roughly half a supergranule), we may be able to probe the transition between
small-scale aligned flows and deeper organizing motions. ARRDI will find application to a variety
of solar phenomena, and with HMI data now becoming available, it will soon allow us to exploit
the capabilities of ring-analysis helioseismology in a manner which was not previously possible.

Bibliography

Abt, H.A. 2000, ApJ, 544, 933
Alecian, E., et al. 2008, A&A, 481, 99A
Aurière, M., et al. 2007, A&A, 475, 1053
Aurière, M., Konstantinova-Antova, P.P., Charbonnel, C., Dintrans, B., Lignières, F., etal., 2008, A&A, 491, 499
Babcock, H.W. 1947, ApJ, 105, 105
Bagnulo, S., Landstreet, J.D., Mason, E., Andretta, V., Silaj, J. & Wade, G.A. 2006,
A&A, 450, 777
Baliunas, S.L., Donahue, R.A., Soon, W., Henry, G.W. 1998, in The Tenth CambridgeWorkshop on Cool Stars, Stellar Systems and the Sun, eds. Donahue, R.A. & Bookbinder,
J.A., ASP Conf. Ser., 154, 153
Barnes, J.R., Lister, T.A., Hilditch, R.W. & Collier Cameron, A. 2004, MNRAS, 348,
1321
Basu, S. & Antia, H.M. 2000, Solar Phys., 192, 469
Basu, S. & Antia, H.M. 2010, ApJ, 717, 488
Basu, S., Antia, H.M. & Bogart, R.S. 2004, ApJ, 610, 1157
Basu, S., Antia, H.M. & Tripathy, S.C. 1999, ApJ, 512, 458
Beck, J.G., Gizon, L. & Duvall, T.L. Jr 2002 ApJ, 575, L47
Bessolaz, N. & Brun, A.S. 2011, ApJ, 728, 115
Birch, A.C., Gizon, L., Hindman, B.W. & Haber, D.A. 2007, ApJ, 662, 730
Birch, A.C., Kosovichev, A.G. & Duvall, T.L., Jr 2004, ApJ, 608, 580
Borra, E.F., Landstreet, J.D., Mestel, L. 1982, ARA&A, 20, 191
Bouret, J.C., et al. 2008, MNRAS, 389, 75
Braithwaite, J. & Spruit, H.C. 2004, Nature, 431, 819

166
Braithwaite, J. & Nordlund, Å 2006, A&A, 450, 1077
Braithwaite, J. 2007, A&A, 469, 275
Braun, D.C. 1995 ApJ, 451, 859
Braun, D.C., Duval Jr, T.L., Labonte, B.J. 1987, ApJ, 319, L27
Braun, D.C., Duvall, T.L., Jr, Labonte, B.J., Jefferies, S.M., Harvey, J.W. & Pomerantz,
M.A. 1992, ApJ, 391, L113
Braun, D.C. & Lindsey, C. 2001, ApJ, 560, L189
Brickhouse, N.S. & Labonte, B.J. 1988, Solar Phys., 115, 43
Brown, B.P., Browning, M.K., Brun, A.S., Miesch, M.S. & Toomre, J. 2010, ApJ, 711,
424
Brown, B.P., Browning, M.K., Brun, A.S., Miesch, M.S. & Toomre, J. 2011, ApJ, in
press
Brown, S.F., Donati, J.-F., Rees, D.E., Semel, M. 1991, A&A, 250, 463
Browning, M.K. 2008, ApJ, 676, 1262
Browning, M.K., Brun, A.S., & Toomre, J. 2004, ApJ, 601, 512
Browning, M.K., Miesch, M.S., Brun, A.S. & Toomre, J. 2006, ApJ, 648, L157
Browning, M.K. & Basri, G. 2007, in Unsolved Problems in Stellar Physics, AIP Conf.
Proc., 948, 157
Brun, A.S. & Toomre, J. 2002, ApJ, 570, 865
Brun, A.S., Browning, M.K. & Toomre, J. 2005, ApJ, 629, 461
Brun, A.S., Miesch, M.S. & Toomre, J. 2004, ApJ, 614, 1073
Brun, A.S. 2007, Astron. Nachr., 328, 1137
Busse, F.H. 2002, Phys. Fluids, 14, 1301
Cally, P.S., Crouch, A.D. & Braun, D.C., 2003, MNRAS, 346, 381
Catala, C., Donati, J.-F., Shkolnik, E., Bohlender, D., Alecian, E. 2007, MNRAS, 374,
L42
Cattaneo, F. & Hughes, D.W. 1996, Phys. Rev. E, 54, 4534
Charbonneau, P. & MacGregor, K.B. 1997, ApJ, 486, 502
Charbonneau, P. & MacGregor, K.B. 2001, ApJ, 559, 1094
Charbonnel, C. & Zahn, J.-P. 2007, A&A, 476, L29
Christensen, U.R. & Aubert, J. 2006, Geophys. J. Int., 166, 97

167
Christensen-Dalsgaard, J., et al. 1996, Science, 272, 1286
Christensen-Dalsgaard, J. 2002, Rev. Mod. Phys., 74, 1073
Chou, D.-Y. & Dai, D.-C. 2001, ApJ, 559, L175
Clune, T.L., Elliott, J.R., Miesch, M.S., Toomre, J., & Glatzmaier, G.A. 1999 Parallel
Comput., 25, 361
Clyne, J., Mininni, P., Norton, A. & Rast, M. 2007, New J. Phys., 9, 301
Craig, I.J.D. & Brown, J.C. 1986, Inverse Problems in Astronomy, (Bristol, England:
Adam Hilger Ltd)
Collier Cameron, A., Robinson, R.D. 1989, MNRAS, 236, 57
Cowling, T.G. 1945, MNRAS, 105,166
DeRosa, M.L., Gilman, P.A. & Toomre, J. 2002, ApJ, 581, 1356
Deubner, F.-L. & Gough, D. 1984, ARA&A, 22, 593
Deutch, A.J. 1956, PASP, 68, 92
Dikpati, M. & Charbonneau, P. 1999, ApJ, 518, 508
Dobler, W., Stix, M. & Brandenburg, A. 2006, ApJ, 638, 336
Donati, J.-F.; Mengel, M.; Carter, B. D.; Marsden, S.; Collier Cameron, A.; Wichmann,R., 2000, MNRAS, 316, 699
Donati, J.-F., Babel, J., Harries, T.J., Howarth, I.D., Petit, P. & Semel, M., 2002,
MNRAS, 333, 55
Donati, J.-F., et al. 2006, MNRAS, 379, 629
Donati, J.-F., et al. 2006b, MNRAS, 365, L6
Donati, J.-F., et al. 2006c, MNRAS, 379, 629
Donati, J.-F., et al. 2008, MNRAS, 386, 1234
Donati, J.-F., et al. 2010, MNRAS, 409, 1347
Donati, J.-F. & Landstreet, J.D. 2009, ARA&A, 47, 333
Donati, J.-F., Morin, J., Delfosse, X., Forveille, T., Farés, R., et al. 2009, in 15th Cambridge Workshop on Cool Stars, Stellar Systems and the Sun, AIP Conf. Proc. 1094,
130
Donati, J.-F., Moutou, C., Farés, R., Bohlender, D., Catala, C., et al. 2008, MNRAS,
385, 1179
Dorren, J.D. & Guinan, E.F. 1982, Astron. J., 87, 1546

168
Duez, V. & Mathis, S. 2009, ApJ, submitted
Duvall, T.L., Jr 1998, in Structure and Dynamics of The Interior of the Sun and Sun-likeStars, (eds.) Korzennik, S., Wilson, A., SOHO 6 GONG 98 Workshop, ESA Conf. Proc.,
418, 581
Duvall, T.L., Jr, D’Silva, S., Jefferies, S.M., Harvey, J.W. & Schou, J. 1996, Nature, 379,
235
Duvall, T.L., Jr & Gizon, L. 2000, Solar Phys., 192, 177
Duvall, T.L., Jr, Jefferies, S.M., Harvey, J.W. & Pomerantz, M.A. 1993, Nature, 362,
430
Duvall, T.L., Jr, Kosovichev, A.G., Scherrer, P.H., Bogart, R.S., Bush, R.I., De Forest,
C., Hoeksema, J.T., et al. 1997, Solar Phys., 170, 63
Featherstone, N.A., Hindman, B.W., Haber, D.A. & Toomre, J. 2004, in Helio- andAsteroseismology: Towards a Golden Future, Editor, D. Danesy, 559, 428
Featherstone, N.A., Hindman, B.W., Haber, D.A. & Toomre, J. 2006, in Proceedings of
SOHO 18/GONG 2006/HELAS I, Beyond the spherical Sun, Editor, M. Thompson, 624,
131
Featherstone, N.A., Hindman, B.W., Thompson, M.J. 2011, in GONG 2010 - SoHO 24:
A new era of seismology of the Sun and solar-like stars, in press
Flowers, E. & Ruderman, M.A. 1977, ApJ, 215, 302
Freyhammer, L.M., et al. 2008, MNRAS, 389, 441
Ghizaru, M., Charbonneau, P. & Smolarkiewicz, P.K. 2010, ApJ, 715, L133
Giles, P.M., 1999, Time-Distance Measurements of Large-Scale Flows in the Solar Convection Zone, Ph.D. Thesis, Stanford Univ.
Giles, P.M., Duvall, T.L., Jr, Scherrer, P.H. & Bogart, R.S. 1997, Nature, 390, 52
Giles, P.M., Duvall, T.L., Jr, T.L. & Scherrer, P.H. 1998, in Structure and Dynamics ofthe Interior of the Sun and Sun-like Stars SOHO 6/GONG 98 Workshop, eds., Korzennik,
S. & Wilson, A., ESA, Noordwijk, Netherlands, 418, 775
Gilman, P.A. & Glatzmaier, G.A. 1980, ApJS, 45, 335
Gizon, L. 2003, Probing Flows in the Upper Solar Convection Zone, Ph.D. thesis, Stanford
Univ.
Gizon, L. & Birch, A.C. 2002, ApJ, 571, 966
Gizon, L. & Birch, A.C. 2005, LRSP, 6, Local Helioseismology
Gizon, L., Birch, A.C. & Spruit, H.C. 2010, ARA&A, 48, 289
Gizon, L., Duvall, T.L., Jr & Larsen, R.M. 2000, J. Astrophys. Astr., 21, 339

169
Gizon, L., Duvall, T.L., Jr & Larsen, R.M. 2001, in Recent Insights into the Physics
of the Sun and Heliosphere: Highlights from SOHO and Other Space Missions, Proc. of
IAU Symposium 203, ed., P. Brekke, F. Bernhard & J.B. Gurman., 189
Gizon, L., et al. 2009, Space Sci. Rev., 144, 249
Gizon, L. & Rempel, M. 2008, Solar Phys., 251, 241
Glatzmaier, G.A., Coe, R.S., Hongre, L. & Roberts, P.H. 1999, Nature, 401, 885
Goldreich, P. & Keeley, D.A. 1977, ApJ, 212, 243
González Hernández, I., Hill, F. & Lindsey, C. 2007, ApJ, 669, 1382
González Hernández, I., Kholikov, S., Hill, F., Howe, R. & Komm, R.W. 2008, Solar
Phys.
González, Hernández, I., Komm, R.W., Hill, F., Howe, R., Corbard, T. & Haber, D.A.
2006, ApJ, 638, 576
González, Hernández, I., Patrón, J., Bogart, R.S. & The SOI Ring Diagram Team 1999,
ApJ, 510, L153
González, Hernández, I., Patrón, J., Cortés, T.R., Bogart, R.S., Hill, F. & Rhodes Jr.,
E.J. 2000, ApJ, 535, 454
Goode, P.R., Gough, D. & Kosovichev, A. 1992, ApJ, 387, 707
Gough, D.O. 1969, J. Atmos. Sci., 26, 448
Gough, D.O. & Toomre, J. 1991, ARA&A, 29, 627
Gregory, S.G., Jardine, M., Gray, C.G. & Donati, J.F. 2010, Rep. Prog. Phys. 73, 1Haber, D.A., Hindman, B.W., Toomre, J., Bogart, R.S. & Hill, F. 2001, in Proceedingsof the SOHO 10/GONG 2000 Workshop: Helio- and asteroseismology at the Dawn of
the Millennium, Ed., A. Wilson, ESA SP-464, 209
Haber, D.A., Hindman, B.W., Toomre, J., Bogart, R.S., Larsen, R.M., & Hill, F. 2002,
ApJ, 570, 855
Haber, D.A., Hindman, B.W., Toomre, J., Bogart, R.S., Thompson, M.J. & Hill, F.
2000, Solar Phys., 192, 335Hale, G.E., Ellerman, F., Nicholson, S.B. & Joy, A.H. 1919, ApJ, 49, 153
Hale, G.E. & Nicholson, S.B. 1925, ApJ, 62, 270
Hindman, B.W., Haber, D.A. & Toomre, J. 2006, ApJ, 653, 725
Hindman, B.W., Haber, D.A. & Toomre, J. 2009, ApJ, 698, 1749
Hathaway, D.H. 1996, ApJ, 460, 1027

170
Harvey, J.W., Hill, F., Hubbard, R., Kennedy, J.R., Leibacher, J.W., et al. 1996, Science,
272, 1284
Henrichs, H.F., et al. 2000, in Magnetic Fields of Chemically Peculiar and Related Stars,
Proc. Int. Meet., eds. YuV Glagolevskij, II Romanyuk, 57
Hindman, B.W., Gizon, L., Duvall, T.L., Jr, Haber, D.A. & Toomre, J. 2004, ApJ, 613,
1253
Howard, R. & LaBonte, B.J. 1980, ApJ, 239, L33
Howe, R., Christensen-Dalsgaard, J., Hill, F., Komm, R.W., Larsen, R.M., Schou, J.,
Thompson, M.J. & Toomre, J. 2000, ApJ, 533, L163
Hussain, G.A., et al. 2009, MNRAS, 398, 189
Hurlburt, N.E. & Rucklidge, A.M. 2000, MNRAS, 314, 793
Hurlburt, N.E. & DeRosa, M. 2008, ApJ, 684, L123
Johns-Krull, C.M., 1997, ApJ, 664, 975
Kochukhov, O., Drake, N.A., Piskunov, N. & de la Reza, R. 2004, A&A, 424, 935
Kochukhov, O., et al. 2004, A&A, 414, 613
Kochukhov, O. & Bagnulo, S. 2006, A&A, 450, 763
Komm, R., Corbard, T., Durney, B.R., González Hernández, I., Hill, F., Howe, R. Toner,
C. 2004, ApJ, 605, 554
Komm, R., Howe, R. & Hill, F. 2009b, Solar Phys., 258, 13
Komm, R., Howe, R., Hill, F., González Hernández, I. & Toner, C. 2005, ApJ, 630, 1184
Komm, R., Howe, R., Hill, F., González Hernández, I. 2009a, Solar Phys., 254, 1
Komm, R.W., Howe, R., Hill, F., Miesch, M.S. Haber, D.A. & Hindman, B.W. 2007,
ApJ, 667, 571
Komm, R., Morita, S., Howe, R. & Hill, F., 2008, ApJ, 672, 1254
Kosovichev, A.G., 1996, ApJ, 461, L55
Kosovichev, A.G. & Duvall, T.L., Jr 1997, in SCORE ’96 Solar Convection and Oscillation and Their Relationship, eds. F.P. Pijpers, J. Christensen-Dalsgaard & C.S. Rosenthal, (Kluwer: Dordrecht), 241
Korzennik, S.G., Rabello-Soares, C. & Schou, J., 2004, ApJ, 602, 481
Kosovichev, A.G. & The HMI Science Team 2007, AN, 328, 339
Kosovichev, A.G. & Schou, J. 1997, ApJ, 482, L207
Krause, F. & Oetken, L. 1976, in Physics of Ap Stars, ed. W.W. Weiss, H. Jenkner, &H.J. Wood (Vienna: Univ. Wien), 29

171
Kuang, W. & Bloxham, J. 1999, J. Comput. Phys., 153, 51
Landstreet, J.D., Bagnulo, S., Andretta, V., Fossati, L., Mason, E., Silaj, J. & Wade,
G.A. 2007, A&A, 470, 685
Landstreet, J.D., et al. 2008, A&A, 481, 465
Leibacher, J.W. & Stein, R.F. 1971, ApJ, 7, 191
Leroy, J.L. 1995, A&AS, 114, 79
Libbrecht, K.G., Popp, B.D., Kaufman, J.M. & Penn, M.J. 1986, Nature, 323, 235
Lindsey, C. & Braun, D.C. 1997, ApJ, 485, 895
Lindsey, C. & Braun, D.C. 2000, Science, 287, 1799
Lindsey, C. & Braun, D.C. 2005a, ApJ, 620, 1107
Lindsey, C. & Braun, D.C. 2005b, ApJ, 620, 1118
Lindsey, C., Braun, D.C., Jefferies, S.M., Woodard, M.F., Fan, Y., Gu, Y. & Redfield,
S. 1996, ApJ, 470, 636
Lisle, J.P., Rast, M.P. & Toomre, J. 2004, ApJ, 608, 1167
MacDonald, J. & Mullan, D.J. 2004, MNRAS, 348, 702
MacGregor, K.B. & Cassinelli, J.P. 2003, ApJ, 586, 480
Markey, P., Tayler, R.J. 1973, MNRAS, 163, 77
Mason, D., Komm, R., Hill, F., Howe, R., Haber, D. & Hindman, B.W. 2006, ApJ, 645,
1543
Mathys, G., Hubrig, S., Landstreet, J.D., Lanz, T. & Manfroid, J. 1997, A&AS 123, 353
Mathys, G. 2001, in Magnetic Fields Across the Hertzsprung-Russel Diagram, ed. G.
Mathys, S. K. Solanki, & D. T. Wickramasinghe (San Francisco: ASP), in ASP Conf.
Ser. 248, 267
Maury, A.C. 1897, Harvard Ann., 28, 96
Maurya, R.A. & Ambastha, A. 2010, ApJ, 714, L196
Mestel, L. 1999, Stellar Magnetism, (Oxford: Clarendon Press)
Miesch, M.S., Brun, A.S. & Toomre, J. 2006, ApJ, 641, 618
Miesch, M.S., Brun, A.S., DeRosa, M.L. & Toomre, J. 2008, ApJ, 673, 557
Morel, P. 1997, A&AS, 124, 597
Moss, D.L. 2004, A&A, 414, 1065

172
Moss, D.L. 2001, in Magnetic Fields Across the Hertzsprung-Russel Diagram, ed. G.
Mathys, S. K. Solanki, & D. T. Wickramasinghe (San Francisco: ASP), in ASP Conf.
Ser. 248, 305
Nagashima, K., Zhao, J., Kosovichev, A.G. & Sekii, T. 2011, ApJ, 726, L17
Neiner, C., et al. 2003, A&A, 406, 1019
Noyes, R. W., Hartmann, L. W., Baliunas, S. L., Duncan, D. K., Vaughan, A. H. 1984,
ApJ, 269, 763
Nye, A., Bruning, D. & Labonte, B.J. 1988, Solar Phys., 115, 251
Olson, P., Christensen, U. & Glatzmaier, G.A. 1999, J. Geophys. Res., 104, 10383
Parker, E.N. 1979, ApJ, 230, 905
Parker, E.N. 1992, ApJ, 390, 290
Parker, E.N. 1993, ApJ, 408, 707
Petit, P., Wade, G.A., Drissen, L., Montmerle, T., Alecian, E. 2008, MNRAS, 387, L23
Piskunov, N. & Kochukhov, O. 2002, A&A, 381, 736
Pizzolato, N., Maggio, A., Micela, G., Sciortino, S., Ventura, P. 2003, A&A, 397, 147
Power, J., Wade, G.A., Aurière, M., Silvester, J., & Hanes, D. 2008, Contrib. Astron.
Obs. Skalnate Pleso, 38, 443
Prendergast, K.H. 1956, ApJ, 123, 498
Preston, G.W. 1967, ApJ, 150, 547
Preston, G.W. 1971, PASP, 83, 571
Radick, R.R., Mihalas, D., Hartmann, L., Worden, S.P., Africano, J.L., Klimke, A. &Tyson, E.T. 1982, PASP, 94, 934
Reiners, A., Basri, G. & Browning, M., 2009, ApJ, 692, 538
Rempel, M. 2011, ApJ, submitted
Rempel, M., Schüssler, M., Cameron, R.H. & Knölker, M., 2009, Science, 325, 171
Richards, M.T. 1990, ApJ, 350, 372
Roberts, P.H. & Glatzmaier, G.A. 2000, Rev. Mod. Phys., 72, 1081
Robinson Jr., R.D. 1980a, ApJ, 239, 961
Robinson, R. D., Worden, S. P. & Harvey, J. W. 1980b, ApJ, 236, L155
Royer, F., Zorec, J. & Gómez, A.E. 2007, A&A, 463, 671
Saar, S.H. 1988, ApJ, 324, 441

173
Saar, S.H. & Linskey, J.L. 1985, ApJ, 299, L47
Sakuraba, A. & Kono, M. 1999, Phys. Earth Planet. Interiors, 111, 105
Scherrer, P.H., Bogart, R.S., Bush, R.I., Hoeksema, J.T., Kosovichev, A.G., et al. 1995,
Solar Phys., 162, 129
Schou, J. & Bogart, R.S. 1998, ApJ, 504, L131
Schou, J. 1999, ApJ, 523, L181
Schou, J. et al. 1998, ApJ, 505, 390
Schüssler, M. & Rempel, M. 2005, A&A, 441, 337
Semel, M. 1989, A&A, 225, 456
Sheeley Jr, N.R. 1969, Solar Physics, 9, 347
Sheeley Jr, N.R. 1972, Solar Physics, 25, 98
Shibahashi, H., 2004, in Multi-Wavelength Investigations of Solar Activity, IAU Symposium, 223, 23, Ed., A.V. Stepanov, E.E. Benevolenskaya & A.G. Kosovichev, Cambridge,
UK: Cambridge University Press.
Shorlin, S.L.S., Wade, G.A., Donati, J.-F., Landstreet, J.D., Petit, P., et al. 2002, A&A,
392, 637
Solanki, S.K. 2003, Astron. Astr. Rev., 11, 153
Spiegel, E.A. & Zahn, J.-P. 1992, A&A, 265, 106
Spruit, H.C. 2002, A&A, 381, 923
Spruit, H.C. 2003, Solar Phys., 213,1
Stibbs, D.W.N. 1950, MNRAS, 110, 395
Strassmeier, K.G., Stepień, K., Henry, G.W. & Hall, D.S. 1999, A&A, 343, 175
Sun, M-T., Chou, D-Y., Lin, C-H. & the TON Team 1997, Solar Phys., 176, 59
Tayler, R.J. 1973, MNRAS, 161, 365
Thompson, M.J., Christensen-Dalsgaard, J., Miesch, M.S. & Toomre, J. 2003, ARA&A,
41, 599
Thompson, M.J., Toomre, J. et al., 1996, Science, 272, 1300
Title, A.M., Tarbell, T.D., Topka, K.P., Ferguson, S.H., Shine, R.A. & the SOUP Team
1989, ApJ, 336, 475
Ulrich, R.K. 1970, ApJ, 162, 993
Ulrich, R.K. 2010, ApJ, 725, 658

174
Valenti, J.A., Marcy, G.W. & Basri, G. 1995, ApJ, 439, 939
Wade, G.A., Alecian, E., Grunhut, J., Catala, C., Bagnulo, S., Folsom, C.P. & Landstreet, J.D. 2009, arXiv:0901.0347
Wilson, O.C. 1978, ApJ, 226, 279
Wittmann, A.D. & Xu, Z.T. 1987, A&AS, 70, 83
Wright, A.A.E. 1973, MNRAS, 162, 339
Zahn, J.-P, Brun, A.S. & Mathis, S. 2007, A&A, 474, 145
Zhao, J. & Kosovichev, A.G. 2003, ApJ, 591, 446
Zhao, J. & Kosovichev, A.G. 2003b, in Local and Global Helioseismology: The Present
and Future, ed., H. Sawaya-Lacoste, ESA Conf. Proc., 517, 417
Zhao, J. & Kosovichev, A.G. 2004, ApJ, 603, 776
Zhao, J., Kosovichev, A.G. & Duvall, T.L., Jr 2001, ApJ, 557, 384
Zhao, J., Kosovichev, A.G., Duvall, T.L., Jr 2004, ApJ, 607, L135
Zhao, J., Kosovichev, A.G. & Sekii, T. 2010, ApJ, 708, 304