ブラックホールと半径
このワークシートはMath by Codeの一部です。
前回までにアインシュタイン方程式に到達しました。
この方程式は、非線形で座標系が任意で解釈の自由度もあるという、
ゆるく複雑さの高いものです。
条件をつけて解釈することで、宇宙にまつわる様々な事象を予測する道具にすることができます。
今回は、ブラックホールの存在予測につながったシュワルツシルト解について学ぼう。
１．前提を明確にする。
条件０・アインシュタインの方程式はRij −1/2gij R＝Tij。
条件１・天体の外部が真空、つまり、Tij=0。
条件２・物質分布は球対称、つまり、球座標系。物質は原点に集中している。
とします。
条件０，１から曲率はゼロとなり、R=Rij＝0
条件２から、球面座標系の一般座標をxi=(c, r, t, p)とする。cはctの略、rは半径、tはθ、pはφ。
共変計量テンソル g_ij=diag(A(r), -B(r), -r2, -r2 sin2t )
としましょう。
２．前提からわかるものを求める
＜リッチテンソルを求める＞
リッチテンソルを求めるには、
リーマンテンソルのRi_jklの２つの軸（1番，3番）を縮約します。
つまり、Rij=Σm Rmimjを計算しましょう。
課題：リッチテンソルを求めるにはどうしたらよいでしょうか。
前もって、リーマンテンソルRijklの辞書を作っておきましょう。
さらにその前に、
計算の高速化のためにg_mnからg^mn=g_mn.inv()を求めておきます。
そして、クリ１とクリ２を同じループで作りながら、これもgamma_ijkとして辞書化しておきます。
そうすれば、リーマンテンソルRijklの辞書作りは、
diff(gamma[i][j][l], xs[k])-diff(gamma[i][j][k], xs[l])に
ループをまわして、
term3 = 0
for m in range(dim):
term3 += gamma[i][k][m] * gamma[m][j][l] - gamma[i][l][m] * gamma[m][j][k]
この合計をさらに加えればできます。
#=====================================
import sympy as sp
from sympy import symbols, sin, diff, Matrix, Function, simplify
# 変数と関数の定義
c, r, t, p = symbols('c r t p')
xs = [c, r, t, p]
dim = len(xs)
A = Function('A')(r)
B = Function('B')(r)
# 計量テンソル g_mn
M = sp.diag(A, -B, -r**2, -r**2 * sin(t)**2)
M_inv = M.inv()
# 辞書を使って gamma[i][j][k] = Γ^i_jk の形で保存
gamma = {}
print("計算中: クリストッフェル記号...")
for i in range(dim):
gamma[i] = {}
for j in range(dim):
gamma[i][j] = {}
for k in range(dim):
# 第1種クリストッフェル記号の計算
res2 = 0
for m in range(dim):
# 第1種 Γ_{m,jk} = 1/2 * (g_{mj,k} + g_{mk,j} - g_{jk,m})
gamma1_mjk = 0.5 * (diff(M[m, j], xs[k]) + diff(M[m, k], xs[j]) - diff(M[j, k], xs[m]))
# 第2種 Γ^i_jk = g^im * Γ_{m,jk}
res2 += M_inv[i, m] * gamma1_mjk
gamma[i][j][k] = simplify(res2)
# リーマンテンソルの辞書を初期化・格納
print("計算中: リーマンテンソル R^i_jkl ...")
riemann = {}
for i in range(dim):
riemann[i] = {}
for j in range(dim):
riemann[i][j] = {}
for k in range(dim):
riemann[i][j][k] = {}
for l in range(dim):
# R^i_{jkl} の計算
term1 = diff(gamma[i][j][l], xs[k])
term2 = diff(gamma[i][j][k], xs[l])
term3 = 0
for m in range(dim):
term3 += gamma[i][k][m] * gamma[m][j][l] - gamma[i][l][m] * gamma[m][j][k]
val = simplify(term1 - term2 + term3)
riemann[i][j][k][l] = val
# if val != 0:
# print(f"R^{xs[i]}_{xs[j]}{xs[k]}{xs[l]} = {val}")
# リッチテンソルの計算
print("\n計算中: リッチテンソル R_ij ...")
ricci = {}
for i in range(dim):
ricci[i] = {}
for j in range(dim):
# R_ij = Σ R^m_imj
res_ricci = 0
for m in range(dim):
res_ricci += riemann[m][i][m][j]
final_val = simplify(res_ricci)
ricci[i][j] = final_val
if final_val != 0:
print(f"R_{xs[i]}{xs[j]} = {final_val}")
#=========================================
[OUT]
計算中: クリストッフェル記号...
計算中: リーマンテンソル R^i_jkl ...
計算中: リッチテンソル R_ij ...
R_cc = 0.5*Derivative(A(r), (r, 2))/B(r) - 0.25*Derivative(A(r), r)*Derivative(B(r), r)/B(r)**2 - 0.25*Derivative(A(r), r)**2/(A(r)*B(r)) + 1.0*Derivative(A(r), r)/(r*B(r))
R_rr = -0.5*Derivative(A(r), (r, 2))/A(r) + 0.25*Derivative(A(r), r)*Derivative(B(r), r)/(A(r)*B(r)) + 0.25*Derivative(A(r), r)**2/A(r)**2 + 1.0*Derivative(B(r), r)/(r*B(r))
R_tt = 0.5*r*Derivative(B(r), r)/B(r)**2 - 0.5*r*Derivative(A(r), r)/(A(r)*B(r)) + 1.0 - 1.0/B(r)
R_pp = (0.5*r*A(r)*Derivative(B(r), r) - 0.5*r*B(r)*Derivative(A(r), r) + 1.0*(B(r) - 1)*A(r)*B(r))*sin(t)**2/(A(r)*B(r)**2)
３．リッチテンソル＝０から解を求める
計算結果を見やすくまとめて＝０としましょう。偏微分はすべてｒに関するものだけです。
最初の2式から、
積の微分公式を思いだそう。(AB)'=0となるから、AB=定数。
時空が先々平坦ならA=1,B=1と仮定してよい。AB=1となるから、B=1/A。
これを第３式に代入しよう。
となる。
'r=１なので、(rA)'-1=0とかけるから、rで積分すると、rA-r=C（積分定数）となり、A=1+C/r。
ニュートン力学の意味で天体の質量をMとして、
たとえば、天体を含む半径ｒの領域の体積積分をしてガウスの定理を使うと、
面積分をしてCを-２GM/c^2と求めることができるから、
A= （B=1/A）となる。
これによって、
共変計量テンソル g_ij=diag( , , -r2, -r2 sin2θ )
が決まるね。
＜シュワルツシルト半径とブラックホール＞
A＝０のとき半径ｒ＝ となる。
これをシュワルツシルト半径と呼ぶ。略してrs
A=0のとき、g_cc＝０となり、g_rrは逆数で発散する。
重力天体の半径rがrsより大きく、真空領域を含まないなら問題ない。
重力天体の半径rがrsより小さいとか、大きくても重力分布がrsの内側に集中しているものが、
ブラックホールだ。
ブラックホールの近辺を想像してみよう。
ｒがrsに近いときどうなるだろうか。
g00の1メモリが無限に小さくなる。進むために時間が膨大にかかり、遠くから見て止まって見える。
grrの1メモリは負の無限大に発散する。半径r方向に空間が無限に引き延ばされる。
rsは事象の境界線（地平線）という名前があり、吸い込まれたら見えなく観測できない。
しかし、吸い込まれる粒子から見れば連続的に移動しているだけだ。急激だけれどね。
＜シュワルツシルト半径の具体化＞
G＝6.7×10-11[m^3/ (kg*s^2)]
c=3×108[m/s] から、
1/c2=1.111*10^-15
とするとき、シュワルツシルト半径rsを求めてみよう。rは実半径です。
太陽の場合
M＝2×1030[kg]、r=6.96×108[m]
rs=(2G*1/c^2)*M=2*6.7*1.111(-11-15)*2×1030=(1.489*10-27)*(2*1030)=2978[m] ＝約3㎞
r/rs=6.96*10**8/2978=230000
今の太陽の半径が23万分の１に中心に固まれば、ブラックホール並みになる。
地球の場合
M＝6×1024[kg]、r=6.38×106[m]
rs=(2G*1/c^2)*M=(1.489*10-27)*(6*1024)=8.9*10-3[m] ＝約9ミリ
3.38*10**6/0.0089=380000000
今の地球の半径が38000万分の１に中心に固まれば、ブラックホール並みになる。