宇宙状態の方程式って何？
このワークシートはMath by Codeの一部です。
アインシュタイン方程式は、非線形で座標系が任意で解釈の自由度もあるという、
ゆるく複雑さの高いものです。
条件をつけて解釈することで、宇宙にまつわる様々な事象を予測する道具にすることができます。
前回は、ブラックホールの存在予測につながったシュワルツシルト解について学びました。
今回は、宇宙項があって、ふくらむアインシュタイン方程式の解をさぐろう。
１．前提を明確にする。
条件１・一様で等方的空間。
条件２・宇宙流体は静止しているので、エネルギーテンソルTij=diag(ρ,p,p,p)
条件３・アインシュタインの方程式は
Rij ー1/2gij R（ーΛgij ）＝8πGTij／c4。
条件４・宇宙は膨張・縮小するとしても一様に伸縮する。
とします。
条件１から、球面座標系をxi=(r, θ, φ)とする。
空間部分の計量テンソル はg_ij=diag(g(r), r2, r2 sin2t ),
g^ij=diag(1/g(r),1/r2, 1/r2 sin2θ )とおけるね。
課題：等方的空間のリッチテンソルをもとめよう。
シュワルツシルト半径を求めたときのように、リッチテンソルを求めるコードを実行すればよいね。import sympy as sp
from sympy import symbols, sin, diff, Matrix, Function, simplify
# 変数と関数の定義
r, t, p = symbols('r t p')
xs = [r, t, p]
dim = len(xs)
g = Function('g')(r)
# 計量テンソル g_mn
M = sp.diag(g, 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_ij ...
R_rr = 1.0*Derivative(g(r), r)/(r*g(r))
R_tt = 0.5*r*Derivative(g(r), r)/g(r)**2 + 1.0 - 1.0/g(r)
R_pp = (0.5*r*Derivative(g(r), r) + 1.0*(g(r) - 1)*g(r))*sin(t)**2/g(r)**2
となる。
つまり、
計量テンソルの対角成分以外は０なので、この３つのテンソルに対応する対角成分をかけた積和が
スカラー曲率になるね。
R=simplify(M_inv[0,0]*ricci[0][0] + M_inv[1,1]*ricci[1][1] + M_inv[2,2]*ricci[2][2])
print(R)
#================================================
[OUT]
2.0*Derivative(g(r), r)/(r*g(r)**2) + 2.0/r**2 - 2.0/(r**2*g(r))
これをさらにカンタンにすると、2.0/r**2{1+ rg'/g2-1/g}
= ＝C（定数）とすると、曲率が一定になる。
条件１からこの曲率が定数ならば、一様で等方的な空間になる。
積分して、g=1/(1-C/6 * r^2+ const/r)
原点を任意にできるようにするためと、Cとrの次元にそろえるために次のようにおこう。
g=1/(1-kr2) ただし、k=0,1,-1のどれかにする。
これで、条件１の一様で等方的な空間の設定になった。
空間部分の計量テンソル は
g_ij=diag(1/(1-k r2), r2, r2 sin2t )
r基底のメモリに着目しよう。
k=0ではメモリは常に１だから、平坦空間になる。
k=-1ではメモリは必ず正だが、原点から遠くなるとメモリが巨大数分の１となり、
微細になる。だから、逆に距離は極度に大きくなり、果てのない開いた宇宙になる。
k=1ではメモリrは１以内でしかありえないために、閉じた球状宇宙となる。
条件１から設定したｋは宇宙の形を決めた。
次は、条件４・宇宙は膨張・縮小するとしても一様に伸縮する。
を考えてみよう。
条件５から、空間の計量の線素をdl2とすると、dlをa倍する相似係数（宇宙のスケール）を考えて、
しかも、これが時間とともに変わる関数a(t)としたい。
その4次元の一般座標の線素は
ds2=(cdt)2-(dl)2=c2dt2-a2(t)[1/(1-k r2) dr2+r2dθ2+ r2 sin2tdφ2]
とおける。
これが、フリードマン・ルメートル・ロバートソン・ウォーカ計量（FLRW計量）です。
課題：FLRW計量のリッチテンソルとスカラー曲率を求めるにはどうしたらよいですか。
変数名がかぶりやすいので、今まで1文字でtheta,phiをt,pとしていたものを略さずにかきます。
その代わり、時間の基底変数はtではなくctにしましょう。
自然単位系ではないので、光速はｃをつかいますが、微分の関係で時間の計量を１にして、
ct自体を１つの変数として微分しておき、最後にtで微分することで、時間微分表示に直しましょう。
微分の連鎖率により、
# dot{a}=da(ct)/dt= da(ct)/dct * dct/dt=a'*c
# ddot{a}=d(a'*c/dt)=d(a'(ct)/dt*c)+a'*dc/dt=a''*c2
# だから、
# a'(ct) = (1/c) * dot{a}, a''(ct) = (1/c^2) * ddot{a}
# 表示を分かりやすくするため、a1=dot{a}, a2=ddot{a} として代入表示
というのがそのやり方です。
import sympy as sp
from sympy import symbols, sin, diff, Function, simplify
# 1. 変数と座標の定義
# 微分は文字通り (ct) で行う
ct, r, theta, phi = symbols('ct r theta phi')
c, k = symbols('c k')
xs = [ct, r, theta, phi]
dim = len(xs)
a = Function('a')(ct)
# 2. 標準的な FLRW 計量テンソル (Signature: + - - -)をそのまま使う。
# ds^2 = 1(dct)^2 - a(ct)^2 [ dr^2/(1-kr^2) + r^2 dθ^2 + r^2 sin^2θ dφ^2 ]
M = sp.diag(1,
-a**2 / (1 - k * r**2),
-a**2 * r**2,
-a**2 * r**2 * sin(theta)**2)
M_inv = M.inv()
# 3. クリストッフェル記号の計算
gamma = {}
for i in range(dim):
gamma[i] = {}
for j in range(dim):
gamma[i][j] = {}
for l in range(dim):
res = 0
for m in range(dim):
g1 = 0.5 * (diff(M[m, j], xs[l]) + diff(M[m, l], xs[j]) - diff(M[j, l], xs[m]))
res += M_inv[i, m] * g1
gamma[i][j][l] = simplify(res)
# 4. リッチテンソル R_ij の計算
print("非ゼロのricciテンソル表示：")
ricci = []
for i in range(dim):
row = []
for j in range(dim):
res_ricci = 0
for m in range(dim):
# R_ij = ∂_m Γ^m_ij - ∂_j Γ^m_im + Γ^m_mn Γ^n_ij - Γ^m_jn Γ^n_im
term1 = diff(gamma[m][i][j], xs[m])
term2 = diff(gamma[m][i][m], xs[j])
term3 = sum(gamma[m][m][n] * gamma[n][i][j] for n in range(dim))
term4 = sum(gamma[m][j][n] * gamma[n][i][m] for n in range(dim))
res_ricci += term1 - term2 + term3 - term4
row.append(simplify(res_ricci))
if res_ricci != 0:
simple_ricci = res_ricci.subs({diff(a, ct, 2): a2/c**2, diff(a, ct): a1/c})
print(f"R_{xs[i]}{xs[j]} = {simple_ricci}")
ricci.append(row)
# 5. スカラー曲率 R = g^ij R_ij の計算
R_val = 0
for i in range(dim):
R_val += M_inv[i, i] * ricci[i][i]
R_final = simplify(R_val)
#print("--- 導出結果 ---")
#print(f"スカラー曲率 R (ct微分版): {R_final}")
# 6. 時間微分 (t) への書き換え
# dot{a}=da(ct)/dt= da(ct)/dct * dct/dt=a'*c
# ddot{a}=d(a'*c/dt)=d(a'(ct)/dt*c)+a'*dc/dt=a''*c2
# だから、
# a'(ct) = (1/c) * dot{a}, a''(ct) = (1/c^2) * ddot{a}
# 表示を分かりやすくするため、a1=dot{a}, a2=ddot{a} として代入表示
a1, a2 = symbols('a1 a2')
R_dot_notation = R_final.subs({diff(a, ct, 2): a2/c**2, diff(a, ct): a1/c})
print(f"スカラー曲率 R (時間微分ドット版): {simplify(R_dot_notation)}")
#=====================================
[OUT]非ゼロのricciテンソル表示：
R_ctct = -3.0*a2/(c**2*a(ct))
R_rr = -2.0*a1**2/(c**2*(k*r**2 - 1)) - 1.0*a2*a(ct)/(c**2*(k*r**2 - 1)) - 2.0*k/(k*r**2 - 1)
R_thetatheta = 2.0*a1**2*r**2/c**2 + 1.0*a2*r**2*a(ct)/c**2 + 2.0*k*r**2 - 1.0*(-tan(theta)**2 - 1)/tan(theta)**2 - 1.0 - 1.0/tan(theta)**2
R_phiphi = 2.0*a1**2*r**2*sin(theta)**2/c**2 + 1.0*a2*r**2*a(ct)*sin(theta)**2/c**2 + 1.0*k*r**2*sin(theta)**2 + 1.0*(k*r**2 - 1)*sin(theta)**2 + 0.5*sin(2*theta)/tan(theta) - 1.0*cos(2*theta)
スカラー曲率 R (時間微分ドット版): 6.0*(-a1**2 - a2*a(ct) - c**2*k)/(c**2*a(ct)**2)
時間微分はダッシュでなくドットで表示しますが、geogebraの表記機能の関係からダッシュにします。
これらを見ると、
時間基底での曲率は「膨張の加速度a''」に連動していることがわかります。
スカラーの曲率は、「膨張の速度a'、その加速度a''」さらに「宇宙の形係数ｋ」の混ざり合い
になってますね。
２．フリードマン方程式から宇宙状態方程式へ
リッチテンソルRttとスカラー曲率Rが求められたから、
アインシュタイン方程式に代入できるね。
条件２・宇宙流体は静止しているので、エネルギーテンソルTij=diag(ρ,p,p,p)
（ρはエネルギー密度、pは運動量。）
条件３・アインシュタインの方程式は
Rij ー1/2gij RーΛgij ＝8πGTij／c4。
条件３の時間座標から、R00-1/2g00R-Λg00 =8πGT00／c4
時間座標ctの成分だけで、代入してみよう。
・g00=M00=1
・条件２の時間座標T00=ρ
条件１，５から得た、
・
・
を入れると、
アインシュタイン方程式は、次のようになるね。
式を整理すると
両辺にc²をかけて、
これがフリードマン方程式だ。
膨張の加速度a''が消えました。
T00=ρエネルギー密度ではなく、T00=ρc²質量密度として、単位系を整理して、
K=c2kとおこう。Kを宇宙曲率と呼ぶことがある。
第１方程式（宇宙のエネルギー）
これが宇宙膨張のフリードマン方程式だ。
これを時間で微分すると、
計算の詳細は省略するが、Tijについても、同様に代入と整理をすることで、
が得られる。第２方程式（宇宙空間の加速度）。
pは運動量、つまり、圧力pだ。圧力ｐのあるもの（光、ガス）が多いほど宇宙の膨張に負の影響、ブレーキがかかるということだね。
この2式の左辺が同じだから、比べることから次の式ができるね。
が得られる。
しかし、圧力pと質量エネルギーρの関係式がこの1個では、宇宙の状態構造は決まらない。
＜宇宙状態方程式＞
圧力pと質量エネルギーρの関係式、
状態方程式p=ωρによって、宇宙状態方程式ができるでしょう。
代表的な例をあげよう。
・p≪ρ（非相対論的物質中心モデル）
光、ガスによる圧力が少なく、塵や星の重力のみで宇宙膨張にブレーキがかかる宇宙。
ω＝０つまり、p=0の場合の宇宙状態で、ρ∝a-3。
K=0、Λ＝０は平坦で非相対論的物質で満ちた宇宙モデル（アインシュタインドジッター宇宙）がある。
・p=ρ/3（相対論的物質モデル）
ω＝1/3の場合。
光、ガスの圧力による体積増加の他に、波長の伸びでも宇宙膨張にブレーキがかかる宇宙。
m=0でもρ∝a-4もありうる。
・p=-ρ（真空モデル）
ω＝-1の場合。ρが一定で膨張しても密度が減らない。
物質がなく宇宙定数Λだけが存在し、指数関数的に膨張する宇宙（ジッター宇宙）
一般には、K、Λともに非ゼロで、非相対論的物質と相対論的物質の両方が存在して存在比をパラメータ化
したものがフリードマン・ルメートル宇宙。これは固定ではなく可変な宇宙像だ。
パラメータそのものを宇宙の成長とともに変えることで、宇宙の歴史を記述しようとするのが、今の標準宇宙論になっているようだ。