Return to EinsteinPy や ctensor でアインシュタイン方程式を計算する

EinsteinPy と SymPy でアインシュタイン方程式を解いてシュバルツシルト解を求める

EinsteinPy を使って球対称な計量から真空のアインシュタイン方程式

$$G^{\mu}_{\ \ \nu} = 0$$

を解き,シュバルツシルト解を求める。

 

ライブラリの import

In [1]:
from sympy.abc import *
from sympy import *
from einsteinpy.symbolic import *
init_printing()

球対称な計量

ランダウ・リフシッツ「場の古典論」の記述にそって(しかし,signature は $(-, +, +, +)$ にして),球対称な計量を以下のようにおく。(ギリシャ文字の $\lambda$ は lambda ではなく,sympy.abc の中で lamda として定義されている。lambda は Python の予約語,変数として使用不可。)

$$ds^2 = -e^{\nu(t,r)} dt^2 + e^{\lambda(t,r)} dr^2 + r^2 \left(d\theta^2 + \sin^2\theta\,d\phi^2 \right)$$

In [2]:
# from sympy.abc import * すると使える
lamda
Out[2]:
$\displaystyle \lambda$
In [3]:
lamda = Function('lamda')(t, r)
nu = Function('nu')(t, r)

Metric = diag(-exp(nu), exp(lamda), r**2, r**2 * sin(theta)**2).tolist()
g = MetricTensor(Metric, [t, r, theta, phi])
g.tensor()
Out[3]:
$\displaystyle \left[\begin{matrix} -e^{\nu{\left(t,r \right)}} & 0 & 0 & 0\\0 & e^{\lambda{\left(t,r \right)}} & 0 & 0\\0 & 0 & r^{2} & 0\\0 & 0 & 0 & r^{2} \sin^{2}{\left(\theta \right)}\end{matrix}\right]$

アインシュタイン・テンソル

$\displaystyle G^{\mu}_{\ \ \nu} = R^{\mu}_{\ \ \nu}  -\frac{1}{2} R \delta^{\mu}_{\ \ \nu} $ = ein とおく。(.change_config('ul') で上付下付に)

In [4]:
ein=EinsteinTensor.from_metric(g).change_config('ul')

$\lambda$ は時間に依存しないこと

$\displaystyle G^{1}_{\ \ 0}$

In [5]:
ein[1,0]
Out[5]:
$\displaystyle \frac{e^{ -\lambda{\left(t,r \right)}} \frac{\partial}{\partial t} \lambda{\left(t,r \right)}}{r}$

$\displaystyle G^{1}_{\ \ 0} = 0$ より

$$\frac{\partial \lambda}{\partial t} = 0, \quad \therefore\ \ \lambda(t, r) \Rightarrow \lambda(r)$$

$\lambda(r)$ として,あらためてアインシュタイン・テンソルを求めてみる。

In [6]:
lamda = Function('lamda')(r)
nu = Function('nu')(t, r)

Metric = diag(-exp(nu), exp(lamda), r**2, r**2 * sin(theta)**2).tolist()
g = MetricTensor(Metric, [t, r, theta, phi])
g.tensor()
Out[6]:
$\displaystyle \left[\begin{matrix} -e^{\nu{\left(t,r \right)}} & 0 & 0 & 0\\0 & e^{\lambda{\left(r \right)}} & 0 & 0\\0 & 0 & r^{2} & 0\\0 & 0 & 0 & r^{2} \sin^{2}{\left(\theta \right)}\end{matrix}\right]$
In [7]:
ein=EinsteinTensor.from_metric(g).change_config('ul')

$\displaystyle G^{0}_{\ \ 0}$

In [8]:
ein[0,0]
Out[8]:
$\displaystyle \frac{\left( -r \frac{d}{d r} \lambda{\left(r \right)}  -e^{\lambda{\left(r \right)}} + 1\right) e^{ -\lambda{\left(r \right)}}}{r^{2}}$

$\displaystyle G^{1}_{\ \ 1}$

In [9]:
ein[1,1]
Out[9]:
$\displaystyle \frac{\left(r \frac{\partial}{\partial r} \nu{\left(t,r \right)}  -e^{\lambda{\left(r \right)}} + 1\right) e^{ -\lambda{\left(r \right)}}}{r^{2}}$

$\displaystyle G^{2}_{\ \ 2}$

In [10]:
ein[2,2]
Out[10]:
$\displaystyle \frac{\left( -r \frac{d}{d r} \lambda{\left(r \right)} \frac{\partial}{\partial r} \nu{\left(t,r \right)} + r \left(\frac{\partial}{\partial r} \nu{\left(t,r \right)}\right)^{2} + 2 r \frac{\partial^{2}}{\partial r^{2}} \nu{\left(t,r \right)}  -2 \frac{d}{d r} \lambda{\left(r \right)} + 2 \frac{\partial}{\partial r} \nu{\left(t,r \right)}\right) e^{ -\lambda{\left(r \right)}}}{4 r}$

$\displaystyle G^{3}_{\ \ 3}$

In [11]:
ein[3,3]
Out[11]:
$\displaystyle \frac{\left( -r \frac{d}{d r} \lambda{\left(r \right)} \frac{\partial}{\partial r} \nu{\left(t,r \right)} + r \left(\frac{\partial}{\partial r} \nu{\left(t,r \right)}\right)^{2} + 2 r \frac{\partial^{2}}{\partial r^{2}} \nu{\left(t,r \right)}  -2 \frac{d}{d r} \lambda{\left(r \right)} + 2 \frac{\partial}{\partial r} \nu{\left(t,r \right)}\right) e^{ -\lambda{\left(r \right)}}}{4 r}$

$\displaystyle G^{2}_{\ \ 2} = \displaystyle G^{3}_{\ \ 3}$ であることを確認。$\displaystyle G^{2}_{\ \ 2}  -\displaystyle G^{3}_{\ \ 3}$ がゼロとなることを示す。

In [12]:
ein[2,2] - ein[3,3]
Out[12]:
$\displaystyle 0$

$\nu =  -\lambda$ とおけること

In [13]:
simplify(ein[1,1]-ein[0,0])
Out[13]:
$\displaystyle \frac{\left(\frac{d}{d r} \lambda{\left(r \right)} + \frac{\partial}{\partial r} \nu{\left(t,r \right)}\right) e^{ -\lambda{\left(r \right)}}}{r}$

$G^1_{\ \ 1}  -G^2_{\ \ 2} = 0$ より,以下の式が得られる。

$$\frac{\partial}{\partial r}\left(\lambda(r) + \nu(t, r) \right) = 0$$

これから,

$$\nu(t, r) =  -\lambda(r) + f(t)$$

となる。

\begin{eqnarray}
\therefore\ \ e^{\nu(t, r)} dt^2 &=& e^{ -\lambda(r) + f(t)} dt^2 \\
&=& e^{ -\lambda(r)} \left( e^{\frac{f(t)}{2}} dt\right)^2
\end{eqnarray}

時間 $t$ のみの任意関数 $f(t)$ の自由度は,$e^{\frac{f(t)}{2}} dt \Rightarrow dt’$ なる新しい時間座標の定義によって吸収できるので,一般性を失うことなく $f(t) = 0$ すなわち

$$\nu(t, r) \Rightarrow  -\lambda(r)$$

とすることができる。

バーコフの定理

ここまでは,球対称真空解は metric が時間によらない,つまり静的であるということを示しているわけで,バーコフの定理の証明になっている。バーコフの定理のもう一つの帰結である漸近的平坦性については,以下で示すように解が $e^{-\mu(r)} = \displaystyle 1  -\frac{r_g}{r}$ となることで $r \rightarrow \infty$ で $e^{-\mu(r)} \rightarrow 1$ となることからわかる。

ということで,あらためて以下のような計量テンソルに対して,アインシュタイン・テンソルを計算してみる。

In [14]:
lamda = Function('lamda')(r)

Metric = diag(-exp(-lamda), exp(lamda), r**2, r**2 * sin(theta)**2).tolist()
g = MetricTensor(Metric, [t, r, theta, phi])
g.tensor()
Out[14]:
$\displaystyle \left[\begin{matrix} -e^{ -\lambda{\left(r \right)}} & 0 & 0 & 0\\0 & e^{\lambda{\left(r \right)}} & 0 & 0\\0 & 0 & r^{2} & 0\\0 & 0 & 0 & r^{2} \sin^{2}{\left(\theta \right)}\end{matrix}\right]$
In [15]:
ein=EinsteinTensor.from_metric(g).change_config('ul')
In [16]:
ein[0,0]
Out[16]:
$\displaystyle \frac{\left( -r \frac{d}{d r} \lambda{\left(r \right)}  -e^{\lambda{\left(r \right)}} + 1\right) e^{ -\lambda{\left(r \right)}}}{r^{2}}$
In [17]:
eq = (ein[0,0] * r**2).expand()
eq 
Out[17]:
$\displaystyle  -r e^{ -\lambda{\left(r \right)}} \frac{d}{d r} \lambda{\left(r \right)}  -1 + e^{ -\lambda{\left(r \right)}}$

微分方程式を解き,$\lambda(r)$ を求める

微分方程式 $\displaystyle  -r \frac{d}{d r} \lambda{\left(r \right)}  -e^{\lambda{\left(r \right)}} + 1 = 0$ を SymPy の dsolve() を使って解く。

In [18]:
sol = dsolve(eq, lamda)
sol
Out[18]:
$\displaystyle  -\lambda{\left(r \right)} + \log{\left(r \right)} + \log{\left(e^{\lambda{\left(r \right)}}  -1 \right)} = C_{1}$

もうちょっとのところなので,

$$e^{-\lambda} \equiv f, \quad \lambda \Rightarrow -\log f$$

とおいて(.subs(lamda, -log(f)))解いてもらう。

In [19]:
print('元の式は')
display(sol)

print('f(r) で書き直すと,解くべき方程式は')
f = Function('f')(r)
eq2 = sol.subs(lamda, -log(f))
display(eq2)

print('f(r) について解いた解は')
sol2 = solve(eq2, f)
sol2
元の式は
$\displaystyle  -\lambda{\left(r \right)} + \log{\left(r \right)} + \log{\left(e^{\lambda{\left(r \right)}}  -1 \right)} = C_{1}$
f(r) で書き直すと,解くべき方程式は
$\displaystyle \log{\left(r \right)} + \log{\left(-1 + \frac{1}{f{\left(r \right)}} \right)} + \log{\left(f{\left(r \right)} \right)} = C_{1}$
f(r) について解いた解は
Out[19]:
$\displaystyle \left[ \frac{r  -e^{C_{1}}}{r}\right]$
In [20]:
sol2[0].expand()
Out[20]:
$\displaystyle 1  -\frac{e^{C_{1}}}{r}$

積分定数 $C_1$ はニュートン近似のときに,

\begin{eqnarray}
g_{00} &=&  -e^{-\lambda(r)} \\
&\simeq&  -\left( 1 + 2 \frac{\phi}{c^2}\right) \\
&=&  -\left( 1  -2 \frac{GM}{r c^2}\right)
\end{eqnarray}

となることから,

$$e^{C_1} = \frac{2 GM}{c^2} \equiv r_g$$

となる。

最終的に

$$ds^2 = -\left(1  -\frac{r_g}{r} \right) dt^2 + \frac{dr^2}{1  -\frac{r_g}{r}} + r^2 \left(d\theta^2 + \sin^2\theta \,d\phi^2 \right)$$