音速抑制法による対流数値計算

音速抑制法による
太陽内部の対流数値計算
横山研究室修士2年
堀田英之
共同研究者:横山央明[1], M. Rempel[2], Y. Fan[2]
[1]東京大学 [2]HAO(High Altitude Observatory)
太陽黒点数11年周期
太陽の黒点数(面積)は11年の周期を持って変動している
90
30
EQ
30
90
1880 1890 1900 1910 1920 1930 1940 1950 1960 1970 1980 1990 2000 2010
1870 1880 1890 1900 1910 1920 1930 1940 1950 1960 1970 1980 1990 2000 2010
Date
Hathaway/web
太陽ダイナモ
太陽の磁場活動を駆動していると考えられているのが太陽ダイナモ
ダイナモとは運動エネルギーを磁場エネルギーに変換すること
Ω効果
ポロイダル方向
α効果
トロイダル方向
(桜井ほか｢太陽｣)
太陽ダイナモはΩ効果とα効果
という2つの機構で成り立って
いると考えられている。
Ω効果：差動回転により磁場を
引き延ばす。ポロイダル磁場
からトロイダル磁場
α効果：いつどこで、効果が
効いているのか諸説あるが、
コリオリ力によって、磁場の方向
を変える。トロイダル磁場から
ポロイダル磁場
修士課程での研究
磁束輸送ダイナモにおける表面乱流拡散の
重要性(Hotta & Yokoyama, 2010a, ApJ, 709,
1009)
太陽大規模磁場の対称性(Hotta & Yokoyama,
2010b, ApJL, 714, L308)
自転角速度の大きい太陽型星の内部角速度
分布(Hotta & Yokoyama, 2010c in prep)
音速抑制法による太陽内部の対流数値計算
(Hotta et al., 2011, in prep)
磁束輸送ダイナモにおける
表面乱流拡散の重要性
Hotta & Yokoyama, 2010a
ApJ
線：
ポロイダル磁場の磁力線
色のコンター：
トロイダル磁場の強さ
今までのモデルでは
極付近に観測と合わない
磁場が生成されて
しまったが、対流層上部
に強い乱流拡散がある
と仮定することでこれを
避けることができる。
太陽大規模磁場の対称性
どうやったら磁場は双極子になるのか
dipole
Hotta & Yokoyama, 2010b, ApJL
表面付近に強い乱流拡散が
あると、双極子磁場になることを
発見。
Symmetric Parameter
quadrupole
diffusivity profile
明るい色:対称な磁場
暗い色:反対称な磁場
回転の速い星
での差動回転
自転角速度が大きくなると
差動回転がTaylor-Proudman
状態に近くなる。
Our sun
Taylor-Proudman
non-Taylor-Proudman
parameter
1
2
Pntp  2 2 
dV
R 0 z
1
2
4
8
16
real solar DR by observation
最終的な動機付け
(博士課程での研究)
太陽内部角速度分布を再現する自己無撞着なモデルの確立
太陽内部の対流を解きつつ、角運動量輸送、放射層での
物理を調べる。Miesch+(2000)などでは高解像度のシミュレーションを
行うことで新たな知見を得ている。
450 nHz
425 nHz
400 nHz
375 nHz
350 nHz
325 nHz
300 nHz
日震学により得られた太陽の差動回転
(Thompson+ 2003)
Taylor-Proudman state
(Rempel 2005)
太陽の対流層を計算するのは
なにが難しいのか？
音速が速いことである。対流層の底で音速が200 km/s。
それに対して、対流の速度は100 m/s。
対流の時間スケールもしくは、11年スケールの現象を知りたいが、
CFL条件により時間ステップは音速で解かなければいけない。
対流層の底でスケールハイトが70000 km 10点で解像したとき
時間ステップは、40秒ほどとなる。
これを100TFLOPSの計算機で11年分計算するのに数10年かかる。
この難しさを回避する何らかの方法をとらなければ、今の計算資源
では太陽対流層を計算して、科学的結果を得ることができない。
アネラスティック近似
100 m/s
200 km/s
velocity
音速   0 v   0
対流速度
他の方法?
アネラスティック近似では、
音速を無限大とする。マッハ数が
1より十分小さい時、良い近似
である。しかし、数値計算では
通信が多く発生し、2000CPU
くらいまでしかスケールしない
(M. Miesch氏談)
通信の少ない方法は逆の方向
に行くことで達成できるのでは?
アネラスティック
近似
M. Miesch氏提供
音速抑制法
RSS(Reduced Sound Speed) technique
以下の方法で、人工的に連続の式を変えることで音速を遅くする

   v 
t
 v 

    2 
t
 
この方法で音速はξ 倍遅くなる。
有限差分法も使える上に、陽解法で全ての方程式が解ける。
しかし、この方法は対流を計算する上で妥当な技法なのだろうか？
本研究では、近似無し、RSS法、anelastic近似での2次元、3次元
対流を比較し、RSS法の妥当性を問う。
対流とは？
粘性熱伝導がなければ、温度勾配が
断熱温度勾配より急な時、対流が起こる。
(Schwartzchild条件)
エントロピーが上に行くに従って、減少して
いるとき、対流が起こるといいかえることも
できる。
エントロピーの勾配を見積もる量として
superadiabaticity という量がある
d (logT )

d (log P)
    ad
superadiabaticityを用いるとエントロピーの勾配は
c p
ds

dz
Hp
方程式系
 0 v 
1
    2  連続の式
t
 
v
p1 1
 ( v   ) v 
 g    Π 運動方程式
t
0 0
s1
1
 ( v  )(s0  s1 ) 
  ( K 0T0s1 )  Π     v
t
 0T0
エネルギー方程式
 1 s1 
p1
背景場は静水圧平衡を仮定




して、そこからの擾乱成分を
 c 
p0
p 
 0
解く。
Reynolds数、Prandtl数、スケールハイトで規格化した計算領域、
superadiabaticityが依存しないパラメーター
数値スキーム
1. 4次精度の中央差分で空間微分
2. 4次精度のRunge-Kuttaで時間積分
3. 人工粘性はRempel(2009)で使われたもの
modified Lax-Wendroff
(CANS) 時間空間2次精度
空間・時間4次精度
スキーム
2.18(Hr)
2次元定常対流
Re=260R Pr=1
8.72 (Hr)
解像度:400x100
音速と典型的対流速度の比は
 (~1/300: bottom)
RMS速度
Horizontal RMS velocity
Height
vrms 
2
v
 dS
S
Vertical RMS velocity
Height
マッハ数
max(v/c)
max(vrms/c)
ξ=1
0.022
0.016
ξ=10
0.28
0.17
ξ=30
0.69
0.54
ξ=50
1.97
1.34
C : reduced sound speed
2.18(Hr)
2次元非定常対流
Re=500R Pr=1
解像度:400x100
音速と典型的対流速度の比は
 (~1/300: bottom)
8.72 (Hr)
full movie
RMS速度
特徴的な速度はあっている、しかし
時間平均が短いため、定常の時に比べて一致は悪い
緑-x方向速度
黒-z方向速度
赤-y方向速度
青-密度
3次元非定常対流
3次元対流の場合。
完全に非定常な振る舞い
を見せる
Re=260 Pr=1
4.36(Hr)
200 grid
4.36(Hr)
200 grid
2.18(Hr)
100 grid
vertical
velocity
entropy
vertical
velocity
RMS速度の時間発展
RMS速度も非線形で非定常な振る舞いをみせる
時間平均したRMS速度
時間平均したRMS速度はξによらない
3次元対流での音速抑制法の妥当性は証明された
太陽内部のsuperadiabaticity
(Stix “The Sun”より)
r/rsun
0.999
0.990
0.949
0.818
0.710
δ(super
-adiabaticity)
5.7x10-2
2.6x10-4
8.1x10-6
4.7x10-7
1.9x10-8
Mach数
の逆数
4.2
62
350
1400
7200
まとめ
• 太陽対流層の計算をする上で、音速が非常に速
いことが大きな困難である
• アネラスティック近似を用いて、今まで計算され
ているが、多くのCPUを使う計算ではスケールし
ないことが知られている
• 音速抑制法を用いればその困難を回避すること
ができる
• これからの目標は音速抑制法を用いてアネラス
ティック近似の研究を再現すること
• 2年後の次世代スパコン完成時に超並列化され
た大規模計算をおこない太陽内部を理解する
おまけのページ
Super artificial viscosity
The diffusive flux is determined as follows
F
i
1
2
1
 c 1 (ur  ul , ui 1  ui )(ur  ul )
2 i 2
2
 (ur  ul ) 

 for (ur  ul )  (ui 1  ui )  0
 (ui 1  ui ) 
otherwise  0
セル境界での値の決め方
ul
ur
i
i+1
セル境界での値の決め方
ul
i
ur
i+1
セル境界での値の決め方
ul
ur
i
i+1
セル境界での値の決め方
ul
ur
i
i+1
ratio of speed of sound and convection
v  Tds
ds  c p
2
T  c / cp
2
s
v  cs 
太陽のpモードが固有振動を
作る。太陽表面の振動を
長い時間観測すると
右図のような図を得ることが
できる。
太陽標準モデルからのどのような
擾乱によって、この図を得ることが
できるか流体の逆問題を解くことに
よって求める。
周波数
What is the helioseismology?
水平方向の波数
inhomogeneous ξ
• 不均一なξは使えないと今では考えている。
 0 v 
0     2 
 
1
1
  2    0 v 
t

• 圧力勾配を変えるのは、ダイナミクスを変え
すぎる。
anelastic approximation
0 v      (Wez )    ( Lez )
anelastic 近似により、速度場が上記のようにポテンシャルで
表せる。これをフーリエ変換により展開して、W、Lの時間発展を
解く。
フーリエ変換時に毎回グローバルの通信が必要になり、
大量のCPUを使う計算では、効率が悪くなる。

Download Report