富岳における太陽内部計算
Hotta & Kusano, 2021, Nature Astronomy accepted堀田英之(千葉大学)
共同研究者：草野完也(名古屋大学)
1

太陽

黒点の面積

NAOJ/JAXA/NASA年

2太陽内部構造

対流層

L = 3.84 ⇥ 1033 erg s 1
11年周期で0.1%の変動あり
輻射エネルギーフラックス

放射層

対流層
計算領域
内部構造(密度、圧力、温度...)などは理論観測により、非常に精密に決定している。
運ぶべきエネルギー量も精密に決定。
初期条件・パラメタの不定性なし。「正しく計算」すれば「正しい答え」が出るはず。

輻射エネルギーフラックス(よく知られた量)3

大規模流れ観測がしっかりしているので、11年周期だけを再現できてもだめ差動回転の観測(R. Howeデータ)

子午面環流
Zhao+2013
SDO

Gizon+2020
SOHO, GONG

差動回転(~2 km/s)の構造はほぼ確定しているが、子午面環流(~10 m/s)は観測装置によって結果が違う(原因不明)→少なくとも差動回転は正しく再現する必要あり
4太陽差動回転理論の問題
太陽自転速度では、赤道を速くできない：Convective conundrum
2倍太陽自転1倍太陽自転格子点数各方向にそれぞれ4倍3倍太陽自転

高解像度化(=現実に近づける)は問題を悪くするセンス5

研究の目的富岳を用いた高解像度計算で熱対流の難問(Convective conundrum)を解決する。
つまり高解像度化で、赤道が速く自転するという実際の太陽の差動回転分布を目指す。
どこまで解像度を上げればいいのかわからないが、富岳でできるところまでを目指す。

6

太陽内部計算で解く方程式系
磁気流体力学の方程式
Magnetohydrodynamics (MHD)@⇢@t=@(⇢v)@t=

ローレンツ力@s
⇢T
@t
@B
@tp
<latexit sha1_base64="ZSr6xrVUPfoDnE9wHLayDckGblc=">AAAEjHichVHLbtNAFL0NAUqAPmCDxGZE2iihtJqEViBQpSpIiB19pa3UqaKxM02t+sV4EiiWf4AfYMEKJBaIz2DDig1i0U9ALIvEhgV3xqbNo5SxbN/XOffcO1boOpGi9HAkdy5//sLF0UuFy1eujo1PTF7biIKOtEXDDtxAblk8Eq7ji4ZylCu2Qim4Z7li09p/pPObXSEjJ/DX1UEodjze9p1dx+YKQ83JkS/MEm3Hj8Uzn0vJD24nBbYruR2zkEvlcJcwuRckJ65KSGmxRGaZzy2XM7sVqLIuIczy4m5SIYwNMPSBTe3f0lOJTFJnZ0maIiGZIcct2iiQlUoYMU2qSTzPQgd5UxrleCLSdfWEVPq8GVLr0XmSYk890eaJlq3T6wPjR6fNbgpJOaPS0tP2FRJhn5Um87jak14seSsZ2keqZ5i2d4LyoMp6utnQVIZmUXdIVCkw4bd6Lq85UaRz1BwybFQzowjZWQ4mvgKDFgRgQwc8EOCDQtsFDhE+21AFCiHGdiDGmETLMXkBCRQQ28EqgRUco/v4baO3nUV99DVnZNA2dnHxlYgkME2/0Q/0iH6mH+l3+vufXLHh0FoO8G+lWBE2x1/dWPv1X5SHfwV7J6gzNSvYhftGq4PaQxPRU9gpvvvy9dHag9XpuETf0R+o/y09pJ9wAr/7036/IlbfGHYfMc/NtJ7p7+N+Y4xb6CVmo7qDjuh5AniRZfV2WseKCUxBEWunEIMPXmx18BqHjY3aXPXu3MJKrbhUz654FG7CLSgj0z1YgiewDA2wcwu57VwrJ/Jj+fn8w/xiWpobyTDXoe/kH/8BXI1Bcw==</latexit>=r · (⇢v)

質量保存則

r · (⇢vv)

rp + ⇢g

重力運動方程式

1+ (r ⇥ B) ⇥ B + 2⇢v ⇥ ⌦4⇡⇢T (v · r)s + Qrad
輻射

コリオリ力

エントロピーの方程式
or エネルギー保存則=r ⇥ (v ⇥ B)=p(⇢, s) 状態方程式(電離の効果を含めるOPAL)

磁場の誘導方程式

これらの方程式を格子法で解く
7

コード開発R2D2コード(Hotta+2019, Hotta+2020)を富岳用にチューニング
京では実行効率20%超程度。
音速抑制法(Hotta+2012)、アルフベン速度抑制(Rempel+2009)などを使って全体通信のないアルゴリズム
京で使っていたコードをそのまま富岳で実行すると実行効率4.5%程度。ある程度のチューニングが必要と考えられる。

8富岳の特徴
値はノードあたり
性能 京：128 GFLOPS/node、A64FX：3 TFLOPS/node
L1$

L1$ BW

L1$ B/F

L2$

L2$ BW

L2$ B/F

京

256 KB

512 GB/s

45 MB

256 GB/s

2

A64FX

3 MB

12 TB/s

432 MB

4 TB/s1.2

また、A64FXは京に比べてレジスタが少なくキャッシュのレイテンシも多い。(数値は公式資料で確認できなかったが、牧野さん資料参照)

9富岳のチューニング方針
1. キャッシュのレイテンシを隠蔽するソフトウェアパイプラインを促進する→レジスタ少ない→最内ループの中の演算量を減らす→ループ分割
2. 低いキャッシュ幅をごまかすためにレジスタを有効に使う→最内ループになるべく多くの物理量を登場させる→ループ内の演算を増やす

10富岳におけるMHD計算の困難
1. 変数が多く、色々な演算に登場
密度、内部エネルギー、速度三成分、磁場三成分、圧力＝8+1成分
ループ分割をするとデータ非局所化の影響が大きい
2. 解像度が重要なので、ステンシルを大きくして多くの要素を参照するがある。R2D2コードでは空間微分は４次の中央差分を利用

L2キャッシュにも3次元データは乗らないのでナイーブにコードを書くと、一度読み込んだデータを有効に活用できない

11

ループ分割とデータ局在化両立のアイディア

微分の各要素を一つのループにまとめる。→結構速くなる12

ループ分割とデータ局在化
ナイーブにループ分割しても、コードが遅くなることがほとんどだった。
MHDには色々な演算があるので、分割・融合を繰り返して、最適な演算順序を考える。

13

ループ一重化

12 openMP並列が基本なので、最外ループの並列化のみでは足りないことが多い。
ループ一重化で並列化効率をアップ。

14

配列の構造(旧バージョンでの計測)
配列の要素の並び順は、チューニングの大きな要素。
四次元配列については全て試した(m,i,j,k)(i,m,j,k)(i,j,m,k)(i,j,k,m)

実行効率

7.3%

8.2%

8.8%

8.9%

実行時間

52.4秒

46.3秒

43.8秒

43.2秒

(i,j,k,m)が一番速かった。京の時は、(m,i,j,k)だった。この前の説明会で聞いた話によると京の時は、このようにしても速くする仕組みがあったとのこと
説明会では(i,j,m,k)が一番速いと言っていた
二次元の一時配列についても(m,i)とするより(i,m)とした方が早い

15

チューニング結果
MPI Processあたり128 64 512 100 step
チューニング前(実行時間)チューニング前(実行効率)チューニング後(実行時間)チューニング後(実行効率)

全体110.1 秒

4.57%

44.5 秒

9.3%

Runge-Kutta
ループ

32.9 秒

4.87%

18.9 秒8.2%

人工粘性68.0 秒

4.26%

19.9 秒

10.5%

輻射輸送

5.27 秒

7.68%

4.2 秒

8.4%

3.8 107 cell update/sec/nodeを達成
4 step Runge-Kutta使っていてCFL numberが2くらいまでいけるので、時間2次精度のMHDスキームと比較するとさらに2倍くらい速い(と思ってもらって良い)
16

MPIチューニング
TofuトポロジーでのMPIスレッドの配置を最適化する。
#PJM –L ”node=3x4x2"
include "mpif-ext.h"call fjmpi_topology_get_dimension(fjmpi_dim,merr)call fjmpi_topology_get_shape(ix0,jx0,kx0,merr)
ノードがどこにあるかしか確かめることができない。
CMGがどこにあるかわわからない。
富岳では1ノードにつき4MPIスレッド使うので、y方向に展開する。
Yin-Yang gridでは、3次元構造をy方向に二つに割って、YinとYangを割り当てているこのようにすれば少なくとも2次元内で相手を探せるので際立った性能の劣化はない(最適ではない)
17

チューニング結果
MPI Processあたり128 64 512 100 step
チューニング前(実行時間)チューニング前(実行効率)チューニング後(実行時間)チューニング後(実行効率)

全体110.1 秒

4.57%

44.5 秒

9.3%

Runge-Kutta
ループ

32.9 秒

4.87%

18.9 秒8.2%

人工粘性68.0 秒

4.26%

19.9 秒

10.5%

輻射輸送

5.27 秒

7.68%

4.2 秒

8.4%

3.8 107 cell update/sec/nodeを達成
4 step Runge-Kutta使っていてCFL numberが2くらいまでいけるので、時間2次精度のMHDスキームと比較するとさらに2倍くらい速い(と思ってもらって良い)18

スケーリング

100万超えてもweak scalingならば、99%を超える並列化効率を示している。大規模計算も滞りなくできそう
19細かい富岳利用Tips
Ø 言語関連マニュアル・言語使用手引書をよく読む
Ø 不老などと一緒に使うのが吉ü 384node=18432core以下混んでる？
Ø on the fly解析推奨
ü 大規模計算は状況を知るだけでも一苦労Ø ストレージ圧倒的に足りない。1ケース100TBくらい利用

20富岳でおこなった計算設定
計算領域：0.71𝑅⊙ < 𝑟 < 0.96𝑅⊙
格子点数：
Low

96 768 1536

Middle192 1536 3072

High384 3072 6144

計算時間：4000日
Highの場合350万ステップほど
3000日ほどでほぼ統計的定常に遷移
3600-4000日でのデータを平均して結果を示す
𝑟 = 0.9𝑅⊙での動径方向速度21熱対流と磁場
Normalized entropy

Magnetic field strength𝑠 − 𝑠 /𝑠"#$

𝐵 [kG]

22

差動回転の角速度依存性
解像度アップ
96 768 1536

192 1536 3072

384 3072 6144

Low caseでは、これまでの計算と同じく太陽回転角速度だと赤道よりも極が速い差動回転が実現してしまっているが、High caseではほぼ太陽と同じ差動回転を実現することができている。
23

熱対流速度とダイナモ生成磁場

解像度が高くなると、熱対流速度は遅くなり、磁場は強くなる。
Low caseの時は、対流層内部の全ての場所で𝐸%&' > 𝐸#()だが、High caseの場合は対流層内部全てで𝐸%&' < 𝐸#()となり状況が全く変わっている。
詳しくは説明しないが、磁場が角運動量を運ぶことがエッセンス

24

さらに高解像度へ
Hotta & Kusano

新しい計算

子午面環流

差動回転

実際の太陽太陽表面近くの構造を再現できていない。
(子午面環流は深い場所では感度なし)

𝑟!"# = 0.96𝑅⊙ , 384×3072×6144
350万時間ステップ→𝑟!"# = 0.99𝑅⊙ , 512×4096×8192
1200万時間ステップ
25

まとめØ 太陽内部の構造はよくわかっているが、数値計算では太陽の差動回転を再現できないという「熱対流の難問」が長らく課題となっていた。
Ø 富岳を用いた大規模計算で解決を試みた
Ø MHD計算を富岳で行うとき、複雑な方程式系・多い変数種類が要因で、ループ分割・データ局所化が困難
Ø 色々な試みの結果、チューニングはある程度の成果をみせ、富岳による大規模計算が可能にØ これまでにない大規模計算で赤道が速い差動回転を実現。「熱対流の難問」解決
Ø 高解像度化はまだまだ価値があり、高解像度化によりさらに観測される太陽の特徴を再現できるフェーズに入ってきている。
26