もっと大規模な計算が必要なのか？�
ニュートリノ輻射輸送、計算科学との連携�

2�

重力崩壊型超新星爆発	

~20Msun	

重力崩壊	

鉄のコア
Fecore	ν� ρc ~1010 g/cm3	� ρc ~1012 g/cm3	

Tc ~1 MeV	

Yec ~0.42 	

電子捕獲	

反応	

ν-閉込め	

Tc ~2 MeV	

Yec ~0.34	

ν
ν
ν
ν反応	

超新星ν
	

ν

1000 kmν

コアバウンス	爆発	� ρc ~3x1014 g/cm3	ν

Tc ~5 MeV	

Yec ~0.31	

中性NS子星	νν νShockwave
�衝撃波�

�

ρc ~4x1014 g/cm3	

Tc ~10 MeV	

Yec ~0.31	

元素合成�

km~ 10	

中性子星	
T~0 MeV, Ye<0.1	

多次元性が重要�
ニュートリノ加熱との組み合わせ�

4�数値シミュレーションによる解明�

古くは	1966年から	• スーパーコンピュータの発展と共に�– ニュートリノ輻射流体計算：近似から厳密へ�
– ニュートリノ・核物理データ：より高い信頼度へ�– 空間次元：球対称から多次元へ������新事実が明らかになっていく��

2000年頃	

１次元: 球対称	

図形：Wikipediaより	

2005年以降	

２次元：軸対称	

2D, 3D超新星: 滝脇氏	ごく最近	

３次元：対称性仮定無し	

5�

多次元計算：非対称性の役割が明らかに�

SN1987A	

対流, 回転, 非等方ν放射, 磁場,... 	

観測による傍証	

流体不安定性とニュートリノ加熱	

Wang (2002)	

ν	ν	

Proto-	NS	

shockwave	

ν-heating�

多次元でのニュートリノ加熱�
Suwa et al. ApJ (2013)	

6�

多次元計算によるEOS影響の研究	 Suwa et al. ApJ (2012)	

Shen EOS vs LS-EOSによる比較 (K, Asymが違う）	

K=180
Asym=29	15Msun	

11Msun	K=290 MeV
Asym=36 MeV	

2D(軸対称）で爆発している例�
バウンス後~500ms経ってから爆発へ�
Acoustic Powered	

Burrows et al. ApJ 640 (2006) 878	

拡散近似	Shen-EOS	SASI + Neutrino-heating	Marek & Janka ApJ 694 (2009)	

ray-by-ray法	

LS-EOS	

課題：異なる計算手法、状態方程式の違い、エネルギーが低い�
8�

3D(対称性無し）での計算が可能になってきた�• ２次元には無い流体不安定性→ニュートリノ加熱に有利？	
• ニュートリノ輻射輸送は近似的扱い	
larger 3D simulations running �@ K-computer by Takiwaki et al.�

太陽質量11倍の星	爆発している！！		

まだエネルギー不足 0.3ベーテ		

もっと重い星の場合は？	

9�

ニュートリノ輻射輸送が鍵�
加熱量：ニュートリノから物質へ�10�

ν-加熱は中間領域で起こる：拡散と自由伝搬の間	ニュートリノ輻射の詳細が必要：� エネルギー・角度分布	 f (Eν ,θν )ex. 拡散近似では不十分	

フェルミ分布は使えない	

freestreaming	 ν	• ν-反応: エネルギー依存性が強い	σ(Eν ) ∝ Eν2outside	R	€shockwave	

核子・原子核・レプトン反応	• ν-加熱率：輻射輸送＋ν反応による�€

heating region	

平均エネルギー・光度: Eν, Lν	

Flux factor:	θν	< µ >=< cosθν >= 0 ~ 12%(MeVLEiν νQν ≈ 110
Xi*' 2s⋅ N & R7 < µ > )€ν-sphere	

Janka A&A (1996)	

€diffusion	ν	ν	center	
11�

爆発の可否への影響は大きい�
10~20 %の加熱率変動でも爆発の可否に影響がありえる	• ν加熱の依存性	

核子吸収による加熱を30, 50, 70%	増加させてみる (from tpb=100 msec)	

30%	

70%	

50%	

• ν光度Lνの3~12%の変化�
�爆発に転ずる(1D,2D)	

• threshold付近�

Eexp	

Exp.	

1051 erg	

No Exp.	

time	Sumiyoshi	

Janka Müller A&A(1996)	

• ニュートリノ輻射輸送による緻密な計算が必要�
12�ν輻射輸送：近似・厳密計算の変遷	
• 厳密計算手法�

1D	– Boltzmann equation	

– Moment formalism, Monte-Carlo	

1Dにおけるニュートリノ放出予測・ミクロ物理の検証	

• 近似計算手法の変遷�

2D	– Light bulb, Leakage scheme	

�放出表面を固定, 放出時間を考慮	

– Diffusion, IDSA, Ray-by-ray	�拡散近似, 1Dを拡大適用�
2D, 3Dにおける爆発メカニズム探索�

3D	

�
13�

ν輻射輸送：近似から厳密計算へ	• 2D, 3D: 近似計算手法�

Ray-by-ray method	

Suwa et al., PASJ (2010)	

– state-of-the-art simulations	– Diffusion (flux limited)	

�中心付近のみ	

– Ray-by-ray (radial transport)	

�動径方向のみ�

	 2D Sn-method	Ott et al. (2008) ApJ	
1D-transport independently	
• 3次元ニュートリノ輻射輸送計算コード�
�- 近似の影響を明らかにする� Sumiyoshi & Yamada, ApJS (2012)	

�- 厳密計算へ置き換える → Grand Challenge�
14�

6次元ボルツマン方程式の解法�
大規模計算で解くための工夫�

15�

3次元輻射輸送計算	

θν	z	• 6次元ニュートリノ分布 (3D+3D)	

θ			fν (r,θ, φ; εν ,θν , φν ; t)€1 ∂fν  
1 ' δfν *+ n ⋅ ∇fν = ) ,c ∂tc ( δt + collisionφν	r	– エネルギー (εν), 伝搬方向 (θν, φν)	• 6次元でのν輻射輸送方程式�

ν	y	x	φ	Boltzmann eq.	

– 左辺：ニュートリノ数の変動�
– 右辺：ニュートリノ反応による変動（衝突項）�

€ν		16�ニュートリノ輻射輸送方程式�
1− µν2 cos φν ∂1− µν2 sin φν ∂fν
1 ∂fν µν ∂ 2+ 2 (r f ν ) +(sin θf ν ) +c ∂t r ∂rr sin θ∂θr sin θ∂φ21−µ1 ∂1 ( δf ν +2ν cos θ ∂+[(1− µν ) f ν ] +(sin φν fν ) = * r ∂µνr sin θ∂φνc ) δt , collision�– 保存形：差分化して解く�(Sn法)�
– 時間に関して陰解法�

€Eν,µν	E´ν,µν´	ν	µν = cosθν€
1 $ δf ν '1= j emission (1− f ν ) −f ν + Cscattering [ ∫ f ν ( E ν- , µν- )dE ν- ]& )c % δt ( collisionλabsoption€– 衝突項：全てのν吸収・放出・散乱過程�
€• 散乱前後のエネルギー・角度ごとに積分計算�
• 状態方程式に依る組成・ニュートリノ反応率、相対論�
Sumiyoshi & Yamada, ApJS (2012)	

17�

衝突項：ニュートリノ反応による寄与�
• 爆発メカニズム解明に必要な基本セット	

Bruenn (1985) 	

+ Shen	• 放出/吸収:�� e- + p ↔ νe + n ��e- + A ↔ νe + A´�� e+ + n ↔ νe + p	

• 散乱:�
� νi + N ↔ νi + N ��νi + A ↔ νi + A	�	• 対生成・消滅:�
�e- + e+ ↔ νi + νi ��
� 3 species:	

���N + N ↔ N + N + νi +νi ��ν , ν , ν 	e•eµ

ここ数年の計算資源で扱うための処方：	

– エネルギーが変化する反応は扱わない（電子散乱）	

– 線形方程式の範囲で取り入れる（対生成・消滅）	

– ローレンツ変換の効果を取り入れる（作業中）	
18�

主な計算負荷：大規模疎行列の解法�
• 線形方程式�Af ν = d• ニュートリノ分布 �€~109	

Nspace=nr x nθ x nφ �Nν=nε x nθν x nφν �Nvector~106 x 103 �
ν分布： 12GB�
行列: 1.5TB�

• 反復法�
BiCGSTAB�
前処理 by 今倉�Imakura et al. JSIAM (2012)	

Kotake et al. PTEP (2012)	

19�

3次元ニュートリノ輻射輸送計算コードの検証�
Sumiyoshi & Yamada ApJS (2012)	

• 解析解との比較	

– ガウス分布拡散、自由伝搬	

– 経路に沿って形式解	

– 平衡解への時間発展	

diffusion	

equlibrium	• 球対称計算コードとの比較	

– 固定分布	• 重力崩壊, バウンス後	– ν密度, フラックス, モーメント	– ν反応率, 平均自由行程	

1D GR ν-radiation hydrodynamicsコードを活用	

Sumiyoshi et al. ApJ (2005)	
20�

2D/3D超新星コアへの応用�
多次元でのν輻射輸送の振舞は？�

21�

3次元ν輻射輸送計算コードの応用�
• MPI並列計算コード	

– KEK, YITP, UT: SR16000, 8ノード	 (256 x 64 x 32, 6 x 12 x 14)	• 物質分布のもとでν輻射輸送を計算できる	– ニュートリノ輻射の定常分布を調べる	

• 6Dボルツマン方程式により時間発展	

• 定常解による流束, 加熱率, モーメント量	• 2D/3D超新星コアでのν輻射の特徴	

– 非動径方向の移流, 近似計算の吟味	• 流体計算コードとの結合（長倉）	

– 2D ν輻射流体コード, 相対論化 → 3Dフル計算へ	22�

３次元超新星コアにおけるニュートリノ輻射輸送	
3D profile by Takiwaki (2012)	

3D（近似）ν輻射流体計算による物質分布	

3D超新星コアでのニュートリノ輻射を調べる	

ニュートリノ密度, フラックス, スペクトル，�加熱率	

ν-heating/cooling rates in 3D	

11.2Msun, 3D	200ms	

νe	

density iso-surface	

Red: heating, Blue: cooling	23�

3次元超新星コアの物質分布�
流体分布 (ρ, T, Ye) 固定, 3D ν輻射輸送を解いて定常解を求める	

entropy	

density	Ye	

entropy	
11.2Msun, 3D	

150 msec after bounce	Rshock~250-400km�

Takiwaki et al. ApJ (2012)	

24�

€

近似ν輻射輸送との比較�
• Full-Boltzmann 6D: �1− µν2 cos φν ∂1− µν2 sin φν ∂fν
1 ∂fν µν ∂ 2+ 2 (r f ν ) +(sin θf ν ) +c ∂t r ∂rr sin θ∂θr sin θ∂φ21−µ1 ∂1 ( δf ν +2ν cos θ ∂+[(1− µν ) f ν ] +(sin φν fν ) = * r ∂µνr sin θ∂φνc ) δt , collision• Ray-by-ray approximation: � 非動径方向の移流をSwitch OFF	
1− µν2 cos φν ∂1− µν2 sin φν ∂fν
1 ∂€fν µν ∂ 2+ 2 (r f ν ) +(sin θf ν ) +c ∂t r ∂rr sin θ∂θr sin θ∂φ21−µ1 ∂1 ( δf ν +2ν cos θ ∂+[(1− µν ) f ν ] +(sin φν fν ) = * r ∂µνr sin θ∂φνc ) δt , collision
Sumiyoshi et al. (2013)	

€25�

3次元超新星でのν輻射輸送: 近似の場合�
Z軸上から見た図	

νe density, flux	

Ray-by-ray: 	

θ, φ-advection OFF	r方向の流束のみ	

Sumiyoshi et al. (2013)	

3D Supernova core 150msec after bounce	26�

3次元超新星でのν輻射輸送: Fullの場合�
Z軸上から見た図	

νe density, flux	

6D-Boltzmann�
非動径方向の	移流も記述できる		φ方向の流束	

Sumiyoshi et al. (2013)	

3D Supernova core 150msec after bounce	

27�

近似計算�との比較�

Ray-by-ray: radial only�

Full-Boltzmann: 6D�

• Ray-by-rayは		

• 空間異方性を	強く出す傾向		

• 非動径方向の	移流を記述できない			νe density, flux	θ方向の流束	

Z軸横から見た図	
Sumiyoshi et al. (2013)	

150msec 	

28�

非動径方向フラックスの寄与：Fth/phi / Fradial�•Full-Boltzmann: 6D�

• Theta direction�

Z	Z	

• Azimuthal direction�

Fr	Fr	

◉	×	Fφ	○	Fθ	R	

Sumiyoshi et al. (2013)	

R	29�ν加熱率の違い�

Qrbr − Q3D
δ=Q3D• Ray-by-ray: radial only� • Full-Boltzmann: 6D� • Deviation�€Sumiyoshi et al. (2013)	

30�

物質分布の角度分布における非一様性�
density	ν-sphere	

Sumiyoshi et al. (2013)	

degeneracy for neutrinos	ν-sphere	Ray-by-rayではLocal fluctuation (hot spot)が強調される	

31�

We need Exa-scale computing�Scaling up to large scales�

32�Scaling of memory and operations: Matrix�
Nν	Nν x Nspace	• Block dense matrix	

for ν-angle changes 		

in collision term	• Block size		

M=nθν x nφν 	M x M	• Memory size		

M2 x nε x Nspace	• Operations		

M3 x nε x Nspace		

Kotake et al. PTEP (2012)	
Estimate by assuming:		

iso-energetic reactions & pre-condition by inversion of blocks	
33�

Necessary computational resources�
• To run 6D Boltzmann solver with high resolution	

Kotake et al. PTEP (2012)	

machine�space�

neutrino�

memory� operations�Current�

256 x 32 x 64�

8 x 12 x 14�

2TB�

K-computer� 512 x 64 x 128�Exascale�

12 x 24 x 20� 200TB�

512 x 128 x 256� 24 x 24 x 24� 3PB�

6 Tera /step�
2 Peta /step�
80 Peta /step�- For 6D Boltzmann + Hydrodynamics with iteration		

Increase by a factor of ten more…	

- Need Exa-flops machines for full 6D simulations	34�まとめと展望�• 超新星爆発メカニズムにおける核物理	– 状態方程式の影響：固さ・柔らかさ, 組成	– ニュートリノ放出量・加熱による影響	

– 超新星ニュートリノ観測で状態方程式を探る	

• 信頼できる原子核データ群の整備が重要	

• 2/3次元超新星での爆発メカニズムの課題	

– ニュートリノ輻射輸送計算の重要性	

– 6Dボルツマン方程式を解く計算コード	

• 3次元ν輻射流体計算コードへ	– エクサスケールへ向けたグランドチャレンジ	35�

原子核・天文・計算科学の共同研究�
• 超新星爆発�– 鈴木英之、中里健一郎�
– 山田章一�– 長倉洋樹、岩上わかな�

• 計算科学�– 松古栄夫�
– 櫻井鉄也、今倉暁�• 天文シミュレーション�
– 固武慶、滝脇知也�
– 関口雄一郎�• 状態方程式テーブル�– H. Shen�– 親松和浩�– 土岐博�

• 状態方程式の拡張�

– 古澤峻�
– 大西明、石塚知香子�
– 鷹野正利、富樫甫�

• ニュートリノ反応�
– 佐藤透、那須翔太、中村聡�

�
新学術領域「計算科学による素核宇融合」	

HPCI戦略P分野５「物質と宇宙の起源と構造」課題３	
Supercomputing resources at KEK, YITP, UT, RCNP, K-Computer	Support from: 	

- Grant-in-Aid for Scientific Research (20105004, 20105005, 22540296, 24244036)	

36�