2022 年第 52 回 天文・天体物理若手夏の学校

星・惑星形成分科会

集録集

■ 謝辞
2022 年度天文・天体物理若手夏の学校は、基礎物理学研究所を始め、国立天文台、光学赤外線
天文連絡会、理論天文学宇宙物理学懇談会、天文教育普及研究会、日本天文学会からのご支援により成り立っております。事務局一同厚く御礼申し上げます。
星・惑星形成分科会2

■ index星惑 a01

玉城 磨生

大マゼラン雲における O 型星アソシエーションと HI 及び CO 分布の比較

星惑 a02

野﨑 信吾

分子雲コアの周辺環境の違いによる質量降着率と星形成効率の変化

星惑 a03

小道 雄斗

星間衝撃波により形成されるフィラメント状分子雲の化学進化星惑 a04

深谷 直史

Corona Australis 分子雲 IRS7 領域における分子雲コアの性質と星団形成星惑 a05

吉田 有宏

分子ガス同位体比の新しい測定手法を用いた TW Hya 円盤における不均一な12

CO/13 CO 比の発見

星惑 a06

安田 郁斗

木星の氷微惑星散乱による雪線内側への水供給についての理論的研究

星惑 a07

宮山 隆志

微惑星衝突に伴う固体蒸発量とその状態方程式依存性星惑 a08

講演キャンセル星惑 a09神原 祐樹

近接遭遇による惑星系の軌道進化と AMD を用いた安定性評価

星惑 a10

河合 優悟

トランジット測光による系外惑星 KELT-9 b のスピン傾斜角の推定星惑 a11

上内 雅暁

SETI につながる月面反射人工電波観測星惑 b01

下川 貴史

原始惑星系円盤の外側の領域における弱い乱流の起源

星惑 b02

瀬野 泉美SPH 法の不連続面における密度計算法の開発

星惑 b03谷安 要

高汎用性 N 体計算コード REBOUND 及びその使用例

星惑 b04植野 雅々

太陽系外惑星 TrES-1b の多波長トランジット観測星惑 b05平野 佑弥

太陽系外惑星の天気の調査

星惑 c01

所司 歩夢

スパースモデリング (SpM) による超解像度画像解析

星惑 c02

赤池 祐貴

若い星状天体の円盤、アウトフロー、磁場の初期進化とダストモデルの違いによる影響

3

——–index へ戻る星惑 a01

大マゼラン雲における O 型星アソシエーションと HI
分布の比較
名古屋大学大学院 理学研究科
玉城 磨生

4

2022 年度 第 52 回 天文・天体物理若手夏の学校

大マゼラン雲における O 型星アソシエーションとHi 分布の比較
玉城 磨生 (名古屋大学大学院 理学研究科)

Abstract
大質量星，特に質量 ≳20M⊙ の O 型星は，紫外線放射や恒星風などにより周囲の星間物質に多大な影響を及ぼす。しかし，その形成機構は明らかではない。大マゼラン雲（LMC）では，小マゼラン雲との潮汐相互作用を起源とする Hi ガスの衝突が大質量星団の形成を誘発したとのシナリオが提案されてい。そこで我々は，Zaritsky et al. (2004) による LMC 全面での UBV 測光データを元に O 型星候補天体の同定を試み，周囲の星間ガスの性質との関係を調べた。Pecaut et al. (2013) の O 型主系列星のカラーのテーブルを参照し，Furuta et al. (2022) で得られた近赤外減光マップから星間赤化を考慮した結果，O 型星候補天体を 1109 個同定した。これら O 型星候補天体の分布は不均一で，グループを形成して存在している傾向が見られた。ここで，それぞれのグループを星団と考え，O 型星候補天体の面密度画像に対して Dendrogram を用いること
で 73 個の星団候補を同定した。1 つの星団は平均で 17 pc の半径を持ち，約 8.3 個の O 型星候補天体が含まれていた。同定された星団候補は Hi や CO の強度が強い領域に多く分布する傾向があり，stellar bar にはほとんど存在しなかった。また，衝突により減速されたと考えられる Hi の中間速度成分 (I 成分) の分布とは良い一致を示した。また，星団候補を構成する O 型星候補天体の数とその周辺の Hi の柱密度，速度の関係から衝突によって生まれる星団の規模は初期条件に依存することを示唆する結果を得た。

1

我々は，Zaritsky et al. (2004) による LMC 全面で

Introductionの UBV 測光データを元に O 型星候補天体の同定を大質量星は超新星爆発，紫外線放射，星風によっ試みた。

て物理的にも化学的にも銀河進化において大きな役割を果たしており，その形成過程の解明は現代の天文学における最も重要な課題の一つであると言える。 2
大マゼラン雲（LMC），小マゼラン雲との潮汐相互Methods/Data

作用を起源とする Hi の衝突が大質量星団の形成を誘本研究では，Pecaut et al. (2013) の O 型主系列星発したとのシナリオが提案されている (Fukui et al.

のカラーのテーブルを参照し，Furuta et al. (2022)

2017; Tsuge et al. 2019)。LMC には潮汐駆動流（L で得られた近赤外減光マップから星間赤化を考慮し，成分），ディスク（D 成分），L 成分と D 成分の間 O 型星候補天体の同定を行なった。得られた O 型星にあり衝突圧縮ガスであると考えられる中間速度成候補天体の位置，測光データ，減光量，質量，絶対等

分（I 成分）の Hi ガス (Tsuge et al. 2020) が存在す

級，見かけの等級をまとめたカタログデータを作成るとした結果も得られた。そのシナリオはどの程度した。分解能である 24 pc pixel−1 の O 型星候補天体

一般性を持つのかは明らかになっておらず，更に研の面密度画像を用いて星団候補の同定を行なった。星究を深めることが重要である。しかし，LMC におけ団候補の同定には，Dendrogram (Rosolowsky et al.

る大質量星のサーベイは限定的である。これまでの2008) を使用した。バックグラウンドの 2σ を星団候

最も包括的なカタログ (Bonanos et al. 2009) は，限補検出の閾値とし，1σ の最小有意性，2 pixel 以上のられた領域に対する分光観測をまとめたもので，不広がりを持つ構造を星団候補として定義した。また，

均一なサンプルであり統計的な研究には不適である。 星団候補の物理量として，ピーク位置，重心の位置，面積，半径，含まれる O 型星候補天体の数を求めた。

52022 年度 第 52 回 天文・天体物理若手夏の学校
ここで得た O 型星候補天体の物理量と LMC におけるオーストラリアコンパクト電波干渉計 (ATCA) と

1 kpc

-66

パークス天文台 (Kim et al. 2003) で得られた LMC
全体の Hi 21cm 線のアーカイブデータを使用し，L

Declination

成分，D 成分，I 成分の Hi データ (Tsuge et al. 2020)を用いた統計的解析を行なった。観測の詳細については Kim et al. (2003) を参照。

-68

-70

3

Results
-65

-72

6:00

Declination

-67

5:00

図 2: 同定された星団候補の分布。背景の積分強度-68

図は Hi の I 成分，灰色のコントアは Hi の I 成分柱
密度 0.5×1021 cm−2 を示す。青色の円は 10 個以上-69の O 型星候補天体から形成される星団候補，赤色の-70

円は 10 個以下の O 型星候補天体から形成される星団候補であり，円の大きさは星団半径を示す。

-71
-72

5:30

Right ascension

-66

6:00

5:30

を占め，50 個以上の O 型星候補天体を含む星団候補

5:00も存在することがわかる。星団半径（r ）は 10–20 pc

Right ascensionの星団候補が約 90%を占める。

図 1: 同定された O 型星候補天体マップ。背景の積分強度図は Hi，青色のクロスが O 型星候補天体を4表す。

Discussion
LMC において母体となる星間ガスの分布と星団の

図 1 は UBV 測光データから同定された 1109 個規模の関係について議論する。星団候補を星団中心の O 型星候補天体の分布を示す。これらの O 型星候から r の円内領域と r から 2r の円環領域に分け，Hi

補天体は不均一で集団的に分布し，グループを形成との比較を行なった。図 4(a)，(b)，(c) は円内領域，

して存在している傾向が見られた。これらは，星団円環領域における Hi のそれぞれ L 成分，D 成分，I

（アソシエーション）を形成していると言える。

成分の柱密度のピークと O 型星候補天体の数の散布

図 2 は，図 1 をもとに Dendrogtam で同定した 73図である。黒色で示した円内領域が，白色で示した個の星団候補のうち楕円近似可能であった 66 個の分円環領域よりグラフの左側にプロットされている傾布を示す。これらは Hi や CO の強度が強い領域に多向があるため，円内領域では円環領域より Hi 柱密度く分布する傾向があり，stellar bar にはほとんど存が減少している領域が多いことがわかる。これによ在しなかった。また，衝突により減速されたと考え

り，円内領域では O 型星候補天体の輻射による Hi のられる Hi の I 成分の分布とは良い一致を示した。

電離などが進んでいると考えられる。つまり円環領
図 3 は星団候補の典型的な物理量を示している。 域では，より星団形成初期の物理状態を示している

2–10 個の O 型星候補天体を含む星団候補が約 90% と考えられる。さらに円内領域，円環領域ともに D

62022 年度 第 52 回 天文・天体物理若手夏の学校
図 2 の青色の円で表される 10 個以上の O 型星候補
天体から構成される星団は 1 つの例外を除いて全て柱
密度 0.5×1021 cm−2 以上の Hi の I 成分が受かってい

freq

(a)
40

る。一方，赤色の円で表される 10 個以下の O 型星候35

補天体から構成される星団は Hi の衝突がまだ生じて30いないとされる Hi ridge 南部や stellar bar にも分布25しているが，その数は僅かである。これらより，より

20

大きな柱密度を持つ雲に −45–−35 km s−1 の適度な15相対速度で衝突することで O 型星を 10 個以上含むよ

10

うな大規模な星団形成が誘発される，0.5×1021 cm−25以上の柱密度を持つ中間速度成分が生じる領域で大規

0

0

10

20

30

40

Number of O-stars

50

模な星団形成が誘発される，大規模な星団形成は Hi

60の衝突による寄与と整合性があることが考えられる。

(b)
30

5

freq

25

Conclusion20

1. Zaritsky et al. (2004) による LMC 全面での

15

UBV 測光データを元に 1109 個の O 型星候補
天体を同定した。

10

2. O 型星候補天体の数密度マップに対して，Den-50drogram を使用することで 73 個の星団候補を同
0

5

10

15

20

radius [pc]

25

30

定した。

35

3. 同定された星団候補は Hi の I 成分の分布と良い一致を示した。
同定された星団候補の典型的物理量。(a) １

4. 今回同定した星団候補は約 90%が 2–10 個の O

つの星団候補に含まれる O 型星候補天体の数分布を型星候補天体を含む，もしくは半径 10–20 pc の示すヒストグラム。(b) 星団候補を楕円近似した際の星団候補であり 50 個以上の O 型星候補天体を

長軸短軸相乗平均の分布を示すヒストグラム。

含む星団候補も存在することが明らかになった。

図 3:5. 星団候補を構成する O 型星候補天体の数と Hi

成分において柱密度のピークと星団を形成する O 型の D 成分柱密度，I 成分柱密度には正の相関が星候補天体の数に正の相関が見られ，I 成分でもやや見られた。

正の相関を示すことが分かる。しかし，L 成分に関6. より大きな柱密度を持つ雲に −45–−35 km s−1

しては無相関であるといえる。
図 4(d) は円環領域内における Hi の D 成分に対の適度な相対速度で衝突することで O 型星をする L 成分の相対速度と星団候補を形成する O 型

10 個以上含むような大規模な星団形成が誘発さ星候補天体の数の散布図である。10 個以上の O 型れる。

星候補天体を含む星団候補の形成は相対速度 −45–

7. 0.5×1021 cm−2 以上の柱密度を持つ中間速度成−35 km s−1 の間で優位であり，衝突速度が遅すぎ分が生じる領域で大規模な星団形成が誘発される，あるいは速すぎる場合は大規模な星団形成は誘る。

発されないと分かる。

7

2022 年度 第 52 回 天文・天体物理若手夏の学校

(a)

60

60
50
Number of O-stars

Number of O-stars

50
40
30
20
10
0

0

0.5 1.0 1.5 2.0 2.5 3.0
HI peak column density [cm-2]

×1021

(c)

20

0

60

0.5 1.0 1.5 2.0 2.5 3.0
HI peak column density [cm-2]

3.5
×1021

(d)

50
Number of O-stars

Number of O-stars

30

0

3.5

50
40
30
20

40
30
20
10

10
0

40

10

60

0

(b)

0.5 1.0 1.5 2.0 2.5 3.0
HI peak column density [cm-2]

0
-75

3.5
×1021

-65
-45
-55
-35
Relative velocity [km s-1]

図 4: Hi の (a)L 成分，(b)D 成分，(c)I 成分の柱密度のピークと星団候補を形成する O 型星候補天体の数の散布図。黒色の点は円内領域，白色の点は円環領域を示す。(d) 円環領域内における Hi の D 成分に対す
る L 成分の相対速度と星団候補を形成する O 型星候補天体の数の散布図。

8. 大規模な星団形成は Hi の衝突による寄与と整合 Zaritsky et al., AJ, 128, 1606, 2004
性がある。
Furuta et al., PASJ, 74, 639, 2022
Tsuge et al, arXiv e-prints, 2020

Reference

Kim et al., ApJS, 148, 473, 2003
Fukui et al., PASJ, 269, 5, 2017
Tsuge et al., Apj, 871, 44, 2019
Pecaut et al., ApJS, 208, 9, 2013
Bonanos et al, AJ, 138, 1003, 2009
Rosolowsky et al., Apj, 679, 1338, 2008

8

——–index へ戻る星惑 a02分子雲コアの周辺環境の違いによる質量降着率と星形
成効率の変化九州大学大学院 理学府地球惑星科学専攻
野崎 信吾

9

2022 年度 第 52 回 天文・天体物理若手夏の学校

分子雲コアの周辺環境の違いによる質量降着率と星形成効率の変化
野崎信吾 (九州大学大学院 理学府地球惑星科学専攻)

Abstract
多くの星形成領域で IMF と CMF は類似しており、孤立した分子雲コアのガスの 30-70%が星に変換されると考えられている。一方、近年の観測から ONC 領域では IMF を CMF と一致させるためには CMF を100%以上増加させる必要があることが分かった。(Takemura et al. 2021a,b) これは分子雲コアがその進化過程で外部から質量を獲得する必要があることを示唆している。 本研究では IMF と CMF の関係を調べるため、3 次元磁気流体シミュレーションコードを用いて星形成過程でのガス雲コアの質量増大に関する数値計算を行った。分子雲コアと同程度、もしくは、より大きな質量の星ができることを再現するため、分子雲
コアの初期密度分布と重力圏の大きさをパラメータとして質量降着率の変化を調べた。 様々なパラメータで長時間の質量降着率を見積もるため中心にシンクを置き、磁場・回転を無視して計算を行った。結果、分子雲コアの周辺密度が比較的高い場合、高い質量降着率が維持されることが分かった。さらにそのような場合の質量降着率は、時間が経つにつれ Bondi 降着率に漸近することが示された。

1とが分かっている。Machida & Matsumoto 2012 でIntroduction

は、3 次元磁気流体数値計算を用いて初期の分子雲コ星はガス密度が高い分子雲コアが重力的に収縮す

アの質量をパラメータとした場合の原始星アウトフ

ることで形成され、様々な質量を持つ。また観測的ローの違いを見積もっている。その結果、原始星かにも星形成領域においてさまざまな質量を持つ星がら放出されるアウトフローは初期の分子雲コア質量

知られている。そのため質量の異なる星が形成段階

が大きくなるほど強くなり、より大きな分子雲コア

でどのように質量を獲得するのかを考えることは、
星形成過程を理解する上で非常に重要である。そしであるほど星形成効率が低くなることがわかっている。またこの計算結果を踏まえると一般的な孤立してその様々な星の質量の頻度分布を表すものとしてた分子雲コアのうち 30-70%が星に進化する。つまりIMF(Initical Mass Function) という関数が考えられている。(Salpeter 1955) この IMF の概形は星形成領域によらずほぼ同じべきを持つような形になることが知られている。また IMF に対応する関数で、星形

星形成効率が 30-70%程度となると結論付けている。
このように理論的にも観測的にも星形成過程において分子雲コアの 30-70%が星に進化する描像が一般的に考えられている。

成領域での分子雲コアの質量頻度分布を表すものとして CMF(Core Mass Function) が考えられている。

Nutter & Ward-Thompson 2007 では OrionA/B 領

一方、Takemura et al. 2021a,b では ONC 領域に
おいて IMF の最大値と CMF の最大値が一致するよ
うな結果が得られている。この領域ではフィラメン

域での prestellar core の mass function と IMF をトに沿って分子雲コアが非常に密集して存在してお比較すると、mass funtion の high mass 側のべきがり、分子雲コア間の距離は非常に近いことが知られSalpeter の IMF のべきに似ていることが分かった。

ている。そして結果からアウトフローなどのフィー

また Andre et al. 2010 で求められた CMF は IMF

ドバックを考慮すると分子雲コアが外部から質量をと概形が似ているため、CMF を 0.3 倍程度すると獲得しながら星に進化することが示唆される。

IMF とよく一致する。このように観測的に CMF とまた星形成過程を考えるうえで、分子雲コアのプ

IMF は星形成過程の上で非常に関係しており、概形

ロファイル、例えばその密度分布はシミュレーショ

が非常に似ていることが分かっている。理論的にもンの初期条件としても用いられるため非常に重要で分子雲コアと星の質量の関係密接に関わっているこ

102022 年度 第 52 回 天文・天体物理若手夏の学校ある。特に孤立した分子雲コアである Bok globule

場と回転を考慮せず計算を行い、用いた流体の基礎

は密度分布や温度分布が観測的に求められ、その密方程式は以下のように与えられる。
∂ρ ∂t + ∇ · (ρv) = 0ρ ∂v + ρ(v · ∇)v = −∇P − ρ∇ϕ度分布は critical Bonnor Ebert 球と非常に良く一致することが分かった。(Alves et al. 2001b) さらに

Nielbock et al. 20012 では分子雲コアより外側の密∂t∇2 ϕ = 4πGρP = P (ρ)

度分布は critical Bonnor-Ebert プロファイルとは一致せず、分子雲自体の密度に接続することを示して(1)いる。また分子雲コア外側の密度分布には多様性があることが確認されている。(Kandori et al. 2005)

計算領域の大きさは 4.5 × 105 au で最小セルサイズ

このような bok globule のような孤立した分子雲コは 28au に設定している。

アだけでなく、さらに分子雲コアが集まったような星形成領域でも、critical Bonnor-Ebert 球のような密2.2度分布を持つ分子雲コアが見つかっている。その分

Parameter Settings

子雲コアの周辺環境は様々であり、分子雲コア周辺の

2 つのパラメータを設定することで様々な分子雲コ密度分布も大きく異なる。Motte & Andre 2001 では

アとその周辺環境を仮定し、それぞれの初期条件の低質量原始星の星周環境に着目して、Perseus cluster

分子雲コアの進化を計算した。まずパラメータとしてや Taurus 領域での様々な分子雲コアの密度分布を

分子雲コアが収縮する範囲 (以降重力圏と呼ぶ) を設調べている。この研究では Taurus 星形成領域での原定した。その範囲の大きさは critical Bonnor-Ebert

始星のエンベロープの密度分布は、Shu の自己相似

Sphere radius(rcri ) で規格化して表し、以下式 (2) のinside-out モデル (Shu 1977) に一致することが確認 ように定義する。
されている。一方で Perseus 領域での class 0 の原始
rg = µrcri星は孤立した分子雲コアよりさらに高密度な critical
Bonnor Ebert 球に似ていて、よりコンパクトである
ことが指摘されている。このようにコアの外側の密µ = const.

(2)

今回は重力圏を rg = 2rcri , 3rcri (以降 2BE、3BE と呼ぶ) として計算を行う。

度分布は天体ごとまたは星形成領域ごとに大きく異もう一つのパラメータは初期の雲の初期密度分布なっていることが分かる。

である。今回は分子雲コア周辺の環境の違いによってこのように星の質量や分子雲コア周辺の密度分布

星形成過程がどのように変化するのかを調べるため、は多様性に富んでいる。また分子雲コア外側の密度

1BE 以内を分子雲コア内部と定義し、そこではモデ

分布も分子雲コアによって大きく異なっている一方

ル全てで同じ critical Bonnor-Ebert 球の密度分布にで、観測での制限により不定性が大きくまだよくわ

従うことを仮定する。そして 1BE より外側を分子雲

かっていないことが多い。そこで分子雲コア周辺の

コア外部と定義し、分子雲コア外部の密度をパラメー

環境が星形成にどのような役割を果たすのかを調べタとする。その初期密度分布を図 1 に示す。分子雲コるために、分子雲コアの初期の周辺密度分布をパラ

ア周辺の密度が critical Bonnor-Ebert 球の力学平衡メータとして星形成過程の数値計算を行う。

解を延長した密度より高くなる場合を比較的密度が高い場合と定義し、図 1 では破線で表している。一

2

Methods

2.1

Simulation Settings

方で分子雲コア周辺の密度が critical Bonnor-Ebert
球の密度より低くなる場合を比較的密度が低い場合と定義し、点線で表している。ここでの分子雲コア

本研究での計算には 3 次元磁気流体シミュレーショ

外部の密度分布は式 (3) を用いて一様密度に漸近するように与えている。

ンコード Nested Grid Code(Machida & Hosokawa

2013) を用いた。また長時間の計算をおこなうため磁11ρ(ξ,0) = ρ(6.45,0) [1 − Atanh{B(ξ − 6.45)}](3)

2022 年度 第 52 回 天文・天体物理若手夏の学校
次に分子雲コアから星に進化する際の質量を獲得した割合として星形成効率 (SFE) を定義し、図 3 に示す。ここでの星形成効率は Mstar /Mcore × 100 のよ
うに定義しており、ここでの Mcore は初期の分子雲
コア内部 (1BE 以内) の質量であり、全てのモデルで同じ値となっている。これより周辺密度が高い場合

図 1: 初期密度分布 (左：2BE、右：3BE)
この式のパラメータ A、B を変えることで様々な分子雲コアの周辺密度分布を作成している。また基本モデル (B1) として、外側から質量が流入することを考慮しない重力圏 1BE の孤立した分子雲コアも比較として計算した。

3

図 3: 重力圏 2BE での星形成効率

Results
ここでは中心数密度 nc が 1.5 × 106 cm−3 で重力は、SFE が 180-270%程度となることが分かる。一方圏 2BE の場合の結果を示す。計算は 0.2Myr 程度行

で周辺密度が低い場合は、SFE が 100%程度となる。

い、その場合の星への質量降着率は図 2 のようにな

これらの結果は、今回の SFE の定義や数値計算する

る。この図では、収縮が開始して 2 万年程度はどのモ際に回転や磁場が考慮されておらずアウトフローな
どのフィードバックが効いてないため、以上のような値になっている。
この結果では、分子雲コア周辺の密度分布により質量降着率や星形成効率が大きく異なることが分かった。そしてその違いは星の点源重力による質量降着率の大きさの違いが影響していると考え、星の点源
重力による質量降着を調べた。点源重力による質量降着として Bondi 降着があり、式 (4) のように与えられる。

παG2 M 2 ρext
dM
=
,
dt
Cs3

図 2: 重力圏 2BE での星への質量降着率5−3γ
 2(γ−1)(4)

しその後は分子雲コア周辺の密度の違いにより、質
量降着率は大きく異なる。周辺密度が低い場合は基本モデル (B1) と同様に質量降着率は上昇せず徐々に着率は上昇もしくは高いまま維持される。

α=2
5 − 3γ

この式での M は点源質量、ρext は周辺密度を表す。

デルでも同じ質量降着率になることが分かる。しか低下する。その一方で周辺密度が高い場合は質量降

ここでは周辺密度 ρext に着目し、時間変化する周辺密度を代入した場合の Bondi 降着を上限値、0.2Myr
程度での周辺密度を代入した場合の Bondi 降着を下限値として Bondi 降着による質量降着率がとる範囲
を図 4 に示す。この図からは分子雲コア周辺の密度が高い場合、質量降着率の高まりが Bondi 降着によ122022 年度 第 52 回 天文・天体物理若手夏の学校

5

Conclusion
分子雲コアが収縮することで星は形成され、一般的な孤立した分子雲コアの 30%程度が星に進化すると考えられている。一方で近年の観測では IMF のピー
クと CMF のピークが一致するような傾向が確認されており、分子雲コアが外部から質量を獲得する必要が示唆されている。本研究では分子雲コア周辺の密度分布をパラメータとして、分子雲コアから星へ
図 4: 重力圏 2BE での星への質量降着率と星の質量の長時間進化の数値計算を行った。その結果、特にによる Bondi 降着 (左：周辺密度が高い場合、右：周分子雲コア周辺の密度が比較的高い場合、周辺から辺密度が低い場合)

質量が流入することにより高い質量降着率が維持されることが分かった。そしてその質量降着の高まり
る質量降着率が取る範囲に含まれることが分かった。 は自己重力収縮による降着が進むと同時に Bondi 降
一方で周辺の密度が低い場合、星の質量降着率は低下

着に似た降着が起こるため発生する可能性があるこ
し、Bondi 降着の範囲も低下することが確認できた。 とを明らかにした。
今後は回転や磁場を考慮した計算を行い、アウト

4フローなどによるフィードバックが質量降着率にどDiscussionのように影響してくるかを調べる予定である。

今回の計算では分子雲コア周辺の密度を変化させることで星への質量降着率や SFE 大きく異なることが分かった。特に SFE はアウトフローなどによるReference

フィードバックが効いていないため分子雲コア周辺

Alves, J. F., Lada, C. J., & Lada, E. A. 2001b, Nature,409, 159の密度が低い場合でも 100%程度となった。そこでAndré, P., Revéret, V., Könyves, V., et al. 2016, A&A,
仮に分子雲コアの 5 割程度のガスが外側に放出され592, A54
たと考えると、分子雲コア周辺の密度が低い場合の Bondi, H. 1952, MNRAS, 112, 195
Kandori, R., Nakajima, Y., Tamura, M., et al. 2005,
SFE は 50%程度である。一方で分子雲コア周辺の密AJ, 130, 2166
度が高い場合に対しても同様な仮定をすると SFE は Kropa, P., 2001, MNRAS, 34, 14
90-135%程度となる。これより仮にアウトフローな Larson, R. B., 1969, MNRAS, 145, 271
Machida, M. N., and Matsumoto, T., 2012, MNRAS,
どによるフィードバックが効くことを考慮しても分421, 588-607
子雲コア周辺の密度が高い場合の SFE は 100%を超 Machida, Masahiro N., and Hosokawa, T., 2013, MNRAS, 431, 1719-1744
えると考えられる。これはつまり Takemura et al. Motte, F., André, P., 2001, A&A, 365, 4402021a,b で報告されたような IMF と CMF の最大値 Nielbock, M., Launhardt, R., Steinacker, J., et al. 2012,
A&A, 547, A11が同程度になる観測結果は、今回の計算のように分Nutter, D., & Ward-Thompson, D. 2007, MNRAS, 374,
子雲コア周辺からの質量降着を考慮すると十分に再4, 1413Palmeirim, P., Andre, P., Kirk, J., et al. 2013, A&A,
現可能であることが分かった。
550, A38星への質量降着率と星の質量による Bondi 降着を Penston, M. V., 1969, MNRAS, 144, 425
比較すると、周辺密度が高い場合の質量降着率の高 Salpeter E. E., 1955, ApJ, 121, 161
Shu F. H., 1977, ApJ, 214, 488
まりと Bondi 降着の範囲はほとんど一致することがTakemura, H., Nakamura, F., Ishii, S., et al. 2020,
確認できる。これにより分子雲コアから星に進化す
PASJ, 73, 3, 487-503
Takemura,
H., Nakamura, F., Kong, S., et al. 2021, ApJ,
るとき、ガスは自己重力収縮により降着すると同時
910, L6に Bondi 降着に似た降着が起こることが示唆される。 Whitworth A., Summers D., 1985, MNRAS, 214, 1

13

——–index へ戻る星惑 a03星間衝撃波により形成されるフィラメント状分子雲の化学進化
東京大学大学大学院 理学系研究科 天文学専攻
小道 雄斗

14

2022 年度 第 52 回 天文・天体物理若手夏の学校

星間衝撃波により形成されるフィラメント状分子雲の化学進化
小道 雄斗 (東京大学大学大学院 理学系研究科 天文学専攻)

Abstract
分子雲では 100 種以上の分子種が見つかっている。分子組成は分子雲の物理進化過程を反映すると考えられており、分子雲の化学組成と物理進化過程の双方の理解が重要である。本研究では、衝撃波圧縮によるフィ
ラメント状分子雲の形成・成長過程における化学進化に着目した。3 次元磁気流体計算によりフィラメント
形成を再現し、フィラメントへの降着を経験する流体素片の物理パラメータ進化を求めた。得られた物理モ
デルを用いて、化学反応ネットワーク計算を行った。その結果、比較的単純な二原子分子や氷組成は観測と整合的であったが、炭素鎖分子の存在度が観測を桁で下回った。これは従来の疑似進化モデルとは異なる傾向であり、分子雲の物理進化過程を考慮することの重要性を示唆する。

1

Introduction
分子雲と呼ばれる主に H2 から成るガス雲は、重力収縮により恒星が形成される領域である。分子雲では 100 種以上の分子種が検出されており、各分子の存在度は分子雲のもつ物理状態を反映する。実際の分子輝線観測でも、同領域の分子雲において輝線ごとに強度分布が大きく異なり、分子雲内での分子組成が不均一であることを示している (e.g., Seo et al.

2019, 図 1)。よって、分子雲の形成・成長過程をその化学組成を含めて理論的側面から理解することは重要である。
分子雲の化学組成については、分子雲の密度や温度、AV (紫外線減光率) といった物理パラメータの典図 1: 分子雲の分子輝線観測例

型値を用いて反応速度式を解く疑似進化モデルで調べられてきた (e.g. Harada et al. 2019)。このモデ(Inutsuka et al. 2015)。近年では、3 次元の磁気流体

ルは分子雲での非平衡な化学組成進化を明らかにし計算によって、星間衝撃波に伴う原子ガスの圧縮にてきた。しかし、物理パラメータは一定で、初期組よるフィラメント形成過程が盛んに研究されている

成は水素のみ分子とし、その他は全て原子と仮定し(e.g. Abe et al. 2021)。しかし、磁気流体計算では、ている。実際には物理パラメータの時間進化と共に、 C+ や CO などのガス冷却に効く原子・分子のみを扱
原子ガスから分子ガスが非平衡に生成されることかう簡易化された化学モデルが用いられており、複雑ら、疑似進化モデルは過度に単純化されている。

な分子の生成過程を考慮するには至っていない。ま一方、分子雲の物理構造に関しては、Herschel 宇宙

た、3 次元の磁気流体計算と現実的な化学反応ネッ

望遠鏡による観測からフィラメント形状の分子雲が普トワークモデルを同時に解くことは、計算コスト上遍的に存在していることが分かっている (e.g. André

困難である。そこで、本研究では、大きく 2 段階に

et al. 2010)。これを受けて、超新星などに伴う衝撃 分けた数値計算を行うことによって、フィラメント
波による星間ガス圧縮過程を経てフィラメント状分 状分子雲の形成・成長過程を経験する流体素片の化子雲が形成されるというシナリオが着目されている

学組成を調べた。まず、3 次元磁気流体計算により、

15

2022 年度 第 52 回 天文・天体物理若手夏の学校
速度でガスが定常に流入するように設定し、y,z 方向には周期境界条件を課した。
計算領域内における流体素片の物理進化をラグラ
ンジュ的に求めるために、各時刻における速度場に沿って運動するテスト粒子を計算領域内に投入する。
テスト粒子は初期に一辺あたり 16 個だけ一様に配置
し、合計 4096 個の粒子を投入した。
計算時の分解能 ∆x = 0.016 pc とした。また、計算には公開コードである Athena++(Stone et al. 2020)を用いた。

図 2: 磁気流体計算で得られたフィラメント状分子雲

2.2Results

計算の結果、x = 0 付近に衝撃波圧縮層が形成さ
星間衝撃波に伴うフィラメント形成過程を再現するれた。初期の圧縮により密度が急増した後、磁気圧(第 2 節)。衝撃波により形成される圧縮層の進化過程 や乱流によって圧縮層が膨張することで密度は一時とフィラメント周りの降着流の解析を行った。次に、 減少した。その後、時刻 2 Myr 以降に自己重力が効
解析の結果得られた各物理パラメータの時間進化をき始め、圧縮層は自己重力崩壊に伴い分裂を開始し、元に、定常流を仮定することで物理モデルを構築し、 数密度 104 cm−3 を超える領域が形成された。特に、それを用いて化学反応ネットワーク計算を行った (第y = 0 付近には十分密度の高いフィラメントが 2 つ
3 節)。以上のシミュレーションを元に、得られた化 形成され、重力不安定性による分裂過程が生じた (図学組成と観測との比較を行う。これにより、フィラ 2)。今後は 2 つのフィラメントに着目する。
メント状分子雲の物理・化学的環境を明らかにする

計算で得られたフィラメント周りの速度場を解析ことが最終目的である。

した結果、フィラメント周りには、(a) 圧縮層周囲の低密度のガスが圧縮層に突入する流れ、(b) 比較的高密度な圧縮層が分裂することに伴う降着流、という 2

2

Physical Models

2.1

Method種類の流れが存在することを確認した。(a) について
は、化学進化が進んでいないと考えられ、原子主体のガスだと考えられる。一方、(b) については、圧縮

星間衝撃波によるフィラメント状分子雲形成を再層内で化学進化が進み、分子主体のガスとなる。以現するために、自己重力と加熱冷却過程を考慮した

上より、2 種類の降着流がフィラメント周りに存在し3 次元の磁気流体計算を行った。加熱冷却過程に関ており、フィラメントの化学組成に影響すると予想

しては、Koyama, & Inutsuka 2002 で提案されていされる。以降は、複雑な分子種の形成に関わると考る加熱冷却関数を用いて評価した。計算領域は一辺えられる、(b) の降着流について調べる。

4 pc の立方体とする。初期の数密度は 50 cm−3の一前述の通り、圧縮層の進化過程は、圧縮層の形成様分布とし、x = 0 平面で衝突するよう、x 軸方向に過程と、自己重力が効き始めることによる分裂過程

相対速度 6 km/s の初期速度を与える。温度は加熱冷の 2 つに分けられる。そこで、後の化学計算で用い却関数から得られる定常状態の温度を採用する。初る物理モデルの構築に当たっては 2 段階に分けて行っ期磁場は、磁場強度 2 µG で x 軸方向から y 軸方向た。まず、圧縮層の形成過程に関しては、最終的にへ 3◦ の方向を向いた一様分布とした。また、初期に

前述の 2 つのフィラメントへ降着するテスト粒子を−3−2− 10 倍程度の亜音速乱流を与える。境界条件は、x 方向には初期条件と同じ
速度分散が音速の 10

特定し、それぞれについて、テスト粒子の記録した
物理パラメータを平均化することで、圧縮層形成時

162022 年度 第 52 回 天文・天体物理若手夏の学校の化学種と約 9000 の反応を扱う。気相、ダスト氷の表面及び内部の 3 相モデルを用いた。ダスト/ガス
比は 1/100 とし、H2, CO, C, N2 の自己遮蔽も考慮する。紫外線強度は G0 = 1 とし、宇宙線電離率は

ζ = 1 × 10−17 s−1 とする。
初期組成は太陽金属量での HI ガスとする。磁気流体計算結果の解析で得られた、各フィラメントへ降着する流体要素のガス密度、ガス温度、ダスト温度、

AV の時間発展を入力することで、各時刻での化学組成を出力する。計算は 2 段階に分けて行った。まず、
上記の初期組成と、圧縮層形成時の物理進化モデル(図 3 上) を用いて化学進化計算を行う。次に、その終状態での化学組成を新たに初期条件として、フィラメントへの定常降着流の物理進化モデル (図 3 下)を用いて計算を行った。

3.2

Results and Discussion

計算から得られた、AV = 5 mag での化学組成について、分子雲 TMC-1、L134N の平均組成 (Terzieva,& Herbst 1998) と比較した結果を図 4 上に示す。結
図 3: フィラメントへ降着する流体素片の経験する

果、両モデルとも CO や CH3 OH などの基本的な分圧縮層形成時の平均的な物理進化 (上) と、フィラメ

子については観測と整合的な結果が得られた。一方、ントへの定常降着流の物理進化 (下)。時間の原点は、 C3 H2 や HC3 N 等の炭素鎖分子については観測を大それぞれ圧縮層への突入時及びフィラメントへの降 幅に下回った。この結果は従来の疑似進化モデルで着開始時に対応する。ここでは片方のフィラメント

得られる化学組成とは異なる傾向であり、分子雲のに関するモデルのみを示した。

物理進化を考慮することの重要性を示唆する。
炭素鎖分子の形成が進まなかった原因を調べるたの平均的な物理パラメータを求めた。AV は圧縮層

めに、宇宙線電離率を ζ = 1 × 10−16 s−1 に上げた場の x 軸方向の柱密度から求めた。次に、フィラメン

合と、モデルの AV を 3 mag だけ加算した場合で計ト形成・分裂過程については、先ほどの (b) の降着算を行った。宇宙線強度が高い場合は CO が作られ流を定常流と仮定することで物理パラメータ進化をやすく、AV が高い場合は紫外線による炭素鎖分子の求めた。ここで、ダスト温度については、(Hocuk et

破壊が抑えられると予想される。先ほどと同じ段階

al. 2017) の式 (8) で評価した。その結果、図 3 に示での組成を図 4 中に示す。その結果、両者とも元のすような物理モデルを作成した。

物理モデルの場合よりも炭素鎖分子の存在度が増加したが、観測の再現には至らなかった。
今回は、1 つの磁気流体計算による結果を用いて

3

Chemical Calculations

3.1

Method計算を行った。しかし、磁場の角度や衝突速度などの初期条件に応じて、圧縮層の物理進化は変化する。
よって、より適切な初期条件を元にした物理モデル

化学計算では、Furuya et al. 2017 で用いられた化学反応ネットワーク計算コードを使用した。945 種

により、炭素鎖分子の組成をより改善できる可能性
がある。

17

2022 年度 第 52 回 天文・天体物理若手夏の学校段階まで引き継がれる可能性がある。

4Conclusion
磁気流体計算と、そのポストプロセスで化学組成計算を行うことにより、これまでの疑似進化モデル
よりも現実的な物理パラメータにより化学組成計算を行った。その結果、単純な二原子分子や氷組成については観測を再現できたものの、炭素鎖分子については観測を大きく下回った。
今後は、圧縮層の形成段階にまず着目し、磁気流体計算における広い初期条件パラメータ範囲を考慮することで、物理進化過程の違いに伴う化学組成の差異について研究を進める。

Acknowledgement
本研究の磁気流体計算については、国立天文台

CfCA の Cray XC50 を使用しました。

Reference
Abe, D., Inoue, T., Inutsuka, S-i., & Matsumoto, T.
2021, ApJ, 916, 83
André, P., Men’shchikov, A., Bontemps, S., et al. 2010,
A& A 518, L102
Furuya, K., Drozdovskaya, M. N., Visser, R., et al 2017,
A& A 599, A40

図 4: 上:AV = 5 mag における、各フィラメントの気相中の分子組成と観測。中:片方のモデルについて、 Harada, N., Nishimura, Y., Watanabe, Y., et al. 2019,
ApJ, 871, 238

各宇宙線電離率と AV 加算時の存在度。下:各モデルでの氷相中の氷組成。灰色線は 1 桁分の誤差を示す。 Hocuk, S., Szuc̋s, L., Caselli, P., et al. 2017, A& A, 604,A58

原始星ではそれまでに形成されたダスト表面上の氷が昇華することから、分子雲における氷組成と原Inutsuka, S-i., Inoue, T., Iwasaki, K., & Hosokawa, T.
2015, A& A, 580, A49

始星の観測で得られる気相中の化学組成を比較する

Koyama, H., & Inutsuka, S-i. 2002, ApJ, 564, L97

ことができると考えられる。そこで、今回の計算でÖberg, K., Boogert, A. C. A., Pontoppidan, K. M., etal. 2011, ApJ 740, 109

得られた氷の組成を、大質量および小質量原始星での典型的な組成 (Öberg et al. 2011) と比較した結果を図 4 下に示す。この結果、両モデルともフィラメ
ントでは氷が十分に生成されており、相対量も観測

Seo, Y. M., Majumdar, L., Goldsmith, P. F., et al. 2019,
ApJ, 871, 134と整合的となった。よって、フィラメントへの降着

Stone, J. M., Tomida, K., White, C. J., & Felker, K. G.
ApJS, 249, 4過程を経験した流体素片での氷組成進化が、原始星

Terzieva, R., & Herbst, E. 1998, ApJ, 501

18

——–index へ戻る星惑 a04

Corona Australis 分子雲 IRS7 領域における分子雲コアの性質と星団形成
名古屋大学大学院 理学研究科
深谷 直史

19

2022 年度 第 52 回 天文・天体物理若手夏の学校

Corona Australis 分子雲 IRS7 領域における分子雲コアの性質と星団形成
深谷 直史 (名古屋大学大学院 理学研究科)

Abstract
小質量星形成のモードの一つである集団的星形成においてどのように高密度コアが形成されるかは十分にはわ
かっていない。我々は星団形成領域における小質量星形成の初期条件の理解を目的として、最近傍 (d = 149 pc)の活発な小質量星団形成領域である Corona Australis 分子雲を対象とした ALMA データの解析を進めている。今回 IRS7 領域の 1.3 mm 連続波と C18 O (J=2–1) 輝線を解析した。その結果、連続波点源を 5 つ同定した。これらのうち 4 つは対応天体が存在する class 0 原始星及び class 0/I 原始星と考えられ、残り 1 つは対応天体がなく新検出である。また、点源を繋ぐ弱く広がった連続波放射も見られた。一方 C18 O の強度分布から、連続波放射に沿ってフィラメント状に広がった構造と、その中に埋もれた数 100 AU 程度の C18 O クラン
プが 10 個程度検出された。これらのクランプおよび連続波点源は 1000 AU ほどの間隔でフィラメント内部に分布していた。C18 O イメージに対してコアを同定し、導出したコアの典型的物理量は柱密度 ∼1023 cm−2 、
密度 ∼107 cm−3 、質量 ∼0.01 M⊙ であり、これらは広がった連続波放射から求めたものとほぼ同程度であった。これらの結果は分子雲中で多数形成された高密度フィラメントが短い波長 (1000 AU 程度) で星形成に直結する微小コアに分裂し星団形成に至ることを示唆する。

1

ら、フィラメントの幅は普遍的に約 0.1 pc であるこ

Introductionと (Arzoumanian et al. 2019)、フィラメント中の分小質量星の星形成は孤立的な星形成と星団を形成

子雲コアの間隔はフィラメントのジーンズ長程度でする集団的な星形成の 2 つのモードがある。孤立した小質量星の星形成については数十年前から多くの観しかし、集団的星形成に関しては孤立的星形成と測的、理論的研究が行われており、その理解が進んできている。一般的に考えられている小質量星の形成れが重力収縮して星形成が進行するというものである。近年の研究から分子雲コアの進化はフィラメント
状分子雲の中で起こることがわかってきている。フィ
ラメントの単位長さあたりの質量（line mass, Mline ）が、臨界値（Mline, crit = 2c2s /G）を超えると重力的に不安定になり、フィラメントが長さ方向に分裂し、星形成が起こる分子雲コアが並んで誕生するという理論的研究が Inutsuka & Miyama (1997) などによっ
し座領域（Aquila Rift）の観測から、星形成が起きの初期条件の違いや、分子雲が分裂し、分子雲コアが密集して形成される仕組みなど未だ明らかになって過程は、分子雲が密度の高い分子雲コアに分裂し、そて示された。観測的研究においても、Herschel のわあること (Zhang et al. 2020) が示唆されている。

いない問題が多い。そこで我々は最近傍 (d = 149 pc)の活発な小質量星団形成領域である Corona Aus-

tralis 分子雲を対象とした ALMA データの解析を進めている。Corona Australis 分子雲の高密度領域のALMA Atacama Compact Array（ACA）による観測 (2018.A.00056.S, 分解能 ∼ 5′′ .8) では、TP array
データを組み合わせることで幅 0.1 pc のフィラメン
ト内部により幅が小さく（≲ 0.02 pc）、高密度なサ
ブフィラメントが埋め込まれていることを確認した（図 1、Nishioka et al. in prep）。この観測では 81 本のサブフィラメントが同定され狭い範囲に密集してているコアの大部分は Mline, crit を超える不安定な存在することがわかった。このように星団形成領域

フィラメントの中にあることが示された (André et

ではこれまでの観測とは違った結果が得られており、al. 2010)。これまでのフィラメントの観測的研究か

集団的星形成の理解にはさらなる星団形成領域の観20

2022 年度 第 52 回 天文・天体物理若手夏の学校
測的研究が必要であると言える。したがって、我々 その中の Corona Australis 分子雲 IRS7 領域に含まは ALMA ACA の観測によって検出された Corona

れる 2 視野のデータを合体させて使用した。解析Australis 分子雲高密度領域のサブフィラメント内部 範囲は 51′′ × 35′′ （7600 AU×5300 AU）である。こをより高い分解能で調べることで、集団的星形成ののデータは ALMA 12 m array で得られたものであ理解を深めることを目的に、ALMA 12 m array にり、観測スペクトルは 1.3 mm 連続波、12 CO(J=2–

よって得られた観測データを解析した。本研究の解1)、13 CO(J=2–1)、 C18 O(J=2–1)、SiO(J=5–4)、

析領域では 2 本のサブフィラメントが含まれ、このH2 CO(J=3(0,3)–2(0,2))、N2 D+ (J=3–2) である。今

解析から集団的星形成の機構を探求する。

回は 1.3 mm 連続波と C18 O についてイメージングと解析を CASA (Common Astronomy Software Ap-

plications) を用いて行った。C18 O については IRS7領域を含むより広域な ALMA ACA 及び TP array
データ (図 1) を組み合わせた。最終的に得られたデータの合成ビームサイズは 1′′ .0 × 0′′ .8（∼134 AU）、ノ
イズレベルは 1.3 mm 連続波で 1.1 mJy/beam、C18 O
で 0.31 K、速度分解能は C18 O で 0.2 km/s となった。

3

Results
解析で得られた 1.3 mm 連続波マップを図 2 に示す。連続波マップから連続波点源を 5 つ同定し、これらについて 2 次元ガウシアンフィットを行い、測定したピーク強度とフラックス密度を表 1 に示す。これらの連続波点源のうち、B9 は Class 0 原始星 (Groppi etal. 2007)、IRS7B、IRS7A、CXO34 は Class 0/I 原
始星 (Peterson et al. 2011) として先行研究で同定されており、Source 1 については今回の解析によって新しく同定されたものであった。また、同領域の Spitzerのデータから B9 と Source 1 については対応する赤外線源が存在しないこともわかった。これらの連続
図 1: 上図：Corona Australis 分子雲高密度領域の波点源は 800–1100 AU の間隔で分布していた。一方

C O ピーク強度図。下図：上図中央部の強度が大き で、連続波点源を繋ぐように広がった ≳2 mJy/beam
い領域周辺を拡大した図。黒コントアは本研究の解 の弱い連続波放射が分布していることも確認した。
18

析領域を示す。

次に、解析で得られた C18 O の積分強度図を図 3に示す。積分強度図の強度分布から 2 本のフィラメ
ント状の分子雲が中央部の南東–北西方向と上部の南

2

北方向に伸びて存在していることを確認した。また、

Observations

強度が大きいクランプ状の構造も 10 個程度存在し、ALMA Science Archive から Cycle 7、Band 6で 観 測 さ れ た デ ー タ (P.I.

Mardones,

project

code:2019.1.01792.S) を使用した。このデータは近これらのクランプや連続波点源、弱く広がった連続波放射は 2 本のフィラメント状分子雲に沿って分布していることが確認できる。

C18 O 輝線強度分布からコアの同定と、局所熱力

傍の様々な星形成領域に存在する原始星からのアウ

トフローの観測的研究のために得られたものであり、 学平衡を仮定し、水素分子の柱密度と密度、質量を21

2022 年度 第 52 回 天文・天体物理若手夏の学校

0.002

Jy/beam
0.004

0.006

C18O Integrated Intensity (K km/s)

5

5

10

15

20

25

-36°57'10"

120"

Declination (J2000)

Declination (J2000)-36°57'10"

0.008 0.010 0.0120.014

3

2

4

30"

40"

20"

30"

40"

1500 AU19h 01m 57s

56s

55s

1500 AU54shm

19 01 57

Right Ascension (J2000)
図 2: IRS7 領域の 1.3 mm 連続波マップ。図中のラ

1

B9

216

330

2

IRS7B

280

338

3

IRS7A

21.2

26.7

4

CXO34

15.2

15.8

5

Source 1

11.3

14.4

56
55
Right Ascension (J2000)s

54s

トアは 1.3 mm 連続波の 3σ 、6σ 、10σ 、20σ 、40σ を示す（σ=1.1 mJy/beam）。

表 1: 同定した連続波点源のデータ
ID Name
Fpeak
Fν

(mJy)

s

図 3: IRS7 領域における C18 O 積分強度図。白コン

ベルは表 1 の同定した連続波点源の ID を示す。

(mJy beam−1 )s

ED

求めた。水素分子の柱密度は励起温度を 25 K とし、A

C18 O と水素分子の存在比は先行研究 (Yonekura etBC

al. 1999) のものを使用して求めた。コアの同定は得られた水素分子の柱密度のマップからピークの半値を境界として、それと連続波点源または弱く広がった
連続波放射の分布が一致するものをコアとして定義した。コアの半径はコアの広がりと同等の円の半径として定義し、コアの密度は得られたコアの質量と

図 4: IRS7 領域の水素分子の柱密度マップ。白コン

半径から求めた。図 4 は C18 O から得られた水素分子

トアは 1.3 mm 連続波の 2σ 、3σ 、10σ 、20σ 、40σ をの柱密度マップと同定したコアをコントアで示したものである。表 2 に求めた各コアの物理量を示した。
今回の方法では連続波点源が付随する星ありコアを示す。黒コントア、灰色コントアはそれぞれ同定した星ありコアと星なしコアを示す。

3 個、連続波点源が付随せず弱く広がった連続波放の方法で広がった連続波放射の強度とフラックス密射のみが付随する星なしコアを 2 個同定した。さら度から水素分子の柱密度と密度、質量を求めた。密に、星なしコアについては Kauﬀmann et al. (2008)22

2022 年度 第 52 回 天文・天体物理若手夏の学校

コア

表 2: 同定したコアの物理量
平均柱密度
密度
質量

半径

(AU)

(cm

599

コア B
コア C

コア A

コア D
コア E

−2

)

(cm

−3

)

連続波点源の付随(M⊙ )∼ 8 × 1022

∼ 8 × 10

6

∼ 0.05あり

355∼ 8 × 1022

∼ 1 × 107

∼ 0.02なし

307∼ 7 × 1022

∼ 1 × 107

∼ 0.01あり

353∼ 8 × 1022

∼ 1 × 107

∼ 0.02あり

267∼ 8 × 1022

∼ 2 × 107

∼ 0.009

なし

度はコアの奥行きをコアの直径と仮定して柱密度か~ 105-6 cm-3ら求めた。広がった連続波放射から求めたコアの物

サブフィラメント（幅 0.02 pc）

理量はコア B で平均柱密度 ∼ 8 × 1022 cm−2 、密度∼ 2 × 107 cm−3 、質量 ∼0.002 M⊙ であり、コア E で
平均柱密度 ∼ 1 × 10

23

cm

−2

、密度 ∼ 1 × 10 cm
7−3~ 0.01 M⦿

、~ 104 cm-3

質量 ∼0.005 M⊙ 、となった。C O から求めた物理18~1000 AU

Herschelのフィラメント幅~0.1 pc

量と比較すると約 2 倍以内の差になっており両者で~ 107 cm-3

おおよそ同等の結果が得られた。
図 5: 集団的星形成の概念図

4

星形成におけるフィラメント中の星形成のスケール

Discussion

が小さくなり、より幅が小さく密集したサブフィラ

サブフィラメント中に連続波点源や C18 O クラン

メント中で星形成が起きる可能性がある。これを検プが平均して約 1000 AU の間隔で存在しており、そ証するためには同定されたサブフィラメントの幅やれらが様々な進化段階を持つ点は、孤立的星形成に

line mass を求め、集団的星形成におけるフィラメン

おけるフィラメントのジーンズ長の間隔で星形成がトの性質を調べるほか、今回同定した星なしコアが進む点と似ている。したがって、集団的星形成領域

星形成し得る分子雲コアであるのか検証する必要がでも孤立的星形成領域と同様にフィラメントの分裂ある。

により分子雲コアの形成が進むが、集団的星形成領域におけるそのスケールは孤立的星形成領域におけるスケールを一桁以上小さくしたものである可能性

Referenceがある（図 5）。一方で、このような状況が重力不安定性によって起こるとは考えにくく、サブフィラメ
ントの増加や分裂は外部からの衝撃波による圧縮や乱流の増幅などによって起こる可能性がある。

5Conclusion

Inutsuka S.-. ichiro., & Miyama S. M., 1997, ApJ, 480,
681.
André P., et al. 2010, A&A, 518, L102.
Arzoumanian D., et al. 2019, A&A, 621, A42.
Zhang G.-Y., et al. 2020, A&A, 642, A76.
Groppi C. E., et al. 2007, ApJ, 670, 489.
Peterson D. E., et al. 2011, ApJS, 194, 43.
Yonekura Y., et al. 1999, PASJ, 51, 911.
Kauﬀmann J., et al. 2008, A&A, 487, 993.

星団形成領域である Corona Australis 分子雲 IRS7領域を解析し、サブフィラメント中に原始星や星なしコアが約 1000 AU の間隔で分布していることを確認した。この結果から集団的星形成の描像は孤立的

23

——–index へ戻る星惑 a05分子ガス同位体比の新しい測定手法を用いた TW Hya
円盤における不均一な 12CO/13CO 比の発見
総合研究大学院大学 物理科学研究科 天文科学専攻
吉田 有宏

24

2022 年度 第 52 回 天文・天体物理若手夏の学校

分子ガス同位体比の新しい測定手法を用いた TW Hya 円盤における不均
一な 12 CO/13 CO 比の発見
吉田 有宏 (総合研究大学院大学 物理科学研究科 天文科学専攻)

Abstract
分子雲から惑星系へと至る物質の起源と進化をたどる上で、同位体組成は有力なトレーサーである。特に、原始惑星系円盤の一酸化炭素ガス同位体組成を測定することには意義がある。しかし、CO 回転輝線の線中心は一般に光学的に厚く、柱密度の導出が困難である。一方で、輝線は熱運動により広がっており、「すそ」の部分は光学的に薄くなる。従って、複数の同位体輝線のすそを観測することで同位体比を求めることができる。我々はこのアイデアを定式化し、原始惑星系円盤における分子ガス同位体比をモデル非依存に測定できる手法を開発した。さらに、本手法を ALMA 望遠鏡による TW Hya 円盤の 12 CO J = 3 − 2, 13 CO J = 3 − 2
輝線の観測アーカイブデータに適用し、12 CO/13 CO 比を測定した。その結果、12 CO/13 CO は半径 100au以内の領域では星間空間の平均値の約 0.3 倍に減少し、半径 100au よりも外側の領域では約 1.4 倍かそれ以上に増大していることがわかった。単一の円盤において 12 CO/13 CO の大きな変化が観測されたのはこれが初めてであり、炭素同位体比が惑星系物質のトレーサーとして機能しうることが示唆された。

1

背景

時間的・空間的に変化しうる。そこで形成される巨大ガス惑星の大気組成は、円盤ガスの組成を反映す

惑星系は原始惑星系円盤で形成される。近年、観ると考えられるので、巨大ガス惑星大気と円盤ガス測技術の発達により、いくつかの円盤では惑星が今の C/O 比を比較することで、巨大ガス惑星の形成過まさに形成されていることの証拠が捉えられている

程を制限できる可能性がある。実際、C/O は太陽組
(e.g., Andrews, 2018, ApJL, 869, L41)。ここで興味 成では 0.6 程度であるが、いくつかの原始惑星系円
深いのは、惑星系と原始惑星系円盤の物質的関係で 盤では 1 よりも大きいということが観測的に示されある。惑星系を構成する物質はどこで形成され、どている (Bergin et al., 2016) ものの、その具体的なのように進化してきたのだろうか。この謎を解くた

値を測定することは容易ではない。同位体比に着目めに強力な道具となるのが同位体比である。例えば、 すると、C/O>1 の場合には C/O の増大に従って、水分子中の重水素比 (D/H 比) は、その水分子が低12C16 O/13 C16 O が小さくなる (Lee et al., in prep)。

温環境で生成された場合に高い値を示す。地球の水

これは一酸化炭素同位体の零点エネルギーの違いにと分子雲中の水分子はいずれも、宇宙での平均組成より起こる熱交換反応

に対して高い D/H 比をもつことが知られており、地13

球の水の少なくとも一部が分子雲に起源を持つことの証拠だと考えられている (Öberg & Bergin, 2021,

C+ +12 CO ⇌12 C+ +13 CO + 35K(1)に 起 因 す る も の で あ る 。こ れ を 用 い る こ と で 、

Phys. Rep., 193, 1)。他にも、酸素、窒素、炭素な

12

どの同位体比が太陽系物質やさまざまな天体で測定

C16 O/13 C16 O から C/O の具体的な値を制限できる可能性がある。

されている。

円盤での CO 同位体組成の測定に関する研究は多同位体比に加えて、惑星系物質のトレーサーとしくない。例えば、Smith et al. (2009, ApJ, 701, 163)て近年注目されているのが、炭素-酸素比 (C/O 比)は一方の円盤がもう一方の中心星と同じ視線方向にである。円盤の化学モデル計算によれば、円盤ガス中の C/O 比はダストへの炭素と酸素の凍結により、あるような特殊な連星系において、CO の吸収線を用いて同位体組成の測定を行っているが、この手法

25

2022 年度 第 52 回 天文・天体物理若手夏の学校
アーによるスペクトルのコンタミネーションを避け

は一般の円盤に関しては適用できない。また、Zhang

et al. (2017, Nat. Astron., 1, 0130) は TW Hya 円 る必要がある本手法の適用に適している。
盤において、12 C18 O, 13 C18 O 輝線を観測し、モデル

観測データは ALMA Science Archive に公開されフィッテングによって円盤内側領域での 12 C/13 C をているものを用いた。角分解能と速度分解能はそれ求めている。しかしながら、この結果は不定性が大

ぞれ 0.5 秒角、0.06 km/s 程度であり、積分時間は 3

きく、またモデル依存性がある。そこで、本研究で時間程度である。
我々はまず、観測ビジビリティから 12 CO J = 3−2,

は、新たに円盤の CO 同位体組成をモデル非依存に13測定する手法を提案する。円盤において CO 同位体
1216

CO J = 3 − 2 のイメージキューブを作成した。さ

組成の測定が難しいのは、 C O 輝線の線中心が一らに、円盤の運動がケプラー回転していると仮定し般に光学的に厚いためである。しかし、輝線はガスて、各半径で、スペクトルを線中心を揃えつつスタッの熱運動による拡がりをもつため、「すそ」の部分にクした。これは S/N を向上するためである。さらに、は光学的に薄い部分が存在する。したがって、「すそ」 得られた各半径のスペクトルから式 (2) を計算し、229/28を二つの同位体種について観測すれば、個々の同位節で述べたように、t12体の柱密度自体は測定できなくとも、それらの比でトした (図 2)。もし、円盤全体が同じ同位体比であれあれば求めることができる。

ば、t12

29/28

/t13 を t12 に対してプロッ/t13 は t12 ≪ 1 である値に漸近するはずである。我々は、12 CO/13 CO 比を星間物質の値に設定2したモデルを観測と同様の条件で解析した。その結手法

果は図 2 の左図に点線で示したが、観測点とは整合

円盤の典型的な条件下では、12 CO J = 3−2 輝線のすそには、光学的厚みが 1 よりも小さくなる速度域が∼ 0.1 km/s 程度存在する。同じ速度域で 13 CO J =3 − 2 も観測し、12 CO J = 3 − 2 との比をとること的でないことがわかる。これは、12 CO/13 CO が円盤の場所によって異なっていることを示唆する。次に、
12

CO/13 CO 比が半径 100 au 以内で 20、100 au 以遠

で 100 に階段関数状にしたモデルを作成し、同様にによって同位体比を求めることができる。より具体解析した (図 2 右図・実線)。この結果は観測とよく的には、まず、一致している。また、比較的 S/N の良い半径 70-110{

Ik (v)
tk ≡ − ln 1 −
B(T )

29/28

au においては、t12}/t13 の t12 ≪ 1 の領域に定数関数をフィットすることで具体的な 12 CO/13 CO の値(2)

を 21 ± 5 と決めることができた。

を計算する。ここで、添字 k は 12 CO と 13 CO のいずれかを表し、I, B(T ) はインテンシティとそれぞれの周波数における温度 T でのプランク関数を示す。 429/28
さらに、t12 /t13 を t12 に対してプロットする。こ29/28の時、t12 ≪ 1 では、t12 /t13 が同位体比でほとん

ど決まる値に収束する。詳しい定式化と手法の模擬観測による検証は Yoshida et al.(2022) を参照されたい。

議論
TW Hya 円盤の 12 CO/13 CO 比は、円盤の半径 100

au 以内で 21、以遠で 100 程度と一つの円盤内で大きく変化していることがわかった。これはいずれも星間物質の平均値 (∼ 69) からは大きく外れた値である。
まず、比較的内側の星間物質よりも非常に小さい値については、C/O > 1 の条件下での同位体交換反

3

観測と結果

応 (式 1) の寄与が考えられる。実際、TW Hya 円盤
では C/O > 1 であることが先行研究で示唆されてい本手法を最初に適用する天体としては、TW Hya

る (Bergin et al., 2016, ApJ, 831, 101) が、今回観周りの円盤を選んだ。この円盤は地球に最も近い円

測された値を説明するには C/O ∼ 2.5 程度が必要で盤であり、最もよく研究されている円盤の一つであある。ダストの進化を考慮した化学反応計算によるる。さらに、天球面からの傾斜角が小さく、速度シ26

2022 年度 第 52 回 天文・天体物理若手夏の学校

29/28

図 1: t12/t13 v.s. t12 。左図は 12 CO/13 CO= 69 のモデル (点線)、右図は 12 CO/13 CO が階段関数状に変化するモデル (実線) の結果もプロットしている。
と、気層の C/O は 10 程度まで上昇しうる (Krijt etの分子ガス同位体比を測定する新手法を開発し、TW

al., 2020, ApJ, 899, 134)。
次に、円盤外側での高い値については、12 CO と

Hya 円盤の 12 CO J = 3 − 2, 13 CO J = 3 − 2 輝線の

13

CO の氷への束縛エネルギーの微小な違いと、COHya 円盤の 12 CO/13 CO は、半径 100 au より内側でガスの大域的な枯渇で説明できることがわかった。実は 20 程度、外側では 100 程度以上であることがわ

際、TW Hya 円盤では CO の水素分子に対する存在

かった。これらの値は、内側では C/O > 1 の環境下

観測アーカイブデータに適用した。その結果、TW度が星間空間での値に対して 1-2 桁程度小さくなっ での同位体交換反応、外側では 1 2CO と 1 3CO の氷ていることが示唆されている。また、この枯渇の程度と観測された12への束縛エネルギーの違いと、CO ガスの枯渇で説13

CO/ CO の値は整合的であった。 明できる。一つの円盤内で大きな CO 同位体比の変最後に、先行研究によれば、半径 100 au 付近には化が観測されたのは初めてであり、炭素同位体比もガスのギャップがあることが示唆されている。具体

惑星系物質のトレーサーとして有用である可能性が的なメカニズムは不明であるが、このギャップによっ 示唆された。
て円盤が内外に分断されていれば、12 CO/13 CO にも何らかの影響を与える可能性がある。

詳細については、Yoshida et al., 2022, ApJ, 932,

126, “A New Method for Direct Measurement of
Isotopologue Ratios in Protoplanetary Disks: A

5

結論

Case Study of the 12 CO/13 CO Ratio in the TW
Hya Disk” をご覧いただきたい。

同位体比は分子雲から惑星系へと至る物質進化の有力なトレーサーである。我々は、原始惑星系円盤

27

——–index へ戻る星惑 a06

木星の氷微惑星散乱による雪線内側への水供給についての理論的研究名古屋大学大学院 理学研究科
安田 郁斗

28

2022 年度 第 52 回 天文・天体物理若手夏の学校

木星の氷微惑星散乱による雪線内側への水供給についての理論的研究
安田 郁斗 (名古屋大学大学院 理学研究科)

Abstract
地球は生命を育む水惑星であるが、海洋の質量は全地球質量の 0.02%程度と、地球の水含有率は非常に低い。つまり地球は氷を集積しない雪線の内側で形成され、形成後に少量の水が供給されたと考えるのが自然である。リュウグウをはじめとする C 型小惑星は水を多く含んでいたと考えられており、はやぶさ 2 ミッションの成果からも、木星や土星よりも外側で形成し、小惑星帯そして、地球近傍へと移動したことが示唆されている。氷微惑星の移動は木星形成後に期待される。木星形成直後には、固体核材料である氷微惑星とガス集積を可能にするガス円盤が残っていると考えられる。従ってガス集積した木星の強い重力により微惑星は
太陽系の内側にも外側にも散乱され、原始惑星系円盤によるガス抵抗で軌道が安定化する。先行研究ではこうしたプロセスで、小惑星帯への微惑星輸送と地球への水供給の可能性が示唆された (e.g., Raymond andIzidoro 2017) が、初期条件や仮定への依存性が強いとも考えられるため、本研究では輸送過程を解析することで、輸送のパラメータ依存性などを知り、輸送メカニズムを明らかにすることを目指す。N 体シミュレー
ションにより、ガス抵抗を考慮した制限三体問題として木星と微惑星の軌道進化を計算し、様々な解析を行う。今回は軌道進化を遠日点距離の変化で見ることで、内側への重力散乱では遠日点距離がほとんど変化しないが、ガス抵抗により遠日点が木星軌道から離れて安定軌道を取ることを明らかにした。こうした軌道の定性的な進化や今後の展望についても議論する。

1にはケンタウルス族など少数の小惑星しかない。一Introduction

方、木星形成には多くの氷微惑星が必要であり、木地球は生命を育む水惑星であるが、海洋の質量は星形成直後には多くの氷微惑星が木星軌道の周りに全地球質量の 0.02%程度と、地球の水含有率は非常残っているだろう。そのため、氷微惑星の移動は木星に低く、地球は主に岩石で形成されている。つまり形成後に期待される。また、木星形成直後にはガス

地球は氷を集積しない雪線の内側で形成され、形成集積を可能にするガス円盤も残存すると考えられる。

後に少量の水が供給されたと考えるのが自然である

従ってガス集積した木星の強い重力により微惑星は(e.g., Morbidelli et al. 2012)。太陽系の母体である

太陽系の内側にも外側にも散乱され、原始惑星系円

原始太陽系星雲では太陽に近い場所では温度が高く、遠ざかるにつれて温度が低くなる。雪線は水が気相盤によるガス抵抗で軌道が安定化することで、氷微惑星が雪線の内側にも安定な軌道をとり得る。こので存在する内側領域と固相で存在する領域の境界でような水の輸送過程について様々な先行研究があるある。
一方、リュウグウをはじめとする C 型小惑星は水を多く含んでいたと考えられている。はやぶさ 2 ミッションの成果から、リュウグウは木星や土星よりもが (e.g., Raymond and Izidoro 2017)、初期条件や仮定への依存性が強いとも考えられる。よって本研究
では水の輸送過程を解析し、輸送メカニズムの解明を目指す。

外側で形成し、小惑星帯に移動し、その後地球近傍へとさらに移動したことが示唆されている。つまり、雪線の外側でできた氷微惑星が雪線内側まで移動が起こり、その一部が C 型小惑星になったことが示唆されている。

木星の質量はおよそ 1030 g であり、100 km サイズの微惑星であれば質量は 1022 g 程度なので、その差は 8 桁ほどもあり、木星形成直後における木星まわりの微惑星の軌道進化について、微惑星同士の重力は木星重力に比べて十分小さく、制限三体問題と現在の木星から海王星までの外惑星の軌道の周り

292022 年度 第 52 回 天文・天体物理若手夏の学校
動方程式はして扱うことができる。また、木星の軌道は離心率が小さいので、このような制限三体問題においては、r̈ = −Gm1

ヤコビ積分と呼ばれる保存量が存在する。そのため、第２節で保存量ヤコビ積分を紹介する。そして、第３、４節で N 体シミュレーションにより、ガス抵抗を考慮した制限三体問題として木星と微惑星の軌道

r1
r2
− Gm2
|r 1 |3|r 2 |3(1)となる。次に回転座標系（図 1 の橙の座標系）に移り、その座標系での運動方程式に速度をかけ、エネ
ルギー積分し、再び座標系を（図 1 の青の座標系に）進化を計算し、ヤコビ積分のような物理量をはじめ戻すことで以下の Ej が保存することが導かれる。
(
)
1 2
GM⊙
GMJ
Ej ≡ v −
+
−h·n
(2)
2
|r|
|r − r J |とする軌道に関する各種パラメータを足掛かりとした解析を行う。第 5 節で、太陽系内側に輸送される
微惑星の軌道進化について、定性的な理解が進んだ点について議論する。

ここで v はデカルト座標系（図 1 の青の座標系）における小天体の速さであり、h は単位質量当たりの2

角運動量ベクトルである。r2 が十分大きいとすると、Jacobi Integral

ヤコビ積分は小天体の軌道長半径 a、軌道離心率 e、
軌道傾斜角 i を用いて次のように書ける。

この節では本研究で重要となるヤコビ積分について簡単に説明する。本研究で扱う太陽・木星・微惑Ej ≈星のように、中心星や惑星の二体と質量が桁で小さ√
1+ a (1 − e2 ) cos i
2a(3)い小天体の三体の運動を考える時（制限三体問題）、 式 (3) はティスランパラメータとも呼ばれる。
小天体についての保存量が存在し、それがヤコビ積分である。図 1 は制限三体問題を考える際の座標と3

Methods
本研究では微惑星や惑星の軌道進化を計算するため、4 次エルミート法 (Makino and Aarseth 1992)
を用いて N 体シミュレーションコードを開発した。
微惑星をテスト粒子として扱うことで、木星質量の惑星と微惑星の軌道を、ガス抵抗を考慮した制限三体問題として計算した。惑星の軌道面は原始惑星系
円盤の赤道面と一致させ、軌道長半径が aJ = 6 AUの円軌道とした。
第２節で示したように（式 (3)）、ヤコビ積分は微惑星の軌道長半径 a、離心率 e、傾斜角 i で決まる

図 1: 制限三体問題における座標系と文字の定義。

ので、これらの初期値を固定して計算を行う。a =文字の定義についてを表す図である。図 1 にあるよaJ − ∆a0 , e = 0.05, i = 0.025 とした。ここで i のうに、中心星の質量を m1 、惑星の質量を m2 、小天

単位はラジアンである。また、∆a0 を 2rH として rH体の質量を m とする。ただし、m1 , m2 ≫ m であは惑星のヒル半径を用いて与える。木星の軌道長半る。中心星と惑星が円軌道で運動しているのを考え、 径 aJ 、木星質量 MJ 、太陽質量 M⊙ を用いて(
)1/3
その軌道面を xy 平面、公転角速度ベクトルを n と
MJ
aJ
(4)
r
=
H
する。小天体の位置ベクトルを r 、中心星から小天
3M⊙
体に向かうベクトルを r 1 、惑星から小天体に向かうベクトルを r とする。この時、小天体についての運 と表せる。今回は rH ∼ 0.4 AU 程度となっている。
2

初期軌道の位相だけを様々に変えて計算を行った。

30

2022 年度 第 52 回 天文・天体物理若手夏の学校
た微惑星を３体だけプロットした結果を示す。

微惑星についての運動方程式は以下のように重力の項 F grav とガス抵抗の項 F drag を考えている。

0.8

(5)
Eccentricityr̈ = F grav + F dragr は微惑星の位置を表す。重力項は次のように太陽
重力と木星重力が含まれる。

GM⊙ r GMJ (r − r J ) GMJ r J
F grav = −
−
−
(6)
|r|3
|r − r J |3|r J |3

0.6
0.4
0.2
0.0

1

2

3

rJ は木星の位置を表し、式 (6) の 3 つ目の項は太陽を原点としているために現れる。また、ガス抵抗力

4

5

6

Semimajor-Axis (AU)

7

8

9

10

図 2: 3 体の微惑星の軌道進化の様子。黄、緑、黒のは次のような表式である (Adachi et al. 1976)。

プロットはそれぞれ異なる微惑星軌 道の時間進化を、
Cd πd2 ρg vrel v rel(7) 赤の点は計算を終了した時点での軌道を表す。灰色
F drag = −2mの領域は小 惑星帯を模しており、破線は近日点距離
ここで Cd は抗力係数、v rel はガス円盤の回転速度にが 1 AU の軌道を、青の実 線はヤコビ積分の値が初対する微惑星の相対速度、m、d はそれぞれ微惑星の期値と等しくなる軌道を表す。
質量、半径である。ガス密度 ρg については
(
)
Σg
z2
図 2 は軌道長半径-離心率空間での 3 体の軌道進化
ρg = √
exp − 2(8)HπHの様子を示しており、6 AU のところに木星が軌道をとした (e.g., Kobayashi 2015)。ガス面密度 Σg とガ 持つ。これと同じ微惑星の軌道進化を遠日点距離-離ススケールハイト H は赤道面上の半径 R を用いて心率空間で描画したものを図 3 に示す。図 2 も図 3 も次の表式を採用した (Raymond and Izidoro 2017)。
)−1
(
R
g/cm2
(9)
Σg = 4000
1 AU
H = 0.05R(10)

Eccentricity

0.8

全ての微惑星の内部密度は 1.5 g/cc とし、半径は

d = 100 km を仮定しているので m ≃ 6 × 10

21

gと

0.6
0.4
0.2
0.0

なっている。微惑星は 500 体を同時に計算し、輸送

1

2

3

4

5

6

Aphelion (AU)

7

8

9

10

過程に注目しているため、内側に輸送され、離心率が 0.1 以下になった微惑星の計算は止めている。そ図 3: 3 体の微惑星の軌道進化の様子を横軸を遠日点

れ以外の微惑星については 10 万年間の軌道進化の計

距離で示す。色や線の対応は 図 2 と同じである。

算を行った。
解析の方法としては、軌道進化の過程について、軌同じ微惑星の結果を示しているが、図 3 では木星重道長半径 a や離心率 e、遠日点距離 p = a(1 + e)、近力の影響について見えやすくなっている。これらの日点距離 q = a(1 − e) のような軌道要素に対する進結果により、ヤコビ積分を保存しながら木星重力に化を見ることで、太陽系の内側、雪線の内側に輸送

より軌道進化する段階とガス抵抗によって離心率がされる微惑星の軌道進化の特徴を見出す。

低下していく段階とがあるのがわかり、今回のパラメータでは小惑星帯に入るような輸送が主要となっていることも見てとれる。

4

Results
今回は計算した微惑星 500 体のうち、現在、雪線として考えられている 2.7 AU より内側に輸送され

31

2022 年度 第 52 回 天文・天体物理若手夏の学校

図 4: 本研究でわかった軌道進化の概略図。

5

ガス抵抗によって遠日点が木星から遠ざかり、内側にDiscussion

移動することで木星重力の影響を受けなくなり、ガ
図 2 から分かるように、ヤコビ積分が一定の線にス抵抗のみによって離心率が減少して安定な軌道に沿って離心率が上昇し、ガス抵抗によって軌道長半

落ち着いていく。さらに、雪線の内側にまで輸送さ径と離心率が減少することで輸送が起こる。しかし、遠日点距離と離心率の変化で見ることによってよりれる際、木星重力散乱による離心率の上昇とガス抵抗による離心率の低下が混在する軌道進化段階では、輸送過程での定性的な軌道進化がわかるようになった。木星重力によって離心率が上昇する段階では、ヤコビ積分を保存しながら進化するが、遠日点距離も離心率がほとんど一定で進化することがわかり、その値が予測できれば主要な輸送先を推定できるという示唆が得られた。離心率が一定となる段階についほとんど変化せずに進化することがわかる。さらに、内側の軌道では木星重力の影響を受けるのが遠日点ての解析は今後も行い、将来的には C 型小惑星の移
住や地球への水供給を可能にする輸送のパラメータであるため、重力散乱による離心率上昇とガス抵抗を決定していく予定である。

による離心率減衰がおおよそ釣り合い、離心率がほぼ一定で進化する段階を経由することもわかった。そして遠日点が木星重力の影響を受けにくい、木星かAcknowledgement

ら 3rH -4rH 程度離れるとガス抵抗だけで進化することになる。つまり、重力散乱での進化からガス抵抗

本講演を行うにあたり、共同研究者である小林浩による進化へ遷移する離心率ほぼ一定の段階が解明

助教をはじめ、理論宇宙物理学研究室の皆様には多でき、あるパラメータに対して釣り合う離心率がわくの助言や指導をいただき大変お世話になりました。

かれば、遠日点が木星から十分離れた点からのガス

また、このような研究発表の機会を設けてくださっ

抵抗のみによる軌道進化だけで、主要な輸送先が予た夏の学校事務局の皆様に感謝申し上げます。

測できることになる。

Reference
6

Conclusion

Adachi I., Hayashi C., Nakazawa K., 1976, PThPh, 56,
1756. doi:10.1143/PTP.56.1756

木星形成後における氷微惑星の雪線の内側への輸

Kobayashi H., 2015, EP&S, 67, 60. doi:10.1186/s40623015-0218-y

送過程を理解するために、ガス抵抗を考慮した制限
三体問題の軌道進化について N 体シミュレーションを用いて調べたところ、遠日点距離の進化が輸送のMakino J., Aarseth S. J., 1992, PASJ, 44, 141
Morbidelli A., Lunine J. I., O’Brien D. P., Raymond S. N., Walsh K. J., 2012, AREPS, 40, 251.
doi:10.1146/annurev-earth-042711-105319

理解に重要であることがわかった。本研究で理解が進んだ遠日点の変化を含む微惑星の軌道進化の様子を図 4 に示す。太陽系の内側への輸送は、木星重力によって遠日点がほとんど変化せず離心率が上昇し、Raymond S. N., Izidoro A., 2017, Icar, 297, 134.
doi:10.1016/j.icarus.2017.06.030

32

——–index へ戻る星惑 a07

微惑星衝突に伴う固体蒸発量とその状態方程式依存性
名古屋大学大学院 理学研究科
宮山 隆志

33

2022 年度 第 52 回 天文・天体物理若手夏の学校

微惑星衝突に伴う固体蒸発量とその状態方程式依存性
宮山隆志 (名古屋大学大学院 理学研究科)

Abstract
惑星は衝突合体を繰り返しながら成長すると考えられているため，惑星形成過程において天体衝突は普遍的な現象である. 中でも，高速度の衝突による温度上昇は物質の蒸発を引き起こす．さらに蒸発が生じる基準は物質によって異なるため，衝突蒸発によって惑星の組成編成を引き起こす現象である. つまり衝突蒸発現象の理解は，多種多様な惑星の起源を解明する物理的素過程として重要である. そこで本研究では数値シミュレー
ションを用いて，天体衝突による蒸発計算を行い，蒸発現象についての理論的な解明を行う．衝突によって生じた圧縮波の伝搬速度が物質の音速を超える様な場合は，衝撃波が発生し，不可逆的な加熱圧縮によりエントロピーの上昇を引き起こす. 衝撃波によって加熱圧縮された物質は熱力学第一法則に従い，等エントロピー
的な膨張によって急速に冷えていく. つまり衝撃圧縮によって，物質の沸点を超えるようなエントロピーの上昇が引き起こされ た場合には物質の蒸発が生じる．したがって，衝突によって生じた衝撃波が引き起こす
エントロピーの増加を求める事によって蒸発量を求めることができる．しかし，その蒸発量は用いる状態方
程式に依存することが知られている．そこで蒸発を適切に評価するために，最低限の物理を失わず，なおかつ高速計算に適用可能な状態方程式が必要とされる．本研究ではまず，惑星衝突分野で広く用いられている

Tillotoson EOS のパラメータ改善を行い，衝突蒸発シミュレーションをおこなった．その上で，より洗練された状態方程式である ANEOS を用いて衝突蒸発シミュレーションを行ない，結果の比較をおこなった．その結果，Tillotson EOS による計算結果が ANEOS による結果と整合的であることが確認できた．

1

限の固体物理が表現できる状態方程式だが，計算時Introduction

間も短い．
惑星形成過程において天体衝突は普遍的な現象でこの 2 つの状態方程式が実装されている，数値衝あり，その衝突速度は被衝突天体の脱出速度以上の突計算コード iSALE-2D (e.g., Amsden et al., 1980;

高速度である．このような高速度での衝突によってIvanov et al., 1997; Wünnemann et al., 2006) を用

衝撃波が生じ，衝撃波による温度上昇は固体物質のいて衝突シミュレーションを行ない，蒸発量の比較

蒸発・溶融を引き起こす．この過程は，主に衝突シをおこなった．衝撃に伴う衝撃波の発生により物質

ミュレーションにより調べられるが，固体蒸発量はの沸点を超えるようなエントロピーの上昇が引き起仮定する状態方程式に強く依存する (e.g., O’Keefe &こされると蒸発が起こる (Ahrens & O’Keefe, 1972).

Ahrens, 1982; Pierazzo & Melosh, 2000; Kraus et

そのため，エントロピーの上昇を精度良く解けるこal., 2011)．

とと，蒸発の条件を適切に評価する状態方程式が重要

本研究では，衝突蒸発を適切に評価するために，最低限の物理を失わず，なおかつ高速計算に適用可能である．簡易な Tillotson EOS でも適切にパラメー
タ設定を行えば，エントロピー上昇は比較的精度良な状態方程式を調査した．そのため，惑星衝突分野で広く用いられている２つの状態方程式を用いた．一く解けた．一方，蒸発判定が簡易なため蒸発量は数倍程度ずれることもあった．衝突速度が音速に比べてつ目の ANEOS (Thompson & Lauson, 1972) は数十分大きい場合は Tillotson EOS の結果は，ANEOS

多くのパラメータを用いて室内実験を詳細に再現での場合の衝突蒸発量とファクター２の精度で整合的きるようにチューニングされているが，長い計算時な結果が得られることが分かった．

間も必要となる．一方，もう一つの Tillotson EOS(Tillotson 1962) では少ないパラメータによって最低342022 年度 第 52 回 天文・天体物理若手夏の学校

2

体と液体の相境界を再定義した上で蒸発・溶融量のMethods

計算を行なった．
ここでは本研究で用いた計算手法と，天体衝突によって生じる衝撃波について説明する．以下では衝突

103

Compressed State

天体をインパクター，被衝突天体をターゲットと呼ぶ．

Internal Energy, E (MJ/kg)

また，本研究では科学利用に供するために公開されている数値衝突計算コードである iSALE-2D(e.g., Ams-

den et al., 1980; Ivanov et al., 1997; Wünnemann etal., 2006) を用い，固体物性を反映する状態方程式は
ANEOS (Thompson & Lauson, 1972) と Tillotson
EOS (Tillotson 1962) を用いた．また本研究では岩石天体同士の衝突を想定し，物質の組成は玄武岩と
した．

2.1

ANEOS

Complete

102 Vaporization
Partial101

Ecv

Vaporization100/ 0 <1

10 10

2000

Eiv
/ 0 >1
4000

6000

Density, (kg/m3)

8000 10000

ANEOS による玄武岩 (Pierazzo et al., 2005) の相図を以下に示す (図 1)．図 1 において赤線で示され

図 2: Tillotson EOS による相図の表現ている曲線は Hugoniot 曲線と呼ばれ，衝撃波面前後の保存則である Rankine-Hugoniot 条件と状態方程ここで ρ0 は標準密度である．このように Tillotson

式から導かれる．衝撃波圧縮はこの曲線上の一点へ

EOS は物質を圧縮状態と膨張状態に分ける．このうの不連続変化で表現される．圧縮後の断熱膨張は等

ち圧縮状態 (固体) は以下の解析表現をとる．
][bρu + Aµ + Bµ2
P = a+
u/(u0 η 2 ) + 1

エントロピー変化であるため，断熱曲線は S-P 平面
上 (図 1) で横軸への垂線で表される．

Pressure (GPa)である．a, b, A, B, ρ0 , u0 は物質ごとに固有のパラメー

Hugoniot

101

S+L

10 1
10 3

Solid

10 5
10 7

ここで u は比内部エネルギー，η = ρ/ρ0 , µ = η − 1

Basalt phase diagram

103

ターである．式 (1) の第 1 項は熱的な圧力を表し，第Super critical fluid

LiquidL+V2，第 3 項は非熱的な圧力を表す．また高エネルギー
極限 (u >> u0 ) では理想気体の状態方程式 P ∼ aρu

Vaporに漸近する．

611Pa

0

2

(1)

4

6

Specific Entropy (kJ/K/Kg)

8

10

2.3
図 1: ANEOS による相図の表現

計算セットアップ

本研究では岩石天体同士を衝突させた．初期条件はインパクターとターゲット共に標準密度 ρ0 = 2.86 ×

103 kg/m3 ，インパクターの衝突速度を 14〜70km/s

2.2の範囲で変化させて蒸発・溶融量を求める．また衝突

Tillotson EOS

シミュレーションによる衝撃圧縮の計算精度は分解能
Tillotson EOS による相図を以下に示す (図 2)． に強く依存する．計算分解能は衝撃波圧縮による圧

Tillotson EOS では ANEOS よりも単純な相図の定力上昇の計算精度が 10 ％以内となるよう (Pierazzo

義が用いられており，液相も定義されていない．本研

et al., 2008)，インパクター半径を 20 メッシュで分

究では ANEOS と整合的な結果が得られるよう，気解している．

35

30

6

25

4

20
15
10

Basalt (vimp = 40 km/s, t/ts = 10.0)
TillotsonEOS

100

ANEOS
Internal Energy, E (MJ/kg)

Internal Energy, E (MJ/kg)

2022 年度 第 52 回 天文・天体物理若手夏の学校

2
0
2
4

5

6

0

8

6

4

2

0

2

Radial distance, R/Rp

4

80

ANEOS
TillotsonEOS

60
40
20
0
1000 2000 3000 4000 5000 6000 7000 8000

6

Density, (kg/m3)

図 3: 衝突速度 40km/s の時の t = 10ts = (2Rp /vimp ) までの時間発展を示している．Rp はインパクターの半径，vimp は衝突速度である．左図は左右に同時刻での Tillotson EOS と ANEOS の内部エネルギー分布を示しており，初期位置 (±Rp , −Rp ) のトレーサ粒子 (白三角) の熱力学状態の時間発展が右図である．右図において，赤線，青線はそれぞれの状態方程式による Hugoniot 曲線を示している．

3

衝突速度と蒸発・溶融量の関係を図 5 に示す．これ

Results
衝突速度 40km/s の計算結果を図 3 に示す．図中，右図から分かるようにパラメータの調整によって定義されておらず，蒸発・溶融量を同時に求めることができなかったが，相境界となる内部エネルギー

Tillotson EOS による衝撃波圧縮と断熱膨張の計算結果が，ANEOS と整合的であることが確認できた．また相図の定義の違いにより Tillotson EOS は ANEOS

まで Tillotson EOS の相図では液相と気相は同時にの値を適切に選ぶことで，Tillotson EOS の簡易な相図を用いても ANEOS と誤差数％の精度で蒸発・溶融量を求めることができた．

よりも早く気体となってしまうため，エジェクタが早く出ることがわかった．しかし，t = 70ts まで計算を行ない十分に膨張をさせた結果，図 4 が示すよ
うにエジェクタの総量は数％の誤差に収まっている
ことが分かった．

図 5: 衝突速度と蒸発・溶融量の関係

図 4: エジェクタ総質量の時間変化

362022 年度 第 52 回 天文・天体物理若手夏の学校

4

Conclusion

Melosh, H. J. (1989) Impact cratering: A Geologic Process. Oxford University Press, 245 pp.

天体衝突によって生じる蒸発・溶融は主に数値シ
ミュレーションによって取り扱われるため，最低限の物理を失わずに，なおかつ高速計算が可能な状態方程式が求められる．そこで本研究では，天体衝突によって生じる蒸発・溶融を適切に評価するための状態方程式を調査するため数値衝突計算コード iSALE-2D を用いてシミュレーションを行なった．その中で非常に単純な解析表現と相図で定義される Tillotson EOSのパラメータを適切に選ぶことで，より洗練された
状態方程式モデルである ANEOS と整合的な結果を得ることができた．

Pierazzo, E., Artemieva, N., & Ivanov, B. 2005, in Large
Meteorite Impacts III (Geological Society of A)
Pierazzo, E., & Melosh, H. J. 2000, Icarus, 145
Pierazzo, E., Artemieva, N., Asphaug, E., et al. 2008,
M&PS, 43, 191
Thompson, S. L., & Lauson, H. S. 1972, Improvementsin the Chart D Radiation-Hydrodynamic Code. III:
Revised Analytic Equations of State, Albuquerque,
New Mexico: Sandia National Laboratory. 1972.
Technical Report SC-RR-71-071
Tillotson, E. (1962), Metallic Equations of State For
Hypervelocity Impact, General Atomic Report GA3216. 1962. Technical Rep

Acknowledgement
iSALE の 開 発 者 で あ る Gareth Collins,

O’keefe,J.D. & Ahrens,T.J. (1982), J.Geophys. Res., 87,
6668

Kai

Wünnemann, K., Collins, G. S. & Melosh, H. J. (2006)
A strain-based porosity model for use in hydrocode.
simulations of impacts and implications for transient
Elbeshausen, Tom Davison の各氏に感謝致します．
crater growth in porous targets. Icarus, 180, 514– 527.

Wünnemann, Boris Ivanov, H. Jay Melosh, Dirk

また，本講演を行うにあたり，共同研究者である小林浩助教をはじめ理論宇宙物理学研究室の皆様には多くの助言をいただき，大変お世話になりました．このような研究発表の場を設けてくださった夏の学校
事務局の皆様に感謝申し上げます．

Reference
Ahrens, T. J., & O Keefe, J. D. 1972, Moon, 4, 2
Ahrens, T. J., & Okeefe, J. D. (1977), in Impact andExplosion Cratering: Planetary and Terrestrial Implications, ed. D. J. Roddy, R. O. Pepin, & R. B.
Merrill, 639–656
Amsden, A., Ruppel, H. & Hirt, C. (1980) SALE: A
simpliﬁed ALE computer program for ﬂuid ﬂow at
all speeds. Los Alamos National Laboratories Report,
LA-8095:101P.
Ivanov, B. A., Deniem, D. & Neukum, G. (1997) Implementation of dynamic strength models into 2-D hydrocodes: Application for atmospheric breakup andimpact cratering. International Journal of Impact Engineering, 20, 411–430.
Kraus, R. G., Senft, L. E., & Stewart, S. T. (2011),
Icarus, 214, 724

37

——–index へ戻る星惑 a08

講演キャンセル38

——–index へ戻る星惑 a09近接遭遇による惑星系の軌道進化と AMD を用いた安定性評価
東京大学大学院 理学系研究科
神原 祐樹

39

未提出

40

——–index へ戻る

星惑 a10

トランジット測光による系外惑星 KELT-9 b のスピン傾斜角の推定東京大学大学院 総合文化研究科
河合 優悟

41

2022 年度 第 52 回 天文・天体物理若手夏の学校

トランジット測光による系外惑星 KELT-9 b のスピン傾斜角の推定
河合 優悟 (東京大学大学院 総合文化研究科)

Abstract
地球や他の太陽系惑星の公転軸は、主星である太陽の自転軸と 10◦ 以内で揃っている。一方で、太陽系外には、主星の自転方向と反対向きの逆行軌道や、主星の北極と南極を通過する極軌道を持つ惑星が続々と見つかってきている。この主星の自転軸に対する惑星の公転軸の傾きの角は「スピン傾斜角」と呼ばれる。近年の多様なスピン傾斜角の発見は、太陽系惑星とは全く異なる惑星の形成や軌道進化の過程の存在を示唆している。よって、これを測定することで「どのような惑星は軌道がズレるのか」を紐解くことができ、惑星の軌道が進化する仕組みの解明につながる。
これまでの研究から、主星の温度が約 6,250K 以下（低温星）の系ではスピン傾斜角は 0 度に近いことがわかっている。これは、低温星は表面に対流層を持ち、その対流層と惑星の潮汐の相互作用により時間とともに軌道が揃ってしまうためであると考えられている。反対に、表面に対流層を持たない 6,250K 以上（高温星）の系では、惑星の形成時の軌道の情報が揃うことなく残っている。
そこで本研究では、トランジット系外惑星探査衛星 TESS の観測データから、高温星まわりの系外惑星

KELT-9 b のスピン傾斜角を 85.1◦ ± 1.2◦ と推定し、その極軌道を確認した。現在までにスピン傾斜角が測定された 57 個の高温星周りの惑星のサンプルからは、スピン傾斜角の分布はランダムではなく、0◦ （巡行軌道）と 90◦（極軌道）周辺に分布することが示唆されているが、本研究はこの仮説を支持する結果である。

1

Introduction
惑星はどのように形成されるのか。天文学はこの謎を解明するため、これまで 5,000 個 以上の惑星を
太陽系の外に発見してきた。観測技術の進歩によりこれらの惑星の詳細な調査が可能となり、惑星の形
成や進化の歴史を解明する手がかりが得られるようになった。中でも、主星の自転軸に対する惑星の公転軸の傾きの角であるスピン傾斜角 ψ は、惑星が形成時に受けた力学的な影響の情報を保存している。この情報から、「どのような惑星は軌道がズレるのか」

図 1: スピン傾斜角の概念図を紐解くことができ、これは惑星の軌道が進化する
仕組みの解明につながる。

で定義される量である（図１）。

定量的には、スピン傾斜角 ψ は、観測者の視点か

これまでの研究から、主星の温度が約 6,250K 以下
ら主星の自転軸がどの程度起き上がっているかを指 （低温星）の系ではスピン傾斜角は 0◦ に近いことがす「自転軸傾斜角 i∗ 」、同様に主星の自転軸がどの程 わかっている。(Albrecht et al, 2021) これは、低温
度起き上がっているかを指す「軌道傾斜角 iorb 」、そ星は表面に対流層を持ち、その対流層と惑星の潮汐して観測者の視線と直交する天球面へのスピン傾斜の相互作用により時間とともに軌道が揃ってしまう

角の射影である「射影スピン傾斜角 λ」によってためであると考えられている。反対に、表面に対流

cosψ = cosi∗ cosiorb + sini∗ siniorb cosλ.

(1)

層を持たない 6,250K 以上（高温星）の系では、惑星42

2022 年度 第 52 回 天文・天体物理若手夏の学校の形成時の軌道の情報が残っている。その例として、
高温星まわりでは、主星の北極と南極を通過する極
軌道や、主星の自転方向と反対の逆行軌道の惑星が観測されている。(eg. Narita et al, 2009)
また、現在までにスピン傾斜角が測定された 57個の高温星周りの惑星のサンプルからは、スピン傾斜角の分布はランダムではなく、0◦ （巡行軌道）と90◦（極軌道）周辺に分布することが示唆されている。
(Albrecht et al, 2021) サンプル数を増やすことで、どのような系は巡行軌道を保ち、どのような系は極軌道へと進化するのかを検証することができる。
本研究では、トランジット系外惑星探査衛星 TESSの観測データから、高温星まわりの系外惑星 KELT-

図 2: 軌道による惑星のトランジットと観測データの

9 b のスピン傾斜角を 85.1◦ ± 1.2◦ と推定した。こ概念図（上：軌道が揃っている場合、下：軌道が 90◦

れは Ahlers et al. (2020) の 87+10−11 と非常によく同ズレている場合）

意する結果である。太陽系惑星のスピン傾斜角は全て 10◦ 程度であるため、KELT-9 b の極軌道は、太のとおりである。星は、自転の遠心力により赤道方向

陽系惑星とは全く異なる形成や軌道進化の過程を示に膨らんでいる。膨らんだ赤道付近は密度が下がり、唆している。本稿では２章で推定に採用した手法と温度が下がる。明るさは温度に比例するため、結果使用したデータについて述べる。その後、３章で解的に星は赤道ほど暗く、極点ほど明るくなる。この主析結果について述べ、4 章で今回の結果が示唆する星の自転と温度の温度の関係は、以下の von Zeipel

KELT-9 b の形成過程について議論する。

の定理 (von Zeipel, 1924) で表される。

(

2

T = Tpole

Methods

g

)β

gpole

(2)

ここで、T と g はそれぞれ主星の任意の地点での表

スピン傾斜角を推定するために、本研究ではトラ

ンジット測光と呼ばれる手法を採用した。トランジッ 面温度と重力、Tpole と gpole は極点の表面温度と重トとは、惑星が主星の前を通ることをさす。（図２） 力である。重力減光の程度はこのように表面重力にこの時主星の明るさを観測（測光）すると、惑星が依存するため、星の自転が速ければ速いほど顕著に隠す分主星は暗くなって見える。この時の観測デー現れる。β は定数であり、0.25 とされることが多い。
この主星の重力減光を組み込んだ惑星のトランジッ

タの形状からスピン傾斜角を推定する。
観測データの形状は惑星の軌道に依存するため、ストのモデル (Barnes, 2009) を用い、マルコフ連鎖モ

ピン傾斜角が推定できる、これは、トランジットをンテカルロ法 (MCMC) を実行することで、モデル通して遮られる主星の光の量が一定ではないことにパラメタの一つであるスピン傾斜角の推定を行った。

よる（星は赤道ほど暗く、極点ほど明るい）。たとえ

この際、Python のオープンソースライブラリである

ば、軌道にズレがないトランジット（暗い赤道のみ

Pytransit (Parviainen, 2015) を使用した。

を通過）では、図２右上のようになる。他方で、軌解析に使用したデータは 2019 年 7 月から 8 月にか道が大きくずれているトランジットでは、図２右下

けてトランジット系外惑星探査衛星 TESS の Sectorのようになる。この観測データの形状の違いからス

14 で観測された KELT-9 b のトランジット 18 回分

ピン傾斜角を推定する。

である。このデータは一般公開されているものであ

なお、主星が星は赤道ほど暗く、極点ほど明るくな

り、doi:10.17909/hfat-m449 にて入手可能である。

る効果は重力減光と呼ばれている。この仕組みは次432022 年度 第 52 回 天文・天体物理若手夏の学校

3

Results

表 1: 先行研究 (Ahlers et al. 2020) と本研究で推定

KELT-9 b は主星である KELT-9 を周期約 1.48 日

された KELT-9 b の主要なパラメタの比較で公転する巨大ガス惑星である。KELT-9 は A 型星

パラメタ

本研究

Ahlers et al. (2020)

で、表面温度は約 10,000K と非常に高く、軌道の状

i∗ (deg)

38.3 ± 2

38+8−7

況は潮汐の影響を受けずに残っていると考えられる。

iorb (deg)88.39+1.01−0.72

86.7 ± 0.3

観測データとモデルを用いた解析結果を図３に示す。

λ(deg)-84.8 ±1.4

-88 ± 15

重力減光を加味したモデル（赤線）は、加味していψ(deg)

85.1 ± 1.2

87+10
−11

ないモデル（黒点線）に比べてより良いフィットになっていることが残差からもわかる。青い点が観測データであり、オレンジ色の点は可視性のためにこのデータを 4 分間隔で bin したものである。

図 4: KELT-9 b の軌道のイメージ図（灰色点線の球体が主星の KELT-9、黒色の点が KELT-9 b。三角はそれぞれ主星の自転方向と惑星の進行方向。）

4Discussion
KELT-9 b のような極軌道を描く惑星は太陽系には存在しないが、太陽系外では発見され始めており、

KELT-9 b が稀有な存在ということではないようである。現在までにスピン傾斜角が測定された 57 個の高温星周りの惑星のサンプルからは、スピン傾斜角の分布はランダムではなく、0◦（巡行軌道）と 90◦（極

図 3: KELT-9 b の測光データの解析結果

軌道）周辺に分布することが示唆されているが、本このデータの解析から推定された KELT-9 b の主研究はこの仮説を支持する結果である。しかし、こ要なパラメタと、その先行研究 (Ahlers et al. 2020)との比較が表１である。本研究では、KELT-9 b のス
ピン傾斜角を 85.1◦ ±1.2◦ と推定した。これは Ahlers
これらのパラメタから観測者視点で KELT-9 の系
星はより明るい極点を先に隠し、その後より暗い赤道付近を通過している。これは左肩下がりで非対称になっている観測データと直感的にも一致する。

ついてはまだ議論が活発に行われている。

KELT-9 b のような主星の周りを短周期で公転す

et al. (2020) と非常によく同意する結果である。
を図示したものが図 4 である。トランジット中に惑

れらの軌道がズレた惑星に関して、その形成過程にる（主星に非常に近い）巨大ガス惑星は、非常に高温になることから「ホットジュピター」呼ばれている。現在の惑星形成理論では木星サイズのような巨大な惑星は主星から 0.1AU 以内で形成できないた
め、KELT-9 b に代表されるようなホットジュピター44

2022 年度 第 52 回 天文・天体物理若手夏の学校
はより外側の軌道で形成された後、内側に移動して

5

Conclusion

くると考えられる。この軌道進化の過程で、ホットジュピターの軌道にズレが生じると予想されている。

(Fabrycky and Tremaine, 2007)
軌道進化の段階でズレを生むシナリオの一つは、第今回は系外惑星 KELT-9 b のスピン傾斜角を85.1◦ ± 1.2◦ と推定した。極軌道が広い宇宙ではありふれているかもしれないということは、太陽系が宇宙の中で標準的な系ではない可能性を示唆してい三の天体の重力的な摂動によるものである。特に軌る。5,000 以上発見されている系外惑星の中で、スピ道が傾いた伴星や惑星が軌道の外側に存在するケーン傾斜角が判明している惑星は一握りであり、その中スは古在-Lidov 機構と呼ばれている（図 5）。古在-でもトランジット測光と重力減光から推定された惑

Lidov 機構ではこれらの第三の天体の重力的な影響

星は 10 個程である。サンプルを増やすことは、太陽によって、内側の惑星のスピン傾斜角と離心率が周期的なトレードオフの関係となり、大きく変化する。

系のイメージにとらわれない惑星の軌道の進化の仕組みを紐解くことにつながる。現在の観測精度で重軌道進化自体のシナリオの一つは high eccentricity

migration と呼ばれており、離心率の高い惑星が主星に近づいた際に、潮汐力によって軌道が円形になる

力減光の効果を捉えられる系が限られているが、今後 James Webb 宇宙望遠鏡や他の望遠鏡によってさらに高精度な観測が可能となれば、さらに多くの惑ことでホットジュピターが形成されると考えられている。古在-Lidov 機構によって離心率とスピン傾斜角が大きくなった際に high eccentricity migration に星のスピン傾斜角をこの手法を用いて測定することが可能である。

よって軌道が円形になれば、KELT-9 b のような軌道のズレたホットジュピターが形成される。この機

Reference

構では、第三の天体の条件（質量など）によっては巡行軌道の惑星に加えて極軌道の惑星も多く生成すAhlers et al. 2020, AJ

ることが可能である。

Albrecht et al. 2021, ApJL
Barnes, J. W. 2009, ApJ
Fabrycky, D and Tremaine, S. 2012, ApJL
Narita et al. 2009, PASJ
Parviainen, H. 2015, MNRAS
Rogers et al. 2012, ApJL
von Zeipel, H. 1924, MNRAS

図 5: 古在-Lidov 機構の概念図
そして、古在-Lidov 機構などは惑星の軌道がズレることを想定しているが、惑星の公転軸ではなく、主星の自転軸がズレる可能性も考えられる。特に対流層を核に持ち、放射層を表面に持つ高温星は、この二つの層の境界で内部重力波が発生し、これが星表面で折り返す。この時この内部重力波が不均一に打ち消されることで光球がランダムに傾くことがある
ことが指摘されている。(Rogers, 2012) KELT-9 のような高温星においては、こうした軌道進化のシナリオも考えられる。

45

——–index へ戻る

星惑 a11SETI につながる月面反射電波観測
熊本大学大学院 自然科学教育部理学専攻物理研究科
学コース
上内 雅暁

46

2022 年度 第 52 回 天文・天体物理若手夏の学校

SETI につながる月面反射電波観測
上内 雅暁 (熊本大学大学院 自然科学教育部理学専攻物理研究科学コース)

Abstract
我々人類のような知的生命は、地球以外の惑星に存在するのだろうか。知的生命が存在すれば、その惑星には文明が存在し人工電波を発信している可能性がある。 地球外知的生命探査 (SETI ) の中で最も 大規模に行われている方法は、電波を地球上の電波望遠鏡で受信し、それを解析することで地球外知的生命から発信されたものを探査するというものである。
しかしながら、地球外文明が放出される人工電波の特徴は未知である。そのため、まずは知的生命が存在する惑星からの人工電波をどのように観測し、解析すれば検出できるかを検証する必要がある。月面反射電波観測はその 1 つで、これは現時点で唯一知的生命の存在が知られている地球からの人工電波のうち、月に反射して戻ってきた月面反射人工電波を電波望遠鏡で観測するという手法である。それを疑似的に知的生命が存在する惑星からの人工電波とみなすことで、知的生命から放射された人工電波の観測につなげていく。
本研究で観測する人工電波には、月面で反射して戻ってきた人工電波と望遠鏡に直接入射する人工電波の 2
種類があるが、両者を区別することは難しい。そのため、月面方向 (on-moon) と月面以外の方向 (oﬀ-moon)から検出される人工電波強度の差を取れば月面反射人工電波のみを検出することができる。今回用いるデ ー
タは水沢 VLBI 観測所で 2 日間観測を行い得られた人工電波である。現在、1 日目のデータを解析中でこの日は K-band で観測した。on-moon と oﬀ-moon を同時に観測し、それぞれを 1 μ秒で時間積分しデータを解析した。その差を図示することで月面反射人工電波の有無を探る。

1

Introduction

1.1

SETI

Interferometry) を実行できるまでタイムラグが生じその信号の位置を特定する手段をどのように実現するかという問題が生じる。そのため、送信帯域幅が狭い

宇宙空間において、我々地球人以外で知的な生命は信号の特定という困難さが伴う。このような、SETI

存在するのだろうか？地球外知的生命探査 (Search forのために必要な周波数や天体の位置に関する膨大なパ

Extra-Terrestrial Intelligence) は略して SETI と呼ば ラメータをサーベイするため現在では SKA(Square
れる。地球外知的生命が存在する惑星では文明が発 Kilometre Array) での観測が必要である。SKA と達しており、そこでは地球のように人工電波を用いては、国際的な取り組みにより建設が進んでいる大型通信を行っていると考えられる。1959 年の Giuseppe

電波干渉計であり、高感度・広視野の観測が可能でCocconi と Philip Mor-rison は地球外文明が他の惑ある。そのため惑星からのわずかな電波観測に大い星との交信を目的に人工電波を発するとすれば、そに貢献できるであろう。[1]れは中性水素の出す 21cm 電波であろうと述べた。これを受けて 1960 年に F.Drake はオズマ計画を実行

1.2

し地球外文明の人工電波の検出を試みたが、見つか人工電波らなかった。その後、様々な SETI の計画が進めら

地球外知的生命が存在する惑星では文明が発達しれたが地球外文明由来の人工電波は検出されていなており、そこでは人工電波が用いられている可能性い。

がある。惑星から宇宙空間に漏れ出る人工電波の種

また、仮に地球外知的生命体からの信号がある方向類として、他の惑星との交信が目的で用いられるもから受信できたとしても、VLBI(Very Long Baselineの、日常生活を送る上で必要不可欠な人工電波の 2472022 年度 第 52 回 天文・天体物理若手夏の学校
種類が考えられる。それを地球上の電波望遠鏡で受2.1.2

信できれば地球外知的生命の存在の証明につながる。
しかし、宇宙空間から飛来してくる電磁波は地球

月面反射人工電波観測

先行研究である Earthshine では地球が反射する太陽光に着目していた。本研究ではそれを応用して、地外文明からの人工電波だけではなく恒星由来のもの

球から溢れ出る電波に着目することにする。月面反や、宇宙空間におけるあらゆる物理現象由来のもの

射人工電波観測は地球外文明からの人工電波の特徴

等様々である。その中から地球外文明からの人工電波のみを検出したいのだが、その特徴について我々は未知である。そのため、まずは地球外文明からの人
工電波の特徴を知るとともに、地球上の電波望遠鏡

や、その観測・解析法を検証する方法の 1 つである。
これは現時点で唯一、我々が知っている知的生命が存在する惑星である地球から漏れ出る人工電波を用いる。地球から漏れ出た人工電波のうち、月で反射

でどのように観測し解析すれば良いのかを模索する必要がある。そこで月面反射人工電波の観測を行う。

して再び地球に戻ってきたものを地球上の電波望遠鏡で観測する。これを疑似的に知的生命が存在する
惑星からの人工電波とみなすことで地球外文明から

2

Methods&Observation

放射された人工電波の観測に繋げていくといった手

2.1

Method

地上の望遠鏡で検出しようという試みは EoR(Epoch

2.1.1

Earthshine法である。このような月面で反射された人工電波をof Reionization) のシグナル検出における正確なキャリブレーションを行う目的でも行われた。[3]

本研究では、地球照 (Earthshine) と呼ばれる現象

また、地球上の電波望遠鏡で検出される人工電波に着目した。[2] これは、地球の表面から反射されたには、月面で反射して戻ってきた人工電波と望遠鏡

太陽光が月面にも反射され再び地球に戻ってくる現に直接入射する人工電波の 2 種類があるが、両者を象である (図 1)。地球照によるシグナルの測定を行区別することは難しい。今回ターゲットとするのはうことで、我々はあたかも地球外から地球を見てい月面から反射されて戻ってくる人工電波のみである

るかのように観測することができる。また、地球照ため、2 種類の人工電波のうち後者を取り除く必要

によるシグナルの解析を行うことで、地球の表面上がある。

の構造を知るだけでなく、地球上の生命からのシグナルの特徴についても知ることができる。この現象を本研究で応用した。

2.2

Observation

今回用いるデータは水沢 VLBI 観測所で 2 日間観測を行い得られた人工電波である。本観測所はマルチビームで 2 方向からの人工電波を同時観測することが可能であるため、1 日目の観測は月面方向 (on-moon) と月面以外の方向 (oﬀ-moon) の 2 方向を同時観測した。前者は月面方向から入射する人工電波の他に、直接電波望遠鏡に入射する人工電波も含まれている。後者はノイズとなる直接望遠鏡に入射する人工電波を観測している。そのため、前者の人工電波の強度から後者からの人工電波の強度を引けば月面から反射された人工電波のみを検出することが
図 1: 地球照のイメージできる。また、2 日目の観測は主に on-moon のみの人工電波を観測し、oﬀ-moon の人工電波は 1 時間ごとに数秒観測した。

48

2022 年度 第 52 回 天文・天体物理若手夏の学校
水沢 VLBI 観測所における観測時間について、1 日ために、周波数ごとにヒストグラムをとった。図 3

目は 2021 年 9 月 28 日 23 時から 9 月 29 日 7 時の 8はそれぞれ積分時間を 1 μ秒とした周波数 92MHz、時間、2 日目は 9 月 29 日 24 時から 9 月 30 日 12 時の360MHz の強度のヒストグラムである。時間 bin 数

12 時間観測を行なった。また観測した周波数帯は、1は 100 万個としているため 1 秒分のデータを見てい日目が K-band(24GHz) で、2 日目は X-band(8GHz)る。図 3 から、on-moon と oﬀ-moon のノイズは χ2である。

分布に従っていることが分かるため、on-moon の方

今回は K-band の観測結果をメインに議論を進めに有意に大きな強度の分布があれば、それはシグナる。また、K-band の帯域にある人工電波には空港

ルの可能性が高い。図のように on-moon だけに現れレーダーが含まれており、これは 1 μ秒のパルス幅ているシグナルがあるが、oﬀ-moon にだけに現れてで約 2〜3 秒ごとに飛来する言われている。月面方向いるシグナルもある。そのため、on-moon だけに現からの人工電波として、これをターゲットとする。

れるシグナルが月面反射人工電波であると断定できなくなった。また、直接望遠鏡に入ってくる人工電
波は on、oﬀ-moon 両方に入ってくるという想定と異3Results&Discussion

3.1

K-bandなる結果が得られた。以上の統計的性質から、今後の研究では on-moon に着目して周期的なシグナルがないか模索することにする。

1 日目の観測で得られた on-moon と oﬀ-moon の電波強度の差をそれぞれスペクトル、ウォーターフォー
ルプロット図で以下に示す (図 2)。図は時間分解能(積分時間) を 1 μ秒とした観測開始から 1 時間 20 分
30 秒後のデータである。直接入ってくる人工電波やノイズとは別に強度が強く検出している電波が求めたい月面反射人工電波ということになる。図 2 のス
ペクトルについて、強度 0 付近で振れているのはノイズだと考えられる。一方で 310MHz 付近等ピークが得られた周波数もあり、これに関しては月面反射
人工電波の可能性がある。しかし、1 μ秒という短い時間でしか見ておらず周期性が分からないため断定はできないため、ウォーターフォールプロット図を用
いて周期性があるかどうか検証することにした。図

3 のウォーターフォールプロット図について、横軸の「bin」は時間の bin 数を表しており 1bin が 1 μ秒分の積分時間を表している。時間 bin 数は 1000 個としているため 1m 秒分のデータを見ていることになる。
しかし、2.2 でも述べたようにターゲットの空港レー

図 2: 観測開始から 1 時間 20 分 30 秒後の K-bandのスペクトル (上) とウォーターフォールプロット図

ダーは 2〜3 秒の周期で検出されるため 1m 秒では非(下)。積分時間はそれぞれ 1 μ秒。

常に短いとわかる。図 3 からもわかるように、1m 秒という短い時間では周期的に現れる大きなシグナルは得られなかった。
また、このように大量のデータを 1 つ 1 つ解析するのは難しい。そのため、on-moon と oﬀ-moon の統計的な性質を見て大きなシグナルがないか検証する

49

2022 年度 第 52 回 天文・天体物理若手夏の学校
はさまざまな時間での解析が必要である。

1 日目、2 日目とも長時間で見た場合の on-moon、oﬀ-moon のスペクトルや統計的な性質などまだ分かっていないことが多いため、今後長時間の解析で研究を進めていき月面反射人工電波検出を断定したい。

図 3: 観測開始から 1 時間 20 分 30 秒後のヒストグ

図 4: 観測開始から 1 時間 30 分 00 秒後の on-moon(青ラム。積分時間は 1 μ秒である。 周波数が 92MHz線)、1 時間 26 分 12 秒後の oﬀ-moon(オレンジ線) のの時 (上 2 つ) と 360MHz の時 (下 2 つ) 赤丸は on-

X-band のスペクトル。積分時間は 10 秒 (1 μ秒分積

moon だけに現れているシグナルで青丸は oﬀ-moon

分した強度を 10 秒分足し上げている) である。

だけに現れているシグナル。

4

Conclusion&Future work

Acknowledgement

今回紹介した 1 日目の K-band の結果については

今回の研究に際し解析に協力していただいた国立

研究当初、様々な時間帯でより多くの観測データを天文台の吉浦さん、観測に協力していただいた水沢

長時間解析する方法が確立できていなかった。その観測所の方々に心より感謝申し上げます。ありがとため、短い時間でしか観測データを見ることができず月面からのシグナルだと断定するに至らなかった。

うございました。

しかし、最近になって解析方法の工夫により多くの観測データを短い計算時間で解析できるようになった。したがって、長時間の解析をすることで大きなReference[1]Andrew P.V.,et al.,2014,Searching for Extraterrestrial Intelligence with the Square Kilometre Array

シグナルが周期的に検出できれば月面反射人工電波の検出の可能性が上がる。

2 日目の観測データに関して、解析途中ではあるが on-moon、oﬀ-moon の電波強度のスペクトルを示す (図 4)。図は観測開始から 1 時間 30 分 00 秒後の[2]Enric Palle & Pilar Montanes-Rodriguez.,2012.,The
Earthshine observations: from climate change to Astrobiology[3]McKinley B.,et al.,2013,The Astronomical Journal,145,23

on-moon と 1 時間 26 分 12 秒後の oﬀ-moon のデータであり、縦軸の強度は 1 μ秒分積分した強度を 10
秒分足しあげたものである。on-moon について強度が強く検出している電波が月面反射人工電波の可能性が高いということになるが、この時間帯では月面
反射人工電波検出に至らなかった。そのため、まず50

——–index へ戻る星惑 b01

原始惑星系円盤の外側の領域における弱い乱流の起源
九州大学大学院 理学府 地球惑星科学専攻
下川 貴史

51

2022 年度 第 52 回 天文・天体物理若手夏の学校

原始惑星系円盤の外側の領域における弱い乱流の起源
下川 貴史 (九州大学大学院

理学府

地球惑星科学専攻)

Abstract
原始惑星系円盤における角運動量輸送の仕組み、そしてこの輸送に伴い強い乱流が発生するかは惑星形成において非常に重要な問題である。近年、ALMA 観測によって円盤の外側ではこれまで行われてきた理論計算から予想される MRI 起源の乱流と比べて弱い乱流 (∼ 0.05 ∼ 0.1cs ) が発生していることが示唆されている。
これらの観測結果は原始惑星系円盤のほとんどが層流にあることを示しているため、これまで理論研究で仮定されていた円盤の初期条件のモデルを見直す必要がある。

Simon et al 2018 では電離度や磁場といった初期条件を変化させて 3 次元の shearing-box でシミュレーショ
ンを行い、観測結果と同様の結果が得られるような初期条件を模索した。このシミュレーション結果によると磁場により強い層流ほとんど近い風が発生している場合であっても局所的に MRI が不安定になる場所が存在し強い乱流が存在することが判明した。また、乱流速度が観測結果の上限値を下回るのは従来の理論研
究で仮定されている電離度及び磁場が低く設定されている時のみであったことを明らかにした。
本講演では上記の Simon et al .2018 をレビューし観測結果およびシミュレーション結果から予想される原始惑星系円盤内の磁場構造及び、それに伴う角運動量の輸送機構についての議論を行う。

1を満たすことが判明した。これは円盤の上層部ではIntroduction

0.1cs ∼ 1.0cs 程度の乱流が発生するという理論研究

原始惑星系円盤が中心星に降着するためには円盤(Simon.et al.2015a) に矛盾する。また、TW-hya にお内のガスが角運動量を失うことが必要である。角運

いても vturb ≤ 0.08cs である (Flaherty et al. 2018)

動量がいかにして輸送されるかは惑星形成を理解すことが知られている。これらの観測結果は円盤内でる上で非常に重要な問題である。こうした角運動量のは強い乱流が発生しておらず、大規模なスケールの層輸送の原因としては磁気回転不安定性 (MRI,Balbus

流が発生していることを示唆していることから、こ& Hawley1998) による乱流、磁気風による角運動量れまで理論研究で仮定されていた円盤の初期条件のの輸送が有力視されている。

モデルを見直す必要がある。

しかしながら、MRI 由来の乱流は原始惑星系円盤では非理想磁気流体効果によりこれらの不安定性は減衰されるため、乱流が起こるのは FUV によりガスが2電離される円盤上層部の表層のみ (Simon etal.2013a)

Methodsである。また、磁気風による角運動量の輸送では

本研究では shearing box 近似を用いて小さな降

Blandford-Payne プロセスと層流トルクにより磁気着円盤の時間進化を非理想磁気流体効果を全て含む

風が駆動されることによりガスが角運動量を失う

3 次元の磁気流体シミュレーションで追った。この際(Suzuki & Inutsuka 2009) ことが示されている。こに用いたコードは ATHENAcode(Stone et al.2008)のように、理論研究では円盤内で乱流と磁気風によるである。そして計算領域は動径、方位角、鉛直方向(x,y,z) 空間において角運動量の輸送が起こるということが示されてきた。 をそれぞれ x,y,z としたときの√しかしながら、HD163296 を ALMA で観測した結 (4H, 8H, 8H)(H = 2cs /Ω) であり、1H あたり 32
果 (Flaherty et al.2015,2017) をモンテカルロ法を用

個の格子を設定した。このような設定で以下の基礎

いて分析すると円盤内の中心から 30au ∼ 100au までの全ての領域で起こった乱流の速度は vturb ≤ 0.05cs

52

2022 年度 第 52 回 天文・天体物理若手夏の学校

図 1: 設定した初期条件
また、N1 = 3.0 × 1021 cm−2 ,N2 = 1.0 × 1024 cm−2方程式に代入して計算を行った。
 ∂ρ+ ∇ · (ρv) = 0


∂t



 ∂ρv
∂t

であり NH は円盤表層の水素の数密度である。
また、円盤表層においては FUV による電離効果を考える必要がある。そのため、このときの両極性拡散+ ∇ · (ρvv + T) = [2v × Ω + 3Ω2 xex − Ω2 zez(T = (P + B 2 /2)I − BB)


 ∂B


= ∇ × (v × B) − ∇ × (ηo J + η (J × B̂) + η J )
A ⊥

H

∂t

効果のエルサッサー数はAmF U V = 2.4 × 1016

(1)

r 3/2
f
−5
10 1au(3)t は時間、ρ は密度、P は圧力、Ω は円盤の角運動量、 ρ0,mid は midplane におけるガス密度で f は全てのvbB は磁場、J は電流、ηO , ηH , ηA はそれぞれ、オー ガス粒子に対するイオンの割合である。
ム散逸、ホール効果、両極性拡散の磁気拡散係数を表 今回のシミュレーションではダストを含まない林円
す。このとき、プラズマベータ β と考慮する電離源

盤モデルを採用し、計算を行った。ダストの影響にを図 1 のように変化させながら得られたシミュレー

関しては Discussion で述べることとする。

ションの結果を ALMA の観測結果と比較した。
基礎方程式で使う磁気拡散係数を求めるためには円盤内の荷電粒子の存在量を求める必要がある。その3

Results

3.1

乱流の速度

ため、化学計算ネットワーク (Bai(2013)) を用いて化学計算を行い、荷電粒子の存在量を求めた。
化学計算の際に考慮する可能性のある電離源は宇宙

表 1 のように初期条件を変化させて計算させると線、X 線、放射性崩壊であり、それらの電離率はそ

乱流速度 δv は図 2、図 3 のような結果が得られた。

れぞれ以下のように表される。

ほとんど場合で |z| ≥ 2H では δv = 0.1cs ∼ cs 程度

ζCR = 1.0 × 10

−17{exp (−Σ/96gcm

−2 −1

s

)

}

の強い乱流が発生することがわかった。また弱い乱

N
N
LX,30
(− H )0.5(− H )0.7 流が発生するかどうかは磁場の強さ (β の大きさ) とζXR = 2.2 [4.0 × 10−11 e N1
+ 2.0 × 10−14 e N2]
rau
電離度に関係があり、弱い磁場 (β が小さい) の時は

ζAl = 1.0 × 10−19

高い電離度が必要という反比例の関係にあることが(2)わかった。また、シミュレーションを行っていないΣ は円盤の面密度であり、LX,30 は X 線の光度を
Lx = 4 × 1029 とおくと LX,30 = LX /1030 である。

β ≤ 103 では磁場が強すぎて MRI 乱流が発生せず層53

2022 年度 第 52 回 天文・天体物理若手夏の学校
流になり (Simon et al.2013a)、β ≫ 106 では MRI

3.2

が活発になるため強い乱流が発生する (Bai & Stone2011) ため、この場合も観測結果を満たさない。こうした結果を考慮すると図 2 のような表を描くことができる。この表から、観測結果を満たすのは赤の破
線と紫の破線に囲まれた領域が観測結果を満たす初期条件であることがわかる。

Zonal ﬂow

ここでは外側領域の midplane 付近でも強い乱流が発生するかどうかを見ていく。今回のシミュレーショ
ン結果では Magnetic wind が卓越しており、大局的には強い乱流が発生していない領域であっても局所的には強い乱流が発生していることが判明した (図4、図 5、図 6)。図 3 は表 1 の R100-B3p-XC-5.25Hを 30 公転周期の時のスナップショットである。この
図を見ると磁場の流束が再分配されているのがわかる。これにより多くの領域で β ≥ 104 になり MRI が活発になっている。103 ≤ β0 ≤ 104 の時は midplaneでも zonal ﬂow が発生しているため、磁場の流束が再分配されて MRI が活発になっている。ここから、層流になっていると考えられている領域であっても局所的には強い乱流が発生していることがわかる。
図 5 は |z| < 0.5H における磁場のエネルギー密度を示している。磁場の強さが時間変動していることが

図 2: 磁場ごとの乱流の強さの分布わかる。図 6 は z ∼ 0.5H − H での 20 から 35 公転
周期の平均の物理量の動径方向ごとの分布を示している。図 5 と図 6 を見比べると時間変化と物理量の変化が一致しているため、ここからも zonal ﬂow が発生していることがわかる。

図 3: 磁場の強さと電離度の関係
この図の色付きの数字は z ≥ 2H での乱流の強さを表している。青色の強さは観測

図 4: 30 公転周期での z 方向の磁場の強さの分布

結果の上限値より大きい乱流の強さ、赤色の数字とその間の観測結果の乱流の下限値となる電離度と磁場の強さの領域の境界、度と磁場の強さの領域の境界、紫色の数字とその間の点線は観測結果の乱流の上限値となる電離度と磁場の領域の境界を表して4Discussion

4.1

円盤の磁場の構造のモデルいる

これまでのシミュレーション結果と観測結果をまとめると原始惑星系円盤では、外側領域では弱い乱54

2022 年度 第 52 回 天文・天体物理若手夏の学校

ダストの影響4.2

ここでは今回のシミュレーションで議論しなかったダストの影響について議論する。ガスダスト比が10−4 として 1.0µm のダストを加えて乱流の大きさを計算しても加えなかった場合と結果がほとんど変わらなかった。また、サブミクロンのダストが多く存在す
る場合は、大きく電離度を下げるため乱流を抑制する

図 5: 時間ごとの磁場の強さの分布可能性があるが、この大きさのダストは coagulation
によってすぐに枯渇するため、問題にならないと判断した。

5

Conclusion
原始惑星系円盤は弱い大規模な磁場構造を持つことにより電離度が低くなっている。これによって外側領域で弱い乱流が発生している。また、外側領域
では乱流が発生していないとされている領域であっても局所的には強い乱流が発生している。

図 6: 25-35 周期での物理量の摂動の平均の分布
流が発生しつつも、内側では降着が進まなくてはな

Reference

らない。したがってこれらを説明する磁場構造のモデルとして図 7 を提唱する。このモデルでは円盤のBai, X.-N. 2013, ApJ, 772, 96

内側領域では円盤風により FUV と X 線による外側

Balbus, S. A., & Hawley, J. F. 1998, RvMP, 70, 1

領域の電離を阻止し、外側領域では大規模な磁場構造による magnetic mirroring 効果などによって宇宙

Flaherty, K. M., Hughes, A. M., Rosenfeld, K. A.,

線による電離を阻止する。このような大規模で弱い
磁場構造とそれに伴う弱い電離効果によって外側領域で弱い乱流が発生する効果が説明できる。

et al. 2015, ApJ, 813, 99
Flaherty, K. M., Hughes, A. M., Rose, S. C., et al.
2017, ApJ, 843, 150
Simon, J. B., Bai, X.-N., Armitage, P. J., Stone, J.
M., & Beckwith, K. 2013a, ApJ, 775, 73
Simon et al.2018, ApJ,865,10
Simon, J. B., Hughes, A. M., Flaherty, K. M., Bai,
X.-N., & Armitage, P. J. 2015a, ApJ, 808, 180
Suzuki, T. K., & Inutsuka, S.-I. 2009, ApJL, 691,
L49

図 7: 原始惑星系円盤の磁場構造

Stone, J. M., Gardiner, T. A., Teuben, P., Hawley,
J. F., & Simon, J. B. 2008, ApJS, 178, 13755

——–index へ戻る星惑 b02

SPH 法の不連続面における密度計算法の開発
名古屋大学大学院 理学研究科
瀬野 泉美

56

2022 年度 第 52 回 天文・天体物理若手夏の学校

SPH 法の不連続面における密度計算法の開発
瀬野 泉美 (名古屋大学大学院 理学研究科)

Abstract
星形成や惑星形成，原始惑星系円盤といった宇宙の多くの現象は流体力学によって記述される。しかし，流体力学の基礎方程式は解析的に解けないため，複雑な流体現象の解析の際は，空間的・時間的に離散化して近似解を求める数値流体計算法がよく用いられている。SPH 法とは，ある時刻で有限の大きさに離散化した流体粒子を流体に乗った座標系で追跡する数値流体計算の手法の一種で，宇宙物理学の分野で非常に幅広く使われている。SPH 法における最大の特徴は，密度を離散化した粒子の質量と Kernel 関数の足し合わせで表現していることである。本研究では Kernel 関数として Gaussian を用いるが，一般に Kernel 関数は離散化した粒子の位置と広がりの情報を持っている。この粒子の広がりはスムージング長と呼ばれ，その大きさ
は粒子間隔程度が良いとされているため粒子の変数となるべきである。しかし，このスムージング長の空間分布が正確に定義されていないため，スムージング長で定義される密度の正確な分布も分かっていない，という問題が存在する。本研究では，スムージング長の新しい定義を導入し，可変スムージング長の考えから拡張される密度定義式のうち単純な２つの方法 Gather 法と Scatter 法について詳細に再検討した。具体的には，手法による差が顕著に現れる密度が急激に変化する不連続面付近での計算を行い，理論的な考察を行う。さらに，今後の発展させるべき方向性についても議論する。

1の２方法があるが，本研究では，粒子法の一種であIntroduction

る Smoothed Particle Hydrodynamics（SPH）法を

星形成や惑星形成，円盤物理といった宇宙の物理扱う。SPH 法とは，ある時刻で流体を有限個の粒子現象の多くは流体力学で記述される。流体力学を記に離散化し，その粒子を流体に乗った座標系で追跡述する方程式は以下 4 式で，dρ= −ρ∇ · v
dt
dv
1
= − ∇pdtρdup=− ∇·vdt
ρp = (γ − 1)ρu

する方法で，離散化した粒子を Kernel 関数で広がりを持たせて定義し，それらの足しあわせによって密(1)

度を場の量として表現していることが最大の特徴である（式 (5)）。なお，本研究では Kernel 関数として(2)

Gaussian（式 (6)）を用いる (d : dimension)。
(3)
ρ(x) =(4)∑j(1) : 連続の式，(2) : 運動方程式，(3) : エネルギー

W (x, h) =方程式，(4) : （理想気体の）状態方程式である。簡m W (x − xj , h)
(

1
√

h π

)d

(

x2
exp − 2
h

(5)
)
(6)

単のため 1 次元の場合を考えると，4 変数 (ρ, v, p, u)

Kernel 関数のうち，h はスムージング長と呼ばれに対して４つの方程式があるため，この方程式系は離散化した粒子の広がりの情報を持っている。広が閉じている。従って解析的に解けるように見えるが， りとしてのスムージング長は，元は連続体であるも現実では厳密解は解かれていない。そこで，時間的・ のを離散化していることから少なくとも近傍の粒子
空間的に離散化することで近似解を求める数値流体とは相互作用すべきである，という考えで経験的に計算がよく用いられる。

粒子間隔程度が良いことが分かっている。すなわち

空間的な離散化方法としては，ある時刻で流体をスムージング長は各粒子の持つ変数であるべきであ

粒子に離散化する粒子法と空間を離散化する格子法

る（以下，可変スムージング長と呼ぶ）。しかし，こ57

2022 年度 第 52 回 天文・天体物理若手夏の学校
ので，i 粒子の密度 ρi を計算する際に h = hi とする
1粒子の分布Gather 法 (8) と h = hj とする Scatter 法 (9) の２方足し合わせ法に拡張される。
∑
ρi =m W (xi − xj , hi )

(8)

j

密度分布

ρi =

x
<latexit sha1_base64="brq0tJOPw4WCXukBViZ0idCCu2o=">AAAChHichVHLSsNAFD1GrbU+WnUjuBFLxYWUiW9cSMGNS9vaB6hIEsc6NE1CMi3W4g/oVnHhSsGF+AF+gBt/wIWfIC4V3LjwNg2ISvWGyZw5c8+9c7i6YwpPMvbUoXR2dYd6wr2Rvv6BwWhsaDjv2VXX4DnDNm23qGseN4XFc1JIkxcdl2sV3eQFvbzavC/UuOsJ29qQdYdvV7SSJfaEoUmi0gc7sThLMj/GfwM1AHEEsW7H7rCFXdgwUEUFHBYkYRMaPPo2oYLBIW4bDeJcQsK/5zhChLRVyuKUoRFbpn+JTpsBa9G5WdPz1QZ1MWm5pBxHgj2yG/bKHtgte2YfbWs1oFON9r0afo/mW+u0663a3NmJHo9m3/9VVWiX2P9S/elJYg9LvhdB3hyfabo0Wvra4flrdjmTaEyyK/ZC/i7ZE7snh1btzbhO88wFIjQg9ec4foP8TFJdSKrpuXhqNRhVGGOYwBTNYxEprGEdOerLcYJTnCkhZVqZVeZbqUpHoBnBt1BWPgGYT4/2</latexit>∑m W (xi − xj , hj )(9)j

SPH Density

図 1: 密度計算の概略図。カラーの点線が 1 粒子の広がりを表しており，それらの足しあわせによって密3

Test Calculation度分布（黒の実線）が記述される。

Gather 法 (8) と Scatter 法 (9) の手法による計算
結果の差が顕著に現れる密度が急激に変化する不連

れまでの抽象的な定義では可変スムージング長の空続面付近での計算を行う。これは，可変スムージン

間分布は明確化されておらず，スムージング長で定義される密度の空間分布も詳細に議論されていない。
本研究では，可変スムージング長の新しい定義を

グ長の抽象的な定義から密度一定の領域ではスムー
ジング長も一定であり，密度が変化する場所ではス
ムージング長も変化することから考えられる。

導入し，可変スムージング長の考えから拡張される
密度定義式のうち単純な２つの方法について詳細に再検討した。具体的には，手法による差が顕著に現3.1

れる不連続面付近での計算を行い，理論的な考察を不連続面とは，密度が急激に変化する領域のこと行う。さらに，今後の発展させるべき方向性についを言い，「岩石と大気」といった多状態を同時に記述ても議論する。

2

Methods

2.1

可変スムージング長の定義

する際には必ず出現する。

3.2

セットアップ

本稿では，不連続面を挟んだ密度比が 2:1 という既知の問題に対して計算を行う。SPH 法では密度比

本研究では，可変スムージング長を式 (7) のよう

2:1 の場合，不連続面右の粒子間隔を左の粒子間隔のに定義する。

∑
1≡
W (xi − xj , hi )
d
hi
j

不連続面

2 倍にすることで表現される（図 2 参照）。
(7)

不連続面
x
<latexit sha1_base64="w37+XIgPoQdBw8hz0n+NBsbM/a4=">AAAChHichVG7SgNBFD1ZNWp8RW0Em2BQLCTM+sZCAoJY5mFUUAm76yQO7ovdTVCDYK2tYmGlYCF+gB9g4w9Y5BPEMoKNhXd3A6Ki3mV27py55945HNXWhesxVo9ILa1t0faOzlhXd09vX7x/YM21Ko7GC5qlW86GqrhcFyYveMLT+YbtcMVQdb6u7i359+tV7rjCMle9A5tvG0rZFCWhKR5B2f1iPMlSLIjEz0RuJkk0I2PF77GFHVjQUIEBDhMe5ToUuPRtQgaDTdg2aoQ5lIngnuMIMeJWqIpThULoHv3LdNpsoiad/Z5uwNZoik7LIWYCo+yJ3bIGe2R37Jm9/9qrBpV6/D6rFszw33pAuxr25nax72Qo//Yvy6Ddw+4n609NHkqYD7QI0mYHiK9SC/nVw4tGfiE3Whtj1+yF9F2xOnsghWb1VbvJ8twlYmSQ/N2On8naZEqeTcnZ6WR6+Ti0qgPDGME4+TGHNFaQQYHmcpziDOdSVJqQpqSZsFSKNO0dxJeQFj8AvzaQeQ==</latexit>

以降，hi ≡ h(xi ) と表記を変える。この表式では，

dx0 = 0.01
<latexit sha1_base64="Mhrz6bL81CZyuBVUcfNczlQw6pQ=">AAACjnichVHLSsNQED2Nr1ofrboR3BSL4qpMRKoIRUEQl9paFaqUJL3V0DQJSVqsRXDtXlwIioIL8QP8ADf+gIt+grhUcOPCaVoQLdUJNzP33Dkz99xRbUN3PaJ6QOrq7untC/aHBgaHhsORkdEt1yo7msholmE5O6riCkM3RcbTPUPs2I5QSqohttXiSuN8uyIcV7fMTa9qi72Ssm/qBV1TPIay+cMcRZNRipOci8TY+RZtD+RWEEPL1q3IA3aRhwUNZZQgYMLj2IACl78sZBBsxvZQY8zhSPfPBY4RYm6ZswRnKIwW+b/Pu2wLNXnfqOn6bI27GLwcZkYxRc90R2/0RPf0Qp8da9Wgco3OvWp+j8Zdq+zVZm1h58Kn4+mPf1kl9h4Ovll/avJQwIKvRWdtto80VGpNfuXo/C29mJqqTdMNvbK+a6rTIys0K+/a7YZIXSDEA5J/j6M92JqNy4m4vDEXW149aY4qiAlMYobnMY9lrGEdGf9Fz3CJKykiJaSktNRMlQKt8Y7hh0lrX+q/kww=</latexit>

hd を離散化した粒子の体積（h−d を粒子の数密度）として定義している。スムージング長をスムージン

x = 0.0<latexit sha1_base64="vLNIG7dTN0hSHzcdTFxfT4/cuPg=">AAACinichVHLSsNAFD2N71q16kZwUywVV2Ui4hOhoIhLX31ALZLEaQ2mSUjSYi2Ca3euBLtScCF+gB/gxh9w0U8QlxXcuPDmAaKlekNm7py55945HNnUVNthrBkSurp7evv6B8KDkaHhkejoWMY2KpbC04qhGVZOlmyuqTpPO6qj8ZxpcaksazwrH6+599kqt2zV0PecmskLZamkq0VVkRyCsiex1RhLsoNonFYvYu2JGCRxBLFlRB+xj0MYUFBBGRw6HMo1SLDpy0MEg0lYAXXCLMpU757jDGHiVqiKU4VE6DGtJTrlA1Sns9vT9tgKTdHot4gZQ4K9sHvWYs/sgb2yz4696pCpR+dZdW+G+9Ya7bLfm5sHIxcTux//ssq0Ozj6Zv2pyUERi54WlbSZHuKqVHx+9fSqtbu8k6hPs1v2RvpuWJM9kUK9+q7cbfOdBsJkkPjbjvYkM5sU55Pi9lw8tXHuW9WPSUxhhvxYQAqb2ELae90lrtEQIsKssCSs+KVCKLB3HD9CWP8C3TORwA==</latexit>

2dx0<latexit sha1_base64="f5INBtOsR1J7NsGkLZat+IeuwTI=">AAACiHichVE9S8NQFD2NX7V+VV0EF7EoTuVGxK9JFMSxVVsFLSVJX+vDNAlJWtQiOLs4qjgpOIg/wB/g4h9w8CeIYwUXB2/Tgmip3vDy7jvvnnvf4eiOKT2f6CWktLV3dHaFuyM9vX39A9HBobRnl1xDpAzbtN1tXfOEKS2R8qVvim3HFVpRN8WWvr9Su98qC9eTtrXpHzoiU9QKlsxLQ/MZSk3nDrKUjcYoTkGMNSdqI4mhEQk7+oBd5GDDQAlFCFjwOTehweNvByoIDmMZVBhzOZPBvcAxIswtcZXgCo3Rff4X+LTTQC0+13p6AdvgKSYvl5ljmKBnuqMqPdE9vdJny14V6Nyj9axKMKP21kPe9Xpv4WQHTkc2Pv5lFXn3sffN+lOTjzzmAy2StTkBUlNp1Pnlo/PqxuL6RGWSbuiN9V3TCz2yQqv8btwmxfoVImyQ+tuO5iQ9HVdn42pyJra0elK3KoxRjGOK/ZjDEtaQQIrnSpzhApdKRCFlTlmolyqhhr3D+BHK8hfHUpHG</latexit>

図 2: セットアップの概要図

グ長で定義しているため，反復計算が必要であるが，この反復計算は各粒子がお互いの体積を調整するためにも必要である。

計算手順としては，スムージング長計算 (7) を行った後で，密度計算 ((8),(9)) を行う。なお，初期のスムージング長は不連続面左の粒子には左の粒子間隔
set up2.2密度定義式の拡張

hi = dx0 を，右の粒子には右の粒子間隔 hi = 2dx0を手で与え，収束するまで反復計算を行った。

元は固定スムージング長で表されていた密度定義
式 (5) は，可変スムージング長の考えから単純なも本計算では，運動方程式とエネルギー方程式を解かないことに注意してもらいたい。

582022 年度 第 52 回 天文・天体物理若手夏の学校

4

密度が大きくなるとスームージング長が小さくな

Results

るような結果が得られた。また，密度の振る舞いが正§3.2 のセットアップより行った計算結果を図 3 及

しく単調減少しているか，という単調性の視点に立つび図 4 に示す。

と Gather 法のみが単調性を満たし，Scatter 法では不連続面付近で値が振動してしまうことが分かった。

Smoothing Length

0.02
0.018

5

0.016
0.014

SPH 法のうち，離散化する粒子の広がりであるス
ムージング長に着目し，可変スムージング長の空間

0.012
0.01
-0.15

Conclusion

-0.1

-0.05

0

0.05

0.1

分布が明確化されていないという問題に対して，ス

0.15

ムージング長の定義式を与えることによって密度がx

図 3: スムージング長計算結果。縦軸はスムージング
長，横軸は位置を表している。

大きいところでスムージング長が小さくなるようなスキームを開発した。
また，可変スムージング長から拡張される密度定義の Gather 法と Scatter 法について手法による差が顕著な不連続面付近でテスト計算を行い，Gather 法

2のみが密度の単調性を満たすことが分かった。

Density

1.8
1.6

6

1.4

Discussion

1.2

単調性を満たすからといって万人が Gather 法を好

1
-0.15

-0.1

-0.05

0

0.05

0.1

0.15

x

んで用いている，というわけではない。Gather 法を用いると，運動方程式やエネルギー方程式の定式化(a) Gather 法

が複雑になる，という問題が存在する。この問題か
ら密度の単調性には目を瞑って Scatter 法を用いる人も少なくない。

Density2
1.8

SPH 法の密度以外の i 番目の粒子の物理量 ⟨f ⟩(xi )1.6

は流体の場の量 f (x) と Kernel 関数の畳み込み積分

1.4で表現され（式 (10)），その勾配についても同様で

1.2
1
-0.15

-0.1

-0.05

0

0.05

0.1

0.15

x

(b) Scatter 法

図 4: 図 3 のスムージング長計算結果を用いたときある（式 (11)）。
∫⟨f ⟩(xi ) ≡ dx f (x)W (x − xi , h)∫⟨∇f ⟩(xi ) ≡ dx ∇f (x)W (x − xi , h)∫
= − dx f (x)∇W (x − xi , h)

(10)

(11)

の密度計算結果。縦軸は密度を表し，横軸は位置をScatter 法の場合は ρ(x) の計算に対して h = hj の表しており，黒点線は理想的な密度分布である。また，(a) : Gather 法，(b) : Scatter 法で密度計算しため式 (11) についてた結果である。

⟨∇f ⟩(xi ) = ∇i ⟨f ⟩(xi )59

2022 年度 第 52 回 天文・天体物理若手夏の学校という変換が可能である。ここで ∇i は xi の微分を

Inutsuka 2002 では，一様質量で定式化されていた表している。対して，Gather 法の場合は h = h(x)ため，非一様な質量を取り扱う際に同様の方程式をなので ∇ → ∇i の変換に関して ∇h を考慮しなけれ用いると不連続面で作用反作用を満たさなくなり，圧ばならず，複雑化してしまう。この畳み込み積分を力のうねりが成長してしまうことが考えられる。

運動方程式・エネルギー方程式にも施し定式化する
ため，勾配項の存在を考慮すると Gather 法では複雑になることは容易に想像できる。

7.4

本研究の到達目標§7.1 ∼ §7.3 といった研究を積み重ねることによって，最終的には「大気を持った原始惑星同士の巨大7

Future Works

7.1

密度計算法の決定衝突」といった多媒体を同時に高精度で記述できる
ような SPH 法の開発を目指したいと考えている。

§6 の Gather 法による定式化も展望の一つであるが，Scatter 法の密度単調性を回復させることができAcknowledgement

れば Gather 法は用いなくて良い。従って，この問題
は Gather 法での定式化と Scatter 法の単調性回復と本講演を行うにあたり、共同研究者である犬塚修一郎教授をはじめ理論宇宙物理学研究室の皆様にはいう 2 方法のアプローチが可能である。

多くの助言をいただき大変お世話になりました。また、安定しない世の状況の中、このような研究発表の7.2

機会を設けてくださった夏の学校事務局の皆様、そ

圧力の接触不連続面問題

して星惑分科会の座長団の皆様に感謝申し上げます。

標準的な SPH 法では，圧力が一定の接触不連続面において力がゼロにならない，という接触不連続面問題が存在する。この問題に対して Inutsuka 2002([1]) では，保存則を考慮しながら運動方程式・エネル

Reference

ギー方程式を再定式化することによって，この問題[2] Inutsuka S. & Imaeda Y., 2001, CFDJを解決した（Godunov SPH 法）。さらに，Kitajima

2022 ([6]) では，Inutsuka 2002 において計算の高速化を目指し，数値積分ではなく補間により積分を評[3] Saitoh T. R. & Makino J., 2013, ApJ, 768, 44.
[4] van Leer, B. 1977, JCP, 23, 276.
[5] van Leer, B. 1979, JCP, 32, 101.

価した結果生じた接触不連続面における圧力の揺ら[6] Kitajima, 2022, 名古屋大学 修士論文

ぎを補正によって解決している。今後は，この方法[7] 富阪幸治 他 編, シリーズ現代の天文学 第 14 巻 シ
ミュレーション天文学, 2015, 日本評論社

を軸として，さらなる高精度化を目指す。

7.3[1] Inutsuka S., 2002, JCP, 179, 238.

[8] 藤井孝藏, 流体力学の数値計算法, 2021, 東京大学出版会

非一様質量を用いた SPH 法の開発

密度が小さくなると粒子間隔が大きくなるため，低解像度になることから異なる 2 つの媒体の相互作用を正確に記述できない，という問題に対して，現在系全体で質量を一様にして離散化しているところを媒体毎に異なる質量で離散化し，初期に等解像度にするような質量の非一様な粒子を扱う SPH 法の開発
を行う。

60

——–index へ戻る星惑 b03高汎用性 N 体計算コード REBOUND とその使用例
京都大学大学院 理学研究科
谷安 要

61

2022 年度 第 52 回 天文・天体物理若手夏の学校

高汎用性 N 体計算コード REBOUND とその使用例
谷安 要 (京都大学大学院 理学研究科)

Abstract
数値計算プログラムにおいては, 非公開のコードである場合論文の再現に時間を割く必要があること, 個人開発では信頼性の担保に問題が生じることの二点が大きな問題として挙げられる. これを解決しうるオープン
ソースコードの一つが REBOUND である. 本発表では REBOUND の利点や詳細について具体例を交えつ
つ紹介を行う.

1

初めにメータを確認する関数も搭載されている. その結果
視覚的にもオプションの確認が容易となり, 複雑な惑星科学において利用される計算手法として,SPH計算をする場合でもスパゲティコードになりにくい.

法や N 体計算が挙げられる. これら手法は歴史的にまた, 太陽系のデータは NASA の Horizons System

様々な発展を見せ, 計算スキームも多様で, その分大から質量と座標, 速度のデータを読み込んで研究に量にプログラムが作成されている.

用いることができる.

但しそこには大きな問題点がある. 論文で使用され二つ目の利点は多機能性である. 関数の引数の種類ているコードが非公開である場合, 論文の再現をするが豊富で, 様々な場面において REBOUND を用いために計算スキームを選択しコードを開発するといることができる. 例として粒子の衝突を判定する探うプロセスが必要になることである. また個人の開発

索アルゴリズムでは, 直接探索の他にツリー法探索

では信頼性の担保が難しくなり, 様々な手法に対応しなど 5 種類のアルゴリズムを選択することができる.

たコード開発も難しくなってしまう.

図 1 の場合直接全ての粒子との距離を計算するため,

従って, オープンソースでかつ信頼性と拡張性のあ計算量としては N 個の粒子系において O(N 2 ) 必要

るコードであるということは計算プログラムにおいであるが, 図 2 のように粒子をそれぞれ格子に分け,ては重要なファクターである. そしてそれを叶える

遠くの粒子との距離などをまとめて計算するツリー

コードの一つが REBOUND(Rain and Liu 2012) で法では O(N 2 ) に計算量を削減することができる.

ある. 本発表では REBOUND を, 自身で行った計算を具体例として用いながら紹介する.

2

REBOUND について
REBOUND は https://rebound.readthedocs.io/en/latest/において公開されている,C 及び Python で利用可能な N 体計算プログラムであり, 約 600 の論文で引用され使用されている.

図 1: 直接探索の例

REBOUND は主に 2 つの利点を持っている. 一つ目は簡易性である.REBOUND はライブラリ形式でプログラムに導入することができる. そのため, 計算
において行うオプションの変更は関数の引数で行うことができる. また, ある時間における各粒子のパラ62

図 2: ツリー法探索の例

2022 年度 第 52 回 天文・天体物理若手夏の学校

3

計算例
ここまでの REBOUND の説明を踏まえ, 実際に計算を行う.
惑星形成におけるトピックの一つとして, 地球外生命が存在するような惑星の理論的探査が挙げられる.
生命存在可能性を示唆する指標の一つとして, ハビタ
ブルゾーン (HZ) が挙げられる.HZ は普通恒星系において液体の水が存在できる領域を示すことが多い.
図 4: 計算結果

図 3 は (Kopparapu et al. 2013) において計算された, 太陽質量に対する HZ である.

つの系については HZ 内に粒子安定領域が存在するため, 惑星ができる可能性があると言うことができる.
これは過去に別の系で同様の計算を行った先行研究(Barbosa et al. 2020) と比較して同様の傾向を持っており, 計算結果としての信頼性も一定あると考えられる.
このように,REBOUND を用いて学術的な計算を行うことができ, 用いたコード自体も数十行で記述ができた. これは, 前項における簡易性によるものだと言うことができる.

図 3: 中心星質量に対するハビタブルゾーン

4

結論
REBOUND は簡易性と多機能性を兼ね備えたコー

この領域内に惑星が存在する場合, 生命存在可能性がドであり, 学術的な計算を一定の信頼性の下で行うこ生じるということである.
今回計算するのは HZ 周辺に巨大惑星が存在するとができるコードである. 実際 HZ 周りの粒子安定性

場合の重力安定性である. これは, 惑星移動によってについて計算した結果, 先行研究と同様の傾向を示す

HZ 付近にホットジュピターのような巨大惑星が移動 結果を得ることができた.
本発表を通じて少しでも REBOUND に対し興味してきた場合, 惑星の素である原始惑星, 微惑星が巨を持っていただき, 研究の際のコードの一つとして一大惑星からの重力の影響を十分に受けるためである.

但し今回は惑星移動を考慮せず, 先に移動が終わっ 考していただけると幸いである.
たとしてその状態での重力安定性のみを検証する.
計算手法としては REBOUND を用い,106 年積分を行う. 中心星及び惑星はテスト粒子に重力を与えるが, テスト粒子は他の粒子に重力相互作用を与えないものとして計算する. 巨大惑星を HZ 付近に持つ恒星
系 3 つについて REBOUND を用いて計算した結果

Reference
H.Rain ,and S.F.Liu,REBOUND:an open-source multipurpose N-body code for collisional dynamics,
A&A,537,128 (2020)
R.K.Kopparapu et al.,Habitable Zones Around Mainsequence Stars: New Estimates, ApJ,765,131(2013)が図 4 の通りである.
黒点で表したテスト粒子は左右の二系では離心率が大きくなりながらも残存しているが, 中心の系ではHZ 外に粒子が安定して存在した. 結果として左右二G.O.Barbosa et al.,Earth-size planet formation in
the habitable zone of circumbinary stars,MNRAS
494,1045-1057(2020)

63

——–index へ戻る星惑 b04太陽系外惑星 TrES-1b の多波長トランジット観測
兵庫県立大学理学研究科光学赤外線天文学研究室
植野 雅々

64

2022 年度 第 52 回 天文・天体物理若手夏の学校

太陽系外惑星 TrES-1b の多波長トランジット観測
植野 雅々(兵庫県立大学理学研究科光学赤外線天文学研究室)

Abstract
1995 年に Mayor と Queloz が 51 Peg の周りを周回する惑星 51 Peg b を発見して以来、様々な方法で系外惑星の研究が続けられており、現在 5000 個以上の系外惑星が発見されている。主星近傍を公転するガス惑星のホットジュピターや地球の数倍の質量を持つ岩石惑星の スーパーアースなど多様な惑星がある。惑星を検出する方法の一つにトランジット法がある。系外惑星が大気を持っている場合、トランジット中に主星の光の一部が惑星の上層大気を通過し、原子・分子によって吸収されることがあり、波長によってわずかにト
ランジットの深さが異なる。本研究では、西はりま天文台なゆた望遠鏡に搭載された近赤外撮像装置 NIC を用いて太陽系外惑星 TrES-1b のトランジットを観測した。近赤外線の波長領域 (0.8〜5.0 µm) には分子の吸収が多く、惑星の大気の成分を調査するには適した波長である。取得した画像は画像処理ソフト IRAF で処理・解析を行い、EXOFAST を用いて光度曲線を求めた。

1

イントロダクション搭載された近赤外撮像装置 NIC を用いて太陽系外惑

星 TrES-1b のトランジットを観測した。NIC は J バ
1995 年に Mayor と Queloz が初めて太陽系外惑星 ンド (中心波長 1.25 µm)、H バンド (1.63 µm)、Ks
51 Peg b を発見して以来、様々な方法で太陽系外惑 バンド (2.15 µm) を同時に観測できる。観測された星の研究が続けられている。主星近傍を公転するガ

トランジットの深さから波長ごとの主星と惑星の半ス惑星のホットジュピターや地球の数倍の質量を持径比を求めた。

つ岩石惑星の スーパーアースなど多様な惑星がある。
惑星を検出する方法の一つにトランジット法がある。
トランジット法とは惑星が主星の前を通過する際の

2

観測・解析

2.1

観測

減光を捉え、惑星を検出する方法である。 主星の明
るさの変化量から、以下の (１) 式を使って主星と惑星の半径比を求めることができる。

∆F
=
F

(

Rp
Rs

観測は西はりま天文台の２ｍなゆた望遠鏡に搭載)2された近赤外撮像装置 NIC で観測を行った。観測の(1)

パラメーターは以下の通りである。

∆F は変化した明るさ、F は主星の明るさ、Rp
は惑星の半径、Rs は主星の半径である。系外惑星が大気を持っている場合、トランジット中に主星の光表 1: 観測のパラメーター
観測日
トランジット予想時間の一部が惑星の上層大気を通過し、原子・分子に吸観測天体収されることがある。惑星の大気による吸収は、波長によってわずかにトランジットの深さを変化させる。吸収が強い波長ほどトランジットの深さが深くなり、主星と惑星の半径比は大きく導出される。近赤外線の波長領域 (0.8〜5.0 µm) には分子の吸収が多く、惑星の大気の成分を調査するには適した波長である。本研究では、西はりま天文台なゆた望遠鏡に

65

2021 年 5 月 30 日

21:50〜26:00

22:34〜25:04
TrES-1（TrES-1b）

バンド

J、H、Ks

露出時間

画像 1 枚につき 30 秒

画像枚数

J,H,Ks バンドそれぞれ 210 枚

ディザリング

半径 10 秒角で 10 カ所

2022 年度 第 52 回 天文・天体物理若手夏の学校

表 2: 観測天体のパラメーター

TrES-1 (主星)
半径 [太陽半径]+0.100.85 −0.05

質量 [太陽質量]0.88 ± 0.07

有効温度 [K]5250

V 等級 [mag]11.79J 等級 [mag]10.294

H 等級 [mag]9.887

K 等級 [mag]スペクトル型

9.819
K0

TrES-1b (惑星)
半径 [木星半径]+0.181.08 −0.04

質量 [木星質量]0.75 ± 0.07

公転周期 [日]3.03
Alonso et al. (2004)

2.2

解析

解析には画像処理ソフト IRAF を用いた。まず、ダーク画像を作り、ダーク引きを行い、トワイライ

図 1: TrES-1b の光度曲線

トフラットを使用してフラット補正を行った。次に、
バットピクセル補正を行い、スカイパターン画像を外の相対フラックスを用いて回帰曲線を求め、回帰取得してスカイ補正を行った。星のない領域から縦曲線上の値を１として規格化を行い、EXOFAST を

縞のパターン画像を作り、縦縞パターンを除去し、最用いてトランジットの光度曲線をフィットした。(図後に宇宙線を除去した。測光は開口測光を用いた。測1)

光半径は J、H バンドは 4.64 秒角、Ks バンドは 4.16

J、H、Ks バンドすべてでトランジットを検出でき

秒角、スカイ内径は 6.88 秒角、スカイ幅は 0.80 秒

た。EXOFAST で得られたトランジットの深さは表角に設定し測光を行った。

3 の通りである。

3

結果

4

解析により得られた測光結果から光度曲線を作成した。対象星のフラックスを標準星のフラックスで考察得られたトランジットの深さから上の (1) 式を用

いて主星と惑星の半径比を求めた。(表 4)割ったものを相対フラックスとした。トランジット
表 4: 主星と惑星の半径比
バンド

表 3: 観測結果

半径比

トランジットの深さJ

0.0233 ± 0.0009

J

0.0233 ± 0.0009

H

0.0231 ± 0.0006

H

0.0231 ± 0.0006

Ks

0.0200 ± 0.0008

Ks

0.0200 ± 0.0008

z (Winn et al. 2007)

0.1369 ± 0.0008

バンド

66

2022 年度 第 52 回 天文・天体物理若手夏の学校
今回の観測では、半径比の値が J バンドと H バン
ドは同程度、Ks バンドで最も小さくなった。1 標準偏差で考えると、J バンドの半径比の最小値は 0.1496、

H バンドの半径比の最小値は 0.1502 であり、Ks バ
ンドの半径比の最大値 0.1445 と比較しても値が大きいので、Ks バンドで最も半径比が小さいと言える。

5

結論
TrES-1b を J バンド、H バンド、Ks バンドの３つの波長で観測を行い、すべてのバンドでトランジッ
トが検出できた。半径比は、J バンドと H バンドが同程度、Ks バンドで最も小さくなった。今後、得られた結果から TrES-1b の大気組成を検討していく。

Reference
Alonso, R., Brown, T. M., Torres, G. et al. 2004, ApJL,
613, L153-L156
Winn, J. N., Holman, M. J. and Roussanova, A. 2007,
ApJL, 657, 1098-1106
多葉田 晴日 2017, 修士論文
平野 佑弥 2021, 修士論文

67

——–index へ戻る星惑 b05太陽系外惑星の天気の調査
兵庫県立大学大学院 光学赤外線天文学研究室
平野 佑弥

68

2022 年度 第 52 回 天文・天体物理若手夏の学校

太陽系外惑星の天気の調査
平野 佑弥 (兵庫県立大学大学院 光学赤外線天文学研究室)

Abstract
太陽系外惑星が晴れた大気を持つ場合、トランジット観測で求められる惑星の半径は波長によって異なる。
近赤外領域には大気分子の振動遷移による大きな吸収が見られ、多波長で太陽系外惑星の半径の波長依存性を調べることで大気組成の推定ができる。一方で惑星大気中に雲や微粒子の存在があれば、分子の吸収の特徴を隠してしまうことがある。しかし惑星の半径が観測された約 4000 天体のうち、大気組成の推定が行われた惑星は約 100 天体しかない。そこで 2020 年から 2021 年の期間に西はりま天文台なゆた望遠鏡に搭載された近赤外撮像装置 NIC を用いて、複数の太陽系外惑星の近赤外トランジット測光観測を行った。2021年に行った Qatar1 b の近赤外トランジット測光観測では、惑星大気モデルと比較することで CO2 を含む晴れた大気を持つ可能性を示した。しかし Covino et al. (2013) の可視光領域 R バンドの観測では作成した大気モデルより約 30%大きな惑星半径を示した。この不一致は観測時期により大気中の微粒子や雲の状態が変わったためとも考えられる。そこで西はりま天文台のなゆた望遠鏡と 60 cm 望遠鏡を用いて長期的かつ高頻度に多波長で観測し惑星半径の変動を捉え太陽系外惑星の天気を観測する。解析の結果、新たに 2 回のトランジットが検出でき、惑星半径の変動を示す結果が得られた。

1

イントロダクション

作成し、観測値との比較を行った。大気モデルと観測した減光率から CH4 の大気を持つと推定した。主1930 年代から太陽系外惑星の探査が本格的に始ま 星近傍を周回する高温のガス惑星であるホットジュ
り、1995 年に 51 Peg という恒星の周りに巨大なガ ピターは H を主成分とする大気を持ち、H O, CH ,2

ス惑星が発見された。2009 年にはケプラー宇宙望遠

2

4

CO2 などの分子が検出されている。これらの分子は鏡が運用を開始した。50 万個以上の恒星を観測し、 近赤外領域で振動遷移により分子固有の大きな吸収

2500 個以上の惑星を発見した。2022 年現在、5000を示す。また惑星大気中に雲や微粒子が存在する場

個以上の太陽系外惑星が発見されており、その多く

合、分子の吸収の特徴を隠してしまうことがあるとはトランジット法で発見されている。トランジット

分かってきた。

法とは、惑星が恒星の前を通過する際のわずかな減

2021 年 9 月 15 日に西はりま天文台なゆた望遠鏡

光を捉え、周期的な明るさの変化により太陽系外惑に搭載された近赤外撮像装置 NIC を用いて、Qatar1星を発見する手法である。Charboneau et al. (1999)

b の近赤外トランジット測光観測を行った。近赤外のと Henry et al. (1999) によりトランジット法で初め

3 バンド (J, H, Ks ) で Qatar1 b のトランジットを検て太陽系外惑星が発見された。

出することができた。その後、惑星大気モデルを作成

トランジット法では惑星に大気が存在する場合、分し、大気モデルと求めた主星と惑星の半径比を比較し子や原子の吸収により波長によって減光率が異なっ た。大気モデルの作成には Planetary Spectrum Genて観測される。そのため多波長でトランジットの減 erator (https://psg.gsfc.nasa.gov/atmosphere.php)
光率の波長依存性を調べることで大気組成の推定がを用いた。近赤外領域のみを考えると Qatar1 b はできる。これまでトランジットが観測された約 4000

H2 99%, CO2 1%の大気モデルが最もよくあった。一

天体のうち約 100 天体で大気組成の推定が行われて方で、Covino et al. (2013) が 2011〜2012 年に観測いる。Swain et al. (2008) はハッブル宇宙望遠鏡を用した可視光領域 R バンドでは、作成した大気モデル
いて HD189733 b の近赤外トランジット観測を行っ より約 30%大きな惑星半径を示した。
た。H2 O と CH4 を含む H2 が支配的な大気モデルを

69

2022 年度 第 52 回 天文・天体物理若手夏の学校
しやすい。また公転周期は約 1.4 日と短く、高頻度にトランジットを観測できる。
表 1: 主星と惑星のパラメータ

Qatar1 (主星)

図 1:Qatar1 b 混 合 大 気 ス ペ ク ト ル (H2 99%,半径 (太陽半径)0.82 ± 0.03

質量 (太陽質量)0.85 ± 0.03

有効温度 (K)log g (cgs)4860 ± 130
4.54 ± 0.02

V 等級 (mag)12.843 ± 0.137 (1)J 等級 (mag)11.00 ± 0.02 (2)H 等級 (mag)10.53 ± 0.02 (2)Ks 等級 (mag)

10.41 ± 0.02 (2)Qatar1 b (太陽系外惑星)

CO2 1%) と観測で求めた主星と惑星の半径比
可視光領域で惑星半径が大きく観測されることは、McCullough et al. (2014) では黒点によるもの、

Kawashima & Ikoma (2019) では大気の渦拡散の効率によるものと考えられると説明しており、解釈は半径 (木星半径)1.16 ± 0.05

質量 (木星質量)1.09 ± 0.08

平衡温度 (K)

1400 ± 40

距離 (AU)0.0234 ± 0.0004

log g (cgs)3.27 ± 0.04

Alsubai et al. (2011) より(1) APASS9 http://www.aavso.org/apass(2) 2MASS https://irsa.ipac.caltech.edu/Missions/2mass.html

様々である。本研究では 2011〜2012 年と 2021 年で大気中の微粒子や雲の状態が変わったことにより引き起こされたと考えた。同じ太陽系外惑星を長期間

西はりま天文台なゆた望遠鏡 (口径 200 cm) に搭

かつ高頻度に多波長で観測することで、惑星の天気載された近赤外撮像装置 NIC(J, H, Ks バンド同時を惑星半径の変動として観測できる可能性がある。そ撮像) と 60 cm 望遠鏡に搭載された冷却 CCD カメ

こで西はりま天文台なゆた望遠鏡 (口径 200 cm) とラ STL-1001(V バンド) を用いて Qatar1 b の近赤外

60 cm 望遠鏡を用いて、同じ太陽系外惑星を長期間

トランジット測光観測を行った。なゆた望遠鏡と 60

かつ高頻度に多波長で観測し、雲や微粒子による惑

cm 望遠鏡は同天文台にあるが、距離が離れているた星半径の変動を捉え太陽系外惑星の天気を観測する。 め、なゆた望遠鏡はリモートで操作し、60 cm 望遠
本発表では 2022 年 5 月 23 日と 2022 年 7 月 6 日に鏡は直接操作することで J, H, Ks , V バンドを同時行った Qatar1 b のトランジットの観測・解析結果をに観測した。

紹介し、本研究と過去に Qatar1 b のトランジットが観測された先行研究の比較を行い、太陽系外惑星の

3

天気の変動の可能性について議論する。

解析
解析には画像解析ソフト IRAF を用いた。最初

2にダーク引き、フラット補正を行った。フラット補観測

正にはトワイライトフラットを使用した。次に検出

Qatar1 は約 60 pc 離れた K 型主系列星である。 器由来のバッドピクセル、宇宙線をそれぞれ補正しQatar1 b は 2011 年に New Mexico Skies Observa-

た。またスカイパターンを取得してスカイ引きを

tory の Qatar Exoplanet Survey によって発見され行った。NIC で取得した画像は検出器の読み出し列

た (Alsubai et al. 2011)。Qatar1 b はホットジュピ

を切り替える際に電圧が変動するためスカイ引き後ターであり減光率も 1%以上あるため地上からも観測には縦縞が目立つ。画像ごとに縦縞パターンは異な702022 年度 第 52 回 天文・天体物理若手夏の学校

図 2: 相対フラックスと光度曲線。H バンドのフラックスは 0.1, Ks バンドのフラックスは-0.2 下げている。

V バンドはフラックス+0.1 上げている。

るので、画像の上半分と下半分のそれぞれ星のない求められた。7 月 6 日のトランジットでは 2021 年 9

領域から縦縞パターンを作成し縦縞を除去した。

月 15 日より J バンドは約 20%, H バンドは約 5%大測光は開口測光を用いた。対象星のフラックス

きく、Ks バンドは約 5%小さな半径比が求められた。

を標準星のフラックスで割ったものを相対フラッ また近赤外領域の半径比は波長によらず、ほぼ一定
クスとした。アパーチャー半径を FWHM の 0.5〜である。V バンドは Covino et al. (2013) と同様に、2.0 倍に設定し、トランジット外の相対フラック

近赤外領域より大きな半径比を示した。

スを求めた。標準偏差が最も小さくなる倍率を測光半径に設定した。次にトランジット外の相対フ

表 2: Qatar1 b 主星と惑星の半径比

ラックスを用いて回帰直線を求め、回帰直線の値

2022 年 5 月 23 日

を 1 として規格化を行った。NIC はディザリング

バンド

観測を行ったため、トランジット外の相対フラッV-

クスをディザリング位置ごとに平均値を求め、対J0.131 ± 0.015

応する同じ位置の画像の相対フラックスを平均値

H

0.156 ± 0.006で割った。最後に 2σ 以上の外れ値は通り除き、

Ks

0.153 ± 0.012

EXOFAST(https://astroutils.astronomy.osu.edu/
exofast/exofast.shtml) を 用 い て 光 度 曲 線 を 作 成

2022 年 7 月 6 日
バンド

した。

4

主星と惑星の半径比結果2022 年 5 月 23 日と 7 月 6 日の観測でトランジッ

トを検出することができた。5 月 23 日のトランジッ
トでは、2021 年 9 月 15 日より J バンド, H バンド
は約 15%大きく、Ks バンドは約 5%大きな半径比が71

主星と惑星の半径比

V0.158 ± 0.001

J

0.140 ± 0.008

H

0.139 ± 0.005

Ks

0.137 ± 0.011

2022 年度 第 52 回 天文・天体物理若手夏の学校

図 3: 本研究で得られた主星と惑星の半径比と先行研究の比較

5

考察

Reference

本研究と過去に Qatar1 b のトランジットが観測された先行研究の比較を行った。可視光領域の半径比Henry et al. 2000, The Astrophysical Journal, Vol. 529,
No. 1, p. L41.

は von Essen et al. (2017) の分光観測も含めて近赤外領域より大きな値を示している。また Ks バンド
で観測を行った Cruz et al. (2016) と 2021 年 9 月 15
日, 2022 年 5 月 23 日の半径比は誤差範囲内で一致している。一方で、2022 年 7 月 6 日の近赤外領域の半径比はほぼ一定であり、2021 年 9 月 15 日の半径比気モデルを作成し、惑星半径の変動を再現できるか検証する。

2021 年 9 月 15 日に加え、新たに 2 回のトランジッ
トを検出し、惑星半径の変動を示す結果が得られた。
時期による比較を行うことで、太陽系外惑星の天気の変動を検出できる可能性がある。

Swain et al. 2008, Nature, Vol. 452, No. 7185, pp.
329–331.
Covino et al. 2013, Astronomy & Astrophysics, Vol. 554,
p. A28.
McCullough et al. 2014, The Astrophysical Journal, Vol.
791, No. 1, p. 55.

の特徴と大きく異なる。今後は雲の割合を変えて大さらに今後も高頻度にトランジットを観測し、観測

Charboneau et al. 1999, The Astrophysical Journal, Vol.
529, No. 1, p. L45.

Kawashima & Ikoma 2019, The Astrophysical Journal,
Vol. 877, No. 2, p. 109.
Alsubai et al. 201,1Monthly Notices of the Royal Astronomical Society, Vol. 417, No. 1, pp. 709–716.
von Essen et al. 2017, Astronomy & Astrophysics, Vol.
603, p. A20.
Cruz et al. 2016, Astronomy & Astrophysics, Vol. 595,
p. A61.

72

——–index へ戻る星惑 c01

スパースモデリング (SpM) による超解像度画像解析
九州大学大学院 理学府 地球惑星進科学専攻
所司 歩夢

73

2022 年度 第 52 回 天文・天体物理若手夏の学校

スパースモデリング（SpM）による超解像度画像解析
所司

歩夢 (九州大学大学院 理学府 地球惑星進科学専攻)

Abstract
Atacama Large Millimeter/submillimeter Array (ALMA) の観測によって、近傍の小質量星形成領域にはダスト半径 r がおよそ 20 − 40au 以下のコンパクトな原始惑星円盤が多く存在することが明らかとなった。

DSHARP や ODISEA の観測によって、比較的大きな原始惑星系円盤にはリング・ギャップ構造が普遍的に存在することが分かってきた。しかし、多数派であるコンパクトな原始惑星系円盤にそのようなサブストラ
クチャーがあるかどうかは空間分解能が足りず、調べられていなかった。Yamaguchi et al.(2021) の新しく開発された画像解析技術であるスパースモデリング（SpM）を T Tau 星系に適用し、原始惑星系円盤の数

au スケールの構造を調べた。
今回は Yamaguchi et al.（2021）のレビューを行う。そして、Cieza et al. (2020) でリング構造の存在が明らかになっている Oph 03 について、低解像度データに対し SpM を適用した結果を紹介し、今後の展望を議論する。

1

不完全ないビジビリティデータを復元することができIntroduction

る。ALMA 望遠鏡の Band4 と Band6 の実観測データ

惑星は恒星周りの原始惑星系円盤 (PPDs) で形成される。PPDs の主成分はガスとダストであり、円盤の構造と進化はコア降着モデルやディスク不安定モデルなどを含む惑星形成過程に密接に関係している。
原始惑星が非常に大きな質量を持つ場合、磁気流体に SpM を適用すると、CLEAN 画像よりも高精度の画像を得ることができた。
（Yamaguchi et al.(2020)）

SpM を使うことで超解像度を実現すると同時に、計算スピードが向上するため、長時間の観測が必要なくなる。

効果や呪力不安定性などにより円盤内にギャップが生まれる。ギャップができる最小質量は円盤の粘性と大きさに依存するため、ギャップを空間的に分解する
ことで成長する惑星の質量を推定することができる。

ALMA 望遠鏡の DSHARP や ODISEA の観測に

今回、Yamaguchi et al.（2021）のレビューを行い、SpM を T Tau 星系に適用した結果を紹介する。そして、Cieza et al.(2020) でリング構造の存在が明らかになっている Oph 03 について、低解像度データ
に対し SpM を適用した結果を紹介し、今後の展望をよって、空間分解能 2 − 5au で半径 r = 50 − 260au

議論する。

の大きな円盤画像が得られ、環状ギャップ構造が明らかになった。小さい円盤は小質量星形成領域で多いことが分かっているが、サブストラクチャーがあるかどうかは空間分解能が足りず調べられていない。

2

Yamaguchi et al.

そのような構造を研究するためには、35 ミリ秒角（mas）の滝分解能が必要であるが、ALMA 望遠鏡のBand6 観測で得られた最も高い分解能は 20 − 50mas
であり、長時間の観測が必要である。

Methods(2021) では T Tau 三重星周辺の PPD に注目している。T Tau は低質量で前主系烈星の三重星である。この星系は北の星 T Tau N と南の連星 T Tau Sa/Sb で構成されており、143.7±1.2pc 離（Gaia Colそのような空間分解が十分でない天体に対して、 れた Taurus 星形成領域に位置している。
Yamaguchi et al.(2021) では新しく開発された画像 laboration et al. 2016,2018)T Tau N は Class II、解析技術がスパースモデリング（SpM）を提案されT Tau Sa/Sb は ClassI に分類されている。T Tauている。SpM による画像解析は空間周波数における Sa と Sb の質量はそれぞれ 2.1M⊙ と 0.4M⊙ である。

74

2022 年度 第 52 回 天文・天体物理若手夏の学校

図 1: SpM を用いた画像解析のフローチャート
則化パラメータを選択するために、統計的方法である(Schaefer et al. (2020))T Tau N は Turus 内の T

Tauri 星系で最も明るい星の一つである。いくつかの CV 法を用いる。10-fold CV 法では、最初にデータ研究で抗精神科モデルと組み合わせた光学スペクトルセットを 10 分割する。10 分割したデータのうち、9タイプを使用して計算されており、放射強度 6.82L⊙

セットを trainingと質量 2.19M⊙ を採用している。以前の ALMA 観

して扱う。 Λl と Λtsv に任意の値を設定し、training

測では T Tau 星系のダスト円盤を対象にしていた

data から SpM 画像解析を用いて、モデル画像を生が、120mas の空間分解能が不十分のため、CLEAN成する。モデル画像に 2 次元フーリエ変換を行い、画像で空間的に分解できなかった。T Tau Sa/Sb 円

model data を得て、CVE の式に残りの validation

盤は単一のガウス分布のように見えた。T Tau N はdata を代入することで計算し、数値的に評価を行う

data、残りを validation data と∼ 180mJy のトータルフラックスを持つ半径 ∼ 20au （Wi はトーレス行列）。Λl と Λtsv の値を変え、同じ作業を繰り返す。CVE が最も小さい値をとるものをの平坦なコンパクト円盤であると判明した。
この T Tau 三重星に適用する SpM 画像解析と交差

適切な値とし、その時再構築したものを SPM 画像と検証 (CV) 法について説明する（図１参照）。画像は輝する。SpM では実際に観測を行なっていないため、度分布の L1 ノルムと Total squared variation(TSV)

ビームサイズがなく、CLEAN 画像の分解能と比較の 2 つの正則化項が使われているコスト関数を最小することができない。この問題は中心点源のトータ限に抑えることによって再構築される。Λl と Λtsv のルフラックスの 5 ∼ 10% の大きさの人工点源の挿入

正則化パラメータの値によって、2 つの正則化項は再し、その大きさを空間分解能とする。人工点源を挿構築された画像のスパース性と滑らかさを調整する

入した SpM 画像をもとにサブストラクチャーの詳細ことができる。画像を再構築するにあたって適切な正に迫っていく。

752022 年度 第 52 回 天文・天体物理若手夏の学校

図 2: T Tau 星系に SpM と CLEAN を適用した結果（Reference[1] の Figure 1 より一部抜粋)

3

れぞれの total ﬂux の大きさは降着円盤モデルによっResults
T Tau 星系に SpM と CLEAN を適用した結果が図て予測される各 SED と大まかに一致している。このような T Tau 星系の解像度の信頼性と所見はより良2 である。SpM 画像の空間分解能は 38mas である。
い角度分解能と感度を備えた将来の観測によって正CLEAN 画像は robust=0.5 の Briggs
weghiting
当性を評価し、サブストラクチャーの存在を固めるで 140 × 100mas となっている。つまり、SpM 画像はことができる。
CLEAN beam size の ∼30% 以下の空間分解能を達成している。SpM 画像では T Tau N 周辺のディスク構造を空間的に分解し、環状ギャップ構造が新たに発見された。CLEAN 画像では空間的に分解できていない T Tau S 星系も２つの放射源に分解することができた。SpM 画像から得られた放射源のトータルフ
ラックスは 5σ を超える CLEAN 画像から得られた値と一致している。SpM は CLEAN 画像よりも観測されたビジビリティデータに合わせて忠実性の高い画像を再現するが、CLEAN モデルと同等の可視性領域で同様の結果を提供する。また、SpM 画像は円盤表面の輝度分布をよりよく再構築し、CLEAN 画像
では円盤の構造を正確に反映していないいくつかの点源の合計で画像を再構築すると考えられる。よって、T Tau 星系の物理的性質を導き出すために SpM
画像を利用する。
図２は T Tau N 円盤のピークを中心に極座標表示
したものと全方位角で平均化された放射強度マップ
を表している。r = 11.6 ± 0.3au の環状ギャップ (円盤半径 24 ± 4au) と r = 15.7 ± 0.1au の外側のピー
クを確認することができた。また、T Tau 星系のそ

図 3: T Tau N の円盤構造を極座標に投影したマップ（Top Panel）と全方位角で平均化された放射強度
マップ (Bottom panel) (Reference[1] の Figure 3)

762022 年度 第 52 回 天文・天体物理若手夏の学校

図 4: Oph 03 の低解像度データを SpM(左) と CLEAN(中央) で解析した結果と高分解能データ CLEAN 画
像 (右)

4

SpM 画像解析についてより深く理解すると同時に、

Futute Prospects

EHT のブラックホール観測の正当性を評価すること
Cieza et al. 2020 で 2 つの環状ギャップ画報告さができる。
れているのが Oph 03 である（図 4 の右)。同じ天体の低解像度データに SpM と CLEAN で画像解析し、それらの結果と高解像度データの CLEAN 画像と比Acknowledgement

較を行った。低解像度データを CLEAN したものは空間分解能 328mas でリング構造を確認できなかっ
たが、SpM 画象は空間分解能 107 × 96mas で 1 つの本研究を進めるにあたり、山口正行 氏（国立天文台)、町田正博 氏、徳田一起

氏、原田直人 氏に多くのご指導をいただきました。この場を借りて感謝

環状ギャップ構造を確認することができた。このこ申し上げます。また、データ解析には国立天文台のとから SpM 解析は CLEAN 解析のビームサイズの多波長データ解析システム (ADC) を使用しました。

30% 以下を達成することができており、サブストラ
クチャーを判定するには 300mas 以下である必要があることが分かった。そして、高解像度データと同

References

様の結果を得るためには空間分解能が 100mas 以下

Cieza, L。A., González-Ruilova, C., Hales, A., et alである必要があると予想される。

2022 年 5 月 12 日に天の川銀河の中心の巨大ブラッ 2020
クホール Sagitarrius A∗ の SpM 画像解析が施され Gaia Collaboration, Brown, A. G. A., Vallenari, A.,た画像が公開された。しかし、ブラックホールシャ et al. 2018, A&A, 616, A1
ドウがリング状であるか、またその正当性について Gaia Collaboration, Prusti, T., de Bruijne, J. H. J.,
賛否両輪がある。そこで ALMA 観測で EHT の観測

et al. 2016, A&A, 595, A1

シミュレーションを行う。既にリング構造の存在がSchaefer, G. H., Beck, T. L., Prato, L., & Simon,

分かっている天体データを用意し、EHT に模して基M. 2020, AJ, 160, 35
Yamaguchi, M., Akiyama, K., Tsukagoshi, T., et al.

線長やアンテナ配置で限定をかける。そのデータにSpM 解析を行うことでリング構造が確認できるか判 2020, ApJ, 895, 84
断するものである。このシミュレーションを通して、 Yamaguchi, M., Tsukagoshi, T., Muto, T., et al.
2021, APJ, 923, 121

77

——–index へ戻る星惑 c02若い星状天体の円盤、アウトフロー、磁場の初期進化とダストモデルの違いによる影響
九州大学大学院 理学府地球惑星科学専攻
赤池 祐貴

78

2022 年度 第 52 回 天文・天体物理若手夏の学校

若い星状天体の円盤、アウトフロー、磁場の初期進化とダストモデルの違いによる影響赤池祐貴 (九州大学大学院 理学府地球惑星科学専攻)

Abstract
今回の講演では、Tsukamoto et al. 2020 をレビューし、今後の課題についての議論をする。
低質量の若い星状天体 (young stellar objects , 以下 YSOs) の形成や時間発展は、３次元の非理想磁気流体
シミュレーションを用いて調べられている。
本研究では、YSOs の質量が 0.1M⊙ に達した、原始星形成後の 104 年までの時間発展による物理現象を調べている。特にダストサイズの違いによる影響に着目し、5 つのモデルを比較していく。結果、どのダストモデルでも、星周円盤は形成されることがわかった。円盤は原始星形成時には約 10AU、104 年後には数十AU にまで達し、円盤質量は中心星質量に匹敵すること、重力不安定性が成長していることがわかった。更
に、小さいダストモデルのシミュレーションにおいて、104 年後 pseudodisk の揺らぎの成長していることがわかった。この揺らぎは、円盤の磁気制動を強める働きを持ち、円盤のサイズは小さくなる。Ion-neutral

drift はダストサイズ ≧ 0.2µm、M ≧ M⊙ でエンベロープ内で発生しうることを見つけた。アウトフローは、ダストのサイズと反相関を示し、大きいダストサイズでは弱いアウトフロー、小さいダストサイズでは強い
アウトフローを出すことを確認した。

1

約 104 年までの初期円盤の進化について、特に、ダIntroduction
原始星や原始惑星円盤が誕生する分子雲コアでは、強力な磁場がかかっている。その磁場は、重力崩壊にストの大きさの違いの影響に着目して研究を行っている。

抵抗するのに十分強くなくても、重力崩壊の段階において重要な役割を果たす。
（Inoue & Inutsuka 2012）

2

Numerical

Methods

また、分子雲コアスケールでのダストの成長の可能性

non-ideal MHD eqation の式として次の４本、に関して、μ m サイズのダスト粒子の散乱光が分子雲

Dv
1
1= − (∇(P + |B|2 − ∇ · (BB))) − ∇ΦDtρ2

コア内の中赤外線放出として解釈できる（Steinacker

et al. 2010; Pagani et al. 2010)。さらに、1 μ m 以
下のサイズのダストが分子雲コア内で成長を果たす(Ormel et al. 2009; Hirashita & Li 2013) ことも報告されている。天体現象との関わりについて、小さいDB= (B · ∇v − B(∇ · v)) − ∇ × (η0 (∇ × B)
Dt
−ηA ((∇ × B) × B̂ × B̂))

サイズのダストによる円盤形成 (Masson et al. 2016;

Wurster et al. 2016; Tomida et al. 2015; Tsukamoto∇2 Φ = 4πGρet al. 2017a; Wurster et al. 2018; Wurster & Bate2019) があり、ダストが天体の減少に関わる事がわ

P = P (ρ) = c2s.iso ρ(1 + (

かる。

ρ 2/3
) )
ρcrit

Tsukamoto et al. 2020 では、ダストのモデルの影 SPMHD 法 (Iwasaki & Inutsuka 2011, 2013) を用い響を調べる研究が必要であると考え、原始星形成後 て方程式を解き、シミュレーションを行った。

792022 年度 第 52 回 天文・天体物理若手夏の学校

図 1: ダストモデル

図 4: model-a03µm x-z

図 2: model-MRN x-z

500AU

500AU

図 5: model-a03µm x-y 200AU

図 3: model-MRN x-y

3

200AU

Results

図 6:
2000AU

シミュレーションで比較したのは、図 1 の model-

MRN と model-a03µm である。各モデルの時間発展について、x-z 図として横から、x-y 図として上から見た様子を表す。また、red arrow は質量流速 whitearrow は磁場を表す。
図２と図４から小さいダスト大きいダスト関わらずローを出している事がわかる。また、図３図５から、大きいダストモデルの方が不安定性としてみられる

spiral arm が発展し続けており、小さいダストモデ
ルの方が不安定性が抑えられているのがわかる。そして、円盤の密度に強い依存がないこともわかる。

アウトフローが出ている事がわかる。さらに、図６
図 7 より、どちらも 1000AU を超える強いアウトフ

model-MRN x-z 図 7: model-a03µm x-z
2000AU

図 8 から、t = 5 × 104 年までは、どのモデルも円盤半径は増加しているが、red,black,yellow line の80

2022 年度 第 52 回 天文・天体物理若手夏の学校

図 10: 本研究と観測結果について

図 8: 円盤半径の時間発展性

4Discussion
本研究で扱った非理想磁気流体効果は円盤形成を起こす有望なメカニズムとなっている。実際に非理想

MHD シミュレーションを用いた先行研究で確認されている (Tsukamoto et al. 2015b)。しかし、(Mellon& Li 2008) による The magnetic braking catastrophe は、原始星進化の初期段階で magnetic braking
により円盤の形成が完全に抑制されると主張しており、長い間、原始星形成の課題となっており、円盤が生き残るかどうかはまだわかっていない。現在の研究で、円盤は数十 AU にまで成長し少なくとも原図 9: アウトフローについて

始星形成後の 104 年は生き残る、また、磁気散逸を考えることで、早い段階での円盤形成に繋がること

小さいダストモデルが円盤半径現象を見せている。 が分かっている（いくつかの観測結果からも原始星
これをシミュレーションの図２と照らし合わせると、 形成初期段階で円盤形成が見られる）。

spiral arm の発展の時間と一致している事がわかる。
つまり、円盤半径と角運動量の時間発展によって揺また、観測結果との比較について、図 10 で星と円らぎが成長できるかどうかは、円盤のサイズ成長に

盤の総質量と円盤の半径についての図がある。赤い影響を与えるとわかった。

点が観測結果で、一致しているものもあるが、少し

図９は、アウトフローのサイズ、質量、線形運動量、 外れているものもある。この観測と理論のズレにつ角運動量の時間発展による変化を表した図である。 いて、図 3 や図 5 の spiral arm を観測したものだと

t > 5x104 で円盤角運動量と同じくらい、もしくは大 すると、観測結果が上に外れていることを説明でききくなり、それと同時期に円盤サイズが減少し始め る。つまり、観測結果と整合性がとれているといえている。アウトフローによる角運動量の減少は、円る。

盤のサイズに役割を持っていると改めて考える事が円盤中の磁場の強さの違いやダストモデルの違いは、できる。

長時間の円盤の進化に影響を及ぼすため、本研究の105 年より長いシミュレーションを実行して確認す
る、他のダストモデルでも検証してみる、事が必要

812022 年度 第 52 回 天文・天体物理若手夏の学校がある。

5

Reference

1.Tsukamoto Y., Machida M. N., Susa H., Nomura
H., Inutsuka S., 2020, ApJ, 896, 158

2.Hirashita, H. & Li, Z. Y. 2013, MNRAS, 434,
L70
3.Inoue, T. & Inutsuka, S. 2012, ApJ, 759, 354.
4.Iwasaki, K. & Inutsuka, S. 2011, MNRAS, 418,
1668
5.Iwasaki, K. & Inutsuka, S. 2013, in Astronomical
Society of the Paciﬁc Conference Series, Vol. 474,
Numerical Model- ing of Space Plasma Flows (ASTRONUM2012), ed. N. V. Pogorelov, E. Audit, &G. P. Zank, 239
6.Masson, J., Chabrier, G., Hennebelle, P., Vaytet,
N., & Commerc,on, B. 2016, A&A, 587, A32
7.Mellon, R. R. & Li, Z.-Y. 2008, ApJ, 681, 1356
8.Ormel, C. W., Paszun, D., Dominik, C., & Tielens, A. G. G. M. 2009, A&A, 502, 845.
9.Steinacker, J., Pagani, L., Bacmann, A., & Guieu,
S. 2010, A&A, 511, A9
10.Tomida, K., Okuzumi, S.,& Machida, M. N.
2015, ApJ, 801, 117
11.Tsukamoto, Y., Iwasaki, K., Okuzumi, S.,
Machida, M. N.,& Inutsuka, S. 2015a, ApJ, 810,
L26 ―. 2015b, MNRAS, 452, 278
12.Wurster, J., Price, D. J., & Bate, M. R. 2016,
MNRAS, 457, 1037
13.Wurster, J. & Bate, M. R. 2019, MNRAS, 486,
2587

Acknowledgement
今回、論文のレビューをするにあたり、指導してくださった九州大学の先生、先輩方、同期の皆さん、
ありがとうございました。

82