Current Research

Numerical Modeling of Solar Dynamo

I. Overview :

A grand challenge in the astrophysics is the origin of self-organizing property of the magnetic field in highly turbulent flows. The solar magnetism is the front line in this area. The solar magnetic field shows a remarkable spatiotemporal coherence though it is generated by highly turbulent convection (Ra ~ 1020) operating within its interior. Our understanding on the solar magnetism has been accelerated over the past decade in response to broadening, deepening and refining of numerical dynamo modeling. However, it is still unclear what dynamo mode is excited in the solar interior and how it regulates magnetic cycle. To unveil the mystery of the solar magnetism, I’m currently working on the global and local numerical dynamo modelings of the solar dynamo with state-of-the-art massively-parallel supercomputers, such as K-computer at RIKEN and CRAY XC 30 at NAOJ.


II. Development of A New Solar Dynamo Simulation Code with Yin-Yang Grid :

Existing simulation models for the global solar dynamo are mostly spectral-based type with the spherical harmonics expansion. However, it is confronted with the parallel scaling difficulty when O(104) or more CPU cores are used. A different approach toward massively-parallel solar dynamo modeling is required for the present peta- or coming exa-scale era. This motivated us to develop a global solar dynamo simulation code based on the grid point-based approach. The Yin-Yang grid (see Kageyama & Sato 2004) is adopted as a coordinate system in our code. Since the Yin-Yang Grid is free from the coordinate singularity and the grid convergence around there, we can avoid severe time-step constraint due to the CFL condition. This newly-developed solar dynamo code opens up the way ensuring compatibility between finer grid spacing and longer time-integration which are required to resolve turbulent dynamos more accurately. With this code, we solved numerically a sun-like MHD dynamo convection system in Masada et al (2013).

III. Theoretical and Numerical Modelings of Flaring Activity on Compact Star

■ Solar-type Model for Magnetar Giant Flare (Masada et al. 2010; Matsumoto et al. 2011)

Soft gamma-ray repeaters (SGRs) are believed as the population of the ultra-strongly magnetized neutron stars, so-called “Magnetar”. The giant flare with enormous energy (~ 1046 erg) and long bursting duration occurs infrequently on SGRs. In Masada et al. (2010), we presented a theoretical model describing the magnetar giant flare based on solar flare/coronal mass ejection theory. The nonlinear expansion of the magnetic field after triggering the solar-type flare on such a compact star was studied by co-authors and I with using a relativistic MHD simulation (Matsumoto, Masada et al. 2011). Theoretical and numerical models are shown in left panel (a) and (b). A key process in our model is the energy transport by the high-energy photon. The photon flux dredges up the baryon from magnetar crust and produces “magnetar prominence”, which is finally erupted due to the loss of equilibrium triggered by the magnetic reconnection at the main burst stage. Our model is getting a lot more attention lately as one of standard scenarios of the magnetar giant flare (see Tiengo et al. 2013, Nature) The predictions of our model will be tested by future observations even though such an event is very rare.
 

III. Sun-like Differential Rotation, Large-scale B-field, and its Quasi-Regular Cycling (Masada et al. 2013) :

The final goal of my challenge is a construction of self-consistent theory accounting for the observed spatiotemporal structures of the “fields” in the Sun. Fundamental large-scale spatiotemporal structures that remain to be explained are the azimuthal average of the azimuthal flow vφ, and the azimuthal average of the azimuthal magnetic component Bφ. The vφis characterized by the conical iso-rotation profile with the accelerated equator in the meridian plane and the thin tachocline layer with a steep radial angular velocity gradient. The Bφis characterized by an anti-symmetric profile about the equator, the polarity reversals with the periodicity of 22 yrs, and the equatorward migration. With using the Yin-Yang solar dynamo code we developed, spherical solar dynamo simulations are performed in Masada et al. (2013).  Outcomes are shown in bottom figures: (a) the rotation profile vφ/r in meridian plane, (b) & (c) the snapshots of Bφ in different phases at the mid radiative zone and (d) the snapshot of Bφ at the mid convection zone. While the details are still far from the solar profiles, the three key features of the solar interior, sun-like vφ, large-scale Bφ, and its quasi-regular polarity reversals, are spontaneously reproduced. The higher resolution simulation with adopting a realistic solar internal structure is the next step of our newly-developed solar dynamo code.
 

Research Statement

 

TIME (t/τcv)

0

200

400

600

800

1000

1200

1400

Radius

0.0

0.4

0.8

1.2

1.6

2.0

(b) r = 0.62Rsun

(c) r = 0.62Rsun

(d) r = 0.85Rsun

(a)

  1. (b)Time-radius diagram of Bx. Since there is no symmetry breaking in horizontal directions in our model, stochastic properties of turbulent convections play an essential role in sustaining the large-scale dynamo.

Schematic View of MRI-dead zone scenario for GRBs

(b)

(a)

(a) Theoretical Model (Masada et al. 2010)

(b) Numerical Model (Matsumoto et al. 2011)

IV. Large-scale Oscillatory Dynamo in Rotating Stratified Convection (Masada & Sano 2014a,b) :

To gain deeper understanding of the fundamental dynamo mode of turbulent convections, and to get a

clue to the dynamo mode operated in the global system, we studied turbulent convective dynamos in

local Cartesian domain in Masada & Sano (2014a). Our simulation model consists of three-layers like 

as the solar interior: bottom radiative zone, middle  convection zone, and top cooling layer. Panel (a)

shows 3D surface visualization of the convective velocity. The time-radius diagram for the horizontally-

averaged horizontal magnetic field Bx is presented in panel (b). The oscillatory large-scale magnetic

field is spontaneously organized in the bulk of convection zone. We found that although the large-scale

dynamo is operated even without the radiative zone, the magnetic cycle period depends on the presence

of the stable layers. We are now studying this large-scale dynamo and its response to the density contrast

with a special emphasis on the effect of the α-effect and its latitudinal dependence.

Accomplishments in Past (Selected Papers)


I. Magneto-Rotational Instability (MRI) in Supernovae and Gamma-Ray Bursts

■ Construction of the Linear Theory of the MRI in Supernovae (Masada et al. 2006, 2007a)

In ultra-dense and hot baryonic matter constituting nascent neutron stars just after massive-stellar core-collapse, the neutrino plays a prominent role in momentum, heat, and leptonic transports. In the works of Masada et al. (2006, 2007a), we developed the linear theory of the MRI (including both axisymmetric and non-axisymmetric modes) with taking account of the density stratification in the proto-neutron star and the neutrino transport effects.

The finding is that while the MRI can grow even in the stellar system consisting of the neutrino- dominated matter, its growth rate is impacted by the neutrino radiation, i.e., the MRI is drastically suppressed by the neutrino viscosity if the magnetic field in the proto-neutron star is weaker than the critical value of B = 1012 G. Therefore, the onset of the MRI depends strongly on the magnetic field strength in the pre-collapse progenitor star. Not only is it important for the magnetorotational mechanism of the supernova explosion, but also it might be important to account for two distinct populations of the neutron star: ordinary pulsar with 1012 G and strongly magnetized one, called as “magnetar”, with 1014-15 G.

■ MRI-Dead Zone Scenario Accounting for Evolution History of GRBs (Masada et al. 2007b)

GRBs are the most energetic events in the universe, and are considered to be powered by a hyperaccretion

of matter onto a stellar-mass black hole formed in the aftermath of stellar core-collapse or compact objects

merger. The linear theory of the MRI with the neutrino transport is applicable to such a hyperaccretion disk.

We discussed, in Masada et al. (2007b), the possible inactivation of the MRI in the hyperaccreting flow.

In it, we discovered inner MRI-dead zone region where is neutrino-opaque and thus stable to the MRI due

to the neutrino viscosity. Since the MRI is expected to be active only in outer neutrino-transparent region,

it is speculated that the efficiency of the angular momentum transport, described by conventional α, due to

the MRI-driven turbulence would vary in the radial direction: “lower α in the inner MRI-dead zone and

higher α in the outer MRI-active zone”.

Under this speculation, the “MRI-dead zone scenario” which can predict the evolution history of the GRB

is proposed. The schematic view of the scenario is shown in right panel. This scenario consists of two stages:

At the first stage, the matter in the outer MRI-active zone is accreted and thus accumulated into the inner

MRI-dead zone. After the dead zone gains a large amount of mass and becomes gravitationally unstable,

a violent mass accretion is triggered by the gravitational torque. This is the subsequent stage. These two

accretion phases are repeated episodically, and then power intermittent jets which are believed as the origin

of internal shocks producing the observed short-term variability in the prompt emission of the GRB. A lot of

important features in timescales and energetics of the observed emission of the GRB are well- predicted by

the activity of the central engine in the framework of my MRI-dead zone scenario.

■ Nonlinear Study of the MRI toward the Application to Supernovae (Masada et al. 2012)

Generally, the MRI-driven turbulence can sustain a turbulent heating much larger than kinematic viscous

heating. Therefore, it can contribute to the enhancement of neutrino-luminosity and then facilitate a neutrino-

driven explosion if it develops in the proto-neutron star. Bearing in mind the application of this MRI-facilitated

mechanism to the supernova explosion, in Masada et al. (2012), we examined quantitatively the turbulent

heating rate sustained by the MRI with local shearing box simulations. Shown in (a) of right panel is volume

rendering visualization of the MRI-turbulence for the models with different shear rate (q), which represents

the degree of differential rotation of the proto-neutron star. The larger the shear rate, the stronger the MRI-driven

turbulence is sustained. The turbulent stress by the MRI is proportional to the square-root of the shear-vorticity

ratio [q/ (2-q)], and then provides a volume-integrated turbulent heating rate as shown (b) in right panel. The

predictor function for the MRI-sustained turbulent heating derived in this work is applicable as a subgrid-scale

turbulent heating model to the global simulation of the supernova explosion.

II. Effects of Magnetic Prandtl number “PrM” on the MRI and Accretion Disks


■ Masada & Sano (2008) - Why is the MRI dependent on PrM ? -

The final goal of the MRI research should be the elucidation of its saturation physics. It cannot be avoided for

it to understand the role of the microscopic diffusivity. This motivated us to determine the impacts of the

microscopic viscosity and resistivity on the MRI. In this work, we found, for the first time, that the nonlinear

behavior of the MRI in the viscous regime is entirely different from that in the resistive regime. Furthermore,

we showed that the characteristics of the linear dispersion relation of the MRI can account nicely for the difference

in the nonlinear nature of the MRI between two regimes. Our theory can predict naturally the PrM-dependence

of the MRI-driven turbulence which was reported in several earlier simulations: As shown in right panel, it is

expected from our theory that the critical Lundquist number to activate the MRI-turbulence depends on PrM

in the regime PrM > 1.


■ Taahashi & Masada (2011) - PrM-dependent Disk Models - .

0.6

“If MRI truly depends on PrM, how it impacts on disk models ?”, that was the motivation of this work.

To answer it, we examined the stability of a geometrically-thin accretion disk model using a “non-standard”

α-prescription of αPrM δ, to thermal, secular and gravitational instabilities. Our finding is that, when

using Spitzer’s microscopic diffusivities, the PrM-dependent disk model becomes more unstable to all

the instabilities than the standard disk model. Shown in right panel is the stability criterion for the thermal

and secular instabilities. The larger δ yields the more unstable disk. The temporal evolution of the PrM-

dependent accretion disk was studied recently by Potter & Balbus (2014), and they showed that various

observations are well-explained by such a PrM-dependent disk model.


IV. The MRI as a Possible Origin of the Sun’s Differential Rotation V.

■ Why does the Sun have conical rotation profile ? (Masada 2011)

The Sun’s differential rotation deduced from helioseismic inversion tells us that our knowledge of angular momentum transports in

the solar interior is still insufficient (See Miesch 2005 for reviews). What is the most puzzling us is why the solar differential rotation in the convection zone has “conical” iso-rotation profile. The solid line shown in left panel is the helioseismic data of the iso-rotation contour of the solar convection zone on the meridian plane. The deviation from the Taylor-Proudman state with cylindrical isotach (∂Ω/∂z = 0) implies that there is an unaccountable latitudinal entropy gradient in the solar interior (see Pedlosky 1987). In Masada (2011), I pointed out a possible role by the MRI-driven turbulence in causing such a latitudinal entropy gradient. The color map on panel (a) demonstrates the region where the MRI-turbulence can be activated (MRI-unstable zone). You can find the MRI-unstable  region is confined in the higher latitude tachocline. In contrast, the bottom panel shows the latitudinal entropy gradient required for the observed rotation profile when assuming the “thermal wind balance” in the solar interior. The MRI-active region at the bottom of the convection zone closely overlaps with the area which requires a steep entropy rise. In my estimation, the turbulent heating sustained by the MRI is large enough to maintain the entropy variation which is required for the observed rotation profile. There is of course a lot of models for explaining the solar conical rotation profile (c.g., Kitchatinov & Rudiger 1995; Rempel 2005). Still no one knows whether my idea can survive against future observational challenges, my finding would be very suggestive and facilitate our understanding of the physics lying behind the solar rotation profile.