はじめに

みなさんこんにちは。物理学科3年の みうら です。いかがお過ごしでしょうか。みなさんもお忙しいでしょうが、私も忙しいです¹。朝起きて ITC-LMS²の更新通知メールで課題が追加されたことを知った時ほど嫌な時はありません。そんな時は洗濯物を強めに洗濯機に入れます。勢いで洗剤を多めに入れないように注意が必要です。

話は変わりますが、5日目には宇宙班の紹介という記事を書いているので、こちらも読んでみてください。

今日は重力崩壊の自己相似的収縮についてお話ししようと思います。これも勉強中ですので、間違いがあったらごめんなさい。そのときは教えてください。

恒星の力学的安定性と断熱指数

以下では、重力崩壊の力学的側面に着目し、ガスを連続体とみなします。すると、恒星のダイナミクスを記述する方程式系は一般に、連続の式

$$\frac{\partial \rho}{\partial t} + \nabla\cdot(\rho \boldsymbol{v}) = 0 \tag{a}$$

Euler方程式

$$\frac{\partial \boldsymbol{v}}{\partial t} + (\boldsymbol{v} \cdot\nabla)\boldsymbol{v} + \frac{1}{\rho}\nabla p + \nabla \phi = 0 \tag{b}$$

エネルギー保存の式

$$\frac{\partial E}{\partial t} + \nabla\cdot(E + p)\boldsymbol{v} + \rho\boldsymbol{v}\cdot(\nabla\phi)=Q \tag{c}$$

です。$\rho~$ は密度、$p~$は圧力、$\boldsymbol{v}~$は速度を表します。また、$E=\frac{1}{2}\rho \boldsymbol{v}^2 + \rho \epsilon~$はガスの全エネルギーの密度を表し（$\epsilon~$は単位質量あたりの内部エネルギー）、$Q~$は光子などガス以外から単位体積単位時間あたりにガスに与えられるエネルギーを表しています。さらに、$\phi$は Poisson 方程式

$$\Delta\phi = 4\pi G\rho \tag{d} $$

を満たしています。

ここからは、星の自転と磁場の効果を無視して、星は球対称であることを仮定します。この仮定のもとで力学平衡状態を考えると、連続の式は自明に満たされます。Euler 方程式からは

$$\frac{dp}{dr} + \rho \frac{d\phi}{dr} = 0 \tag{e} $$

が得られ、Poisson 方程式は

$$\frac{1}{r^2}\frac{d}{dr}\left(r^2\frac{d\phi}{dr}\right) = 4\pi G\rho \tag{f}$$

となります。エネルギー方程式は今回はもう忘れていいです。

式$~(\mathrm{e})~$は圧力勾配力と重力のつり合いを表しています。この式を満たす力学平衡状態が動径方向の小さな擾乱に対して安定かどうかを、次元解析で見ていきましょう。圧縮が素早く起こると思って、断熱と仮定します。半径の変化

$$r \rightarrow r - \delta r \tag{g}$$

による重力ポテンシャルの変化$~\delta\phi~$は

$$\delta\phi \sim -G \frac{m_r}{r^2}\delta r \tag{h}$$

程度です。ここで$~m_r~$は半径$~r~$の内側にあるガスの質量です。一方で、圧力の変化$~\delta p~$は、収縮に伴う密度の上昇

$$\delta\rho \sim \delta\left(\frac{m_r}{r^3}\right) \sim \frac{3m_r}{r^4}\delta r \sim \frac{3\rho}{r} \delta r \tag{i}$$

によって生じて、

$$\delta p = \left(\frac{\partial p}{\partial \rho}\right)_s \delta \rho \tag{j}$$

で与えられます。断熱指数

$$\gamma = \left(\frac{\partial \ln p}{\partial \ln \rho}\right)_s \tag{k}$$

を用いると、

$$\delta p \sim \gamma \frac{p}{\rho}\delta \rho \sim 3\gamma\frac{p}{r}\delta r \tag{l}$$

とかけます。以上を用いて$~\left|\frac{dp}{dr}\right|\sim\frac{p}{r}$と$~\rho \left|\frac{d\phi}{dr}\right| \sim \rho\frac{|\phi|}{r}~$の変化の絶対値を比較します。

$$\left|\delta\left(\frac{dp}{dr}\right)\right| \sim \delta\left(\frac{p}{r}\right) \sim (3\gamma + 1)\frac{p}{r}\frac{\delta r}{r} \tag{m}$$

$$\left|\delta\left(\rho\frac{d\phi}{dr}\right)\right|\sim \left|\delta\rho\left(\frac{Gm_r}{r^2}\right)+\rho \delta\left(\frac{Gm_r}{r^2}\right)\right|\sim5\rho\frac{Gm_r}{r^2}\frac{\delta r}{r} \tag{n}$$

これらと式$~(\mathrm{e})~$の近似式$~\frac{p}{r}\sim\frac{\rho Gm_r}{r^2}$を合わせると

$$3\gamma + 1 = 5 ~~~~\therefore ~\underline{\gamma = \frac{4}{3}~}(=\gamma_\mathrm{c}) \tag{o}$$

を境に様子が異なることがわかります。$\gamma > \gamma_\mathrm{c}~$の場合は圧力の増加が重力の増加に勝るので、最初の擾乱を打ち消す方向に動くことになり、平衡状態は安定であるといえます。一方で$~\gamma < \gamma_\mathrm{c}~$の場合は圧力の増加よりも重力の増加が勝るので、さらに収縮が進み、平衡状態は不安定であるといえます。

以上より、ガスの熱力学状態が、断熱指数を通して力学平衡状態の安定性を決めていることがわかりました。

重力崩壊の自己相似的収縮

ここでは重力崩壊の力学的側面に着目し、電子捕獲反応などは直接扱わないことにします。先ほどと同様に星の自転、磁場の効果は無視し、星は球対称であるとします。さらに、ポリトロープの状態方程式

$$p = K\rho^\gamma \tag{1}$$

を仮定します。$K~$は定数で、$\gamma~$は断熱指数です³。上で見たように$~\gamma < \frac{4}{3}$のとき収縮が進むので、ここでは$~\gamma \lesssim \frac{4}{3}$とします。

以上の仮定のもとで、基礎方程式

$$\frac{\partial \rho}{\partial t} + \frac{1}{r^2}\frac{\partial (r^2\rho v)}{\partial r} = 0 \tag{2}$$

$$\frac{\partial(\rho v)}{\partial t} + \frac{1}{r^2}\frac{\partial(r^2\rho v^2)}{\partial r} + \frac{\partial p}{\partial r} + \rho\frac{\partial \phi}{\partial r} = 0 \tag{3}$$

$$\frac{1}{r^2}\frac{d}{dr}\left(r^2\frac{d\phi}{dr}\right) = 4\pi G\rho \tag{4}$$

の解を求めます。ポリトロープの状態方程式を仮定したのでエネルギー保存の式は要りません。ここで$~v~$は動径方向の速度です。

ここでは解に自己相似性を仮定します。重要な点は、上の方程式系$~(1)\sim(4)~$には次元を持った定数が$~K,~G~$の二つしかないという点です。半径$~r~$と時間$~t~$を組み合わせて無次元量$~X~$が

$$X = K^{-\frac{1}{2}}G^\frac{\gamma - 1}{2}r(-t)^{\gamma - 2} \tag{5}$$

のようにつくれます⁴。この$~X~$のみに依存する解が存在することを期待します。収縮する過程を考えており、一点に収縮する時間を$~t=0~$にとるために$~t~$の前にマイナスがついています。

各物理量の次元を考慮して、次のような形を仮定します：

$$\rho = G^{-1}(-t)^{-2}D(X) \tag{6}$$

$$v = K^\frac{1}{2}G^\frac{1-\gamma}{2}(-t)^{1-\gamma}V(X) \tag{7}$$

$$m = K^\frac{3}{2}G^\frac{1-3\gamma}{2}(-t)^{4-3\gamma}M(X) \tag{8}$$

ここでそれぞれの大文字の量は無次元の変数で$~X~$のみに依存します。また$~M(X)~$は

$$M(X) = 4\pi \int_0^X dx~x^2 D(X) \tag{9}$$

で定義される量です。

式$~(6)\sim(9)$を見ると、物理量の$~t~$依存性は$~X~$からくるものと物理量の次元を表す因子からくるものがあることがわかります。式$~(6)~$をパッとみると$~t\rightarrow0~$のとき密度$~\rho~$が発散してしまいますが、$~D(X)~$が$~D(X)\varpropto X^{-2/(2-\gamma)}$という$~X~$依存性を持っていれば$~t~$依存性が打ち消しあい、このような特異な振る舞いをしないことがわかります。このような議論を式$~(6)\sim(8)~$について行うと

$$D(X) \varpropto X^{-\frac{2}{2-\gamma}}~~~(X \gg 1)\tag{10}$$

$$V(X) \varpropto X^\frac{1-\gamma}{2-\gamma}~~~(X \gg 1) \tag{11}$$

$$M(X) \varpropto X^\frac{4-3\gamma}{2-\gamma}~~~(X \gg 1) \tag{12}$$

とわかります。

式$~(1),~(6)\sim(8)~$を連続の式$~(2)~$と Euler 方程式$~(3)~$に代入することで

$${V + (2 - \gamma)X}\frac{D'}{D}+V' = -2 - \frac{2V}{X} \tag{13}$$

$$\gamma D^{\gamma - 2}D' + {V + (2-\gamma)X}V' = -\frac{M}{X^2}-(\gamma-1)V \tag{14}$$

が得られます。記号$~'~$は$~X~$での微分を表します。

ここで、収縮速度を homologous な収縮⁵をする部分とそれ以外の部分に分け、その部分を$~U~$とします：

$$V = (\gamma - 2)X + U \tag{15}$$

実際、$~X=~$一定となる点$~r(t)$ の速度は

$$\left.\frac{dr(t)}{dt}\right|_X = K^\frac{1}{2}G^\frac{1-\gamma}{2}(-t)^{1-\gamma}(\gamma-2)X \tag{16}$$

であり、これと式$~(7)$を比較すれば、$~U~$がその差に対応していることがわかります。この式$~(15)~$を式$~(13),~(14)~$に代入すると

$$U\frac{D'}{D} + U'= 4-3\gamma -\frac{2U}{X} \tag{17}$$

$$\gamma D^{\gamma-2}D' + UU' = -\frac{M}{X^2} + (\gamma-1)(2-\gamma)X+(3 - 2\gamma)U \tag{18}$$

となります。式$~(17)~$は積分できて

$$4\pi X^2DU = (4-3\gamma)M \tag{19}$$

となります。$X \ll 1$の極限では$~M(X)\sim \frac{4\pi}{3}DX^3~$なので、

$$U = \left(\frac{4}{3} - \gamma\right)X ~~~(X \ll 1)\tag{20}$$

$$V = -\frac{2}{3}X ~~~(X \ll 1) \tag{21}$$

が得られます。ここで無次元の音速

$$C_\mathrm{s} = \gamma^\frac{1}{2}D^\frac{\gamma-1}{2} \tag{22}$$

を導入します。この$~C\mathrm{s}~$は$~X\rightarrow 0~$で一定値に近づくので、$~X\ll 1~$では$~C_\mathrm{s} > V~$となり、この領域の収縮速度は亜音速です。一方で、$X \gg 1~$では式$~(11),~(15)~$から$~U\rightarrow \infty~$であり、式$~(10)~$から

$$C_\mathrm{s} \varpropto X^\frac{1-\gamma}{2-\gamma}~~~(X \gg 1) \tag{23}$$

となるので、この領域では$~U~$が超音速です。実際に比較すべきは$~v~$と真の音速$~c_\mathrm{s}~$ですが、この意味でも内部コアは亜音速、外部コアは超音速になります。

内部コアで亜音速、外部コアで超音速ということは、解は遷音速点を滑らかに通過するように決めなければなりません。この（$U~$の）遷音速点は式$~(17),~(18)~$の特異点になっています。そのような望ましい解は式$~(9),~(19)~$を用いて式$~(18)~$を数値積分すれば得られます。無次元の中心密度を試行パラメータとして解を探すらしいです[3]。得られた速度分布の図を載せたかったのですが、何かと面倒なので、気になった方は[1]か[3]を見てください。

反跳と衝撃波

重力崩壊する大質量星のコアは、密度が原子核密度$~\rho_0 \sim 10^{14}~\mathrm{g~cm^{-3}}~$を超えると急激に硬くなり、断熱指数が$~\gamma \gtrsim 2~$となるため、重力崩壊に対する安定性を回復します。それによって急激に収縮が止まり、そして膨張に転じます。これを反跳といいます。上で見たように、内部コアは亜音速で収縮しているので、反跳は内部コアで一斉に起こります。一方で外部コアは超音速で収縮しているので、反跳の情報は外部コアには届きません。その結果、膨張しようとする内部コアと収縮しようとする外部コアがその境界でぶつかり合い、衝撃波が形成されます。衝撃波は不連続面で、物質が衝撃波面を通過すると収縮が瞬時に膨張に転じます。内部コアがピストンのようにはたらき、衝撃波は外部コア内を外に向かって進んでいきます。そしてこの衝撃波が外部コア、さらに星の外層部を伝播し星の表面に達したとき、光学観測的意味で超新星爆発が起こったことになります。しかし衝撃波は外部コア内で失速し（停滞衝撃波になる）、そのまま超新星爆発を起こすことはできないと考えられています。停滞衝撃波の復活機構としてはニュートリノ加熱などがあります。

ここで述べたことはまだ詳しく勉強していません。これから勉強したいと思っていますが、如何せん忙しいので... ⁶。

おわりに

アドベントカレンダーということで何かクリスマス要素を探したところ、参考文献[1]の本の第１版第１刷発行日が12月25日でした。めでたいですね。また、ガンマ線バーストという現象は１日に１回程度の頻度で観測されているのですが、2010年のクリスマスに観測されたものは特に異常な現象だったらしく、クリスマスガンマ線バーストと呼ばれているそうです。おしゃれですね。

今日はここまでです。アドベントカレンダーはまだまだ続きますので、お楽しみください。

参考文献

[1] :「超新星 新天文学ライブラリー 第4巻」 山田章一 (日本評論社, 2016)
[2] : 「宇宙物理学：星・銀河・宇宙論」 高原文郎 (朝倉書店, 2015)
[3] : A. Yahil 1983, Astrophys. J. 265, 1047
[4] : nature ダイジェスト 「不可解な『クリスマスガンマ線バースト』」 https://www.natureasia.com/ja-jp/ndigest/v9/n3/不可解な「クリスマスガンマ線バースト」/43027

1. 中性子星の密度と比べればスカスカなスケジュールです。 ↩

2. Information Technology Center - Learning Management System の略で、教員が講義資料や課題を配布したり学生が課題を提出したりするシステムです。ちなみに授業がオンラインになる前は "ITC-LMS" という名前をちゃんと把握している教員はいなかった気がします。 ↩

3. 参考文献[1]には、「断熱指数に対応するもの」と書いてあり、どう違うのか調べたのですがよくわかりませんでした。 ↩

4. 実際には初期中心密度など他にも次元を持つ量は存在しますが、これらの影響は一時的なものに過ぎないらしいです。 ↩

5. 速度が半径に比例した崩壊のことを homologous な崩壊といいます。 ↩

6. 中性子星の密度と比べればスカスカなスケジュールです。ちなみにこのコア崩壊型超新星爆発ののちに形成されるのが中性子星です。 ↩