如何用Mathematica做張量計算

05-30

如何用Mathematica做張量計算

來自專欄天文筆記

Mathematica可謂是符號計算的神器。算什麼東西的時候，可以先用mathematica算一下，如果可以算出來，再用手算，這樣比較省時間。

上一篇文章

李剛：Gravity by Hartle 筆記：真空中的場方程?

zhuanlan.zhihu.com

中我們學習了真空中的場方程。對於四維時空來說，其度規有16個分量，克里斯托夫符號有64個分量，黎曼張量有256個分量，里奇張量有16個分量，計算量不可謂不大。

手算是不可能手算的，這輩子是不可能手算的，作業又不會，只能用mathematica才能維持得了及格這樣子。。

所以我們還是來看看怎麼用mathematica計算張量吧。。

我們先聲明好坐標的變數和度規，比如用 $(t,r, heta,phi)$ 和施瓦西度規來做個例子

v={t,r,theta, phi}L1={(1-2M/r), (1-2M/r)^(-1), r^2, r^2*(Sin[theta])^2}len=Dimensions[v][[1]];(*number of dimension*)g=DiagonalMatrix[ L1 ];(*metric*)Print["The metric is"]Print[g//MatrixForm](*其中 {a,b,c,e..}組成了一個表格，表格可以擴充為張量*)

有了度規，我們還要算其逆，可以這樣

invg=Simplify[Inverse[ g ]];(*inverse metric*)Print["The inverse metric is"];Print[Simplify[invg]//MatrixForm];

然後我們要用 $Gamma^delta_{etagamma}=g^{alphadelta}left(frac{partial g_{alphaeta}}{partial x^{gamma}}+frac{partial g_{alphagamma}}{partial x^eta}-frac{partial g_{etagamma}}{partial x^alpha} ight)$ 來算克里斯托夫符號

dg1=Outer[D, g, v]; (*Partial devirative of the metric with respect to v*)dg2=Transpose[dg1, {1,3,2}];(*change indices orders*)dg3=Transpose[dg1, {2,3,1}];G=Simplify[(1/2)invg . (dg1+dg2-dg3)];(*christoffel symbols*)

其中Outer[D, g, v]可以計算右邊括弧里第一項 $frac{partial g_{alphaeta}}{partial x^{gamma}}$ ，然後我們用Transpose交換角標的順序，於是就有了 $frac{partial g_{alphagamma}}{partial x^eta}$ 和 $frac{partial g_{etagamma}}{partial x^alpha}$ 。Simplify能化簡表達式。

然後我們就能算那256個量的黎曼張量了

$R^{alpha}_{~etagammadelta}=frac{partial Gamma^alpha_{etadelta}}{partial x^gamma}-frac{partial Gamma^alpha_{etagamma}}{partial x^delta}+Gamma^alpha_{gammaepsilon}Gamma^{epsilon}_{etadelta}-Gamma^alpha_{deltaepsilon}Gamma^{epsilon}_{etagamma}$

dG=Outer[ D, G, v ];(*Partial devirative of the Christoffel symbols with respect to v*)(*Riemann tensor*)Rie1=dG;Rie2=Transpose[Rie1,{1,2,4,3}];Rie3=Table[ Sum[G[[i,k,l]]G[[l,j,h]],{l, len}] ,{i,len},{k,len},{j,len},{h,len} ];Rie4=Table[ Sum[G[[i,h,l]]G[[l,j,k]],{l, len}] ,{i,len},{k,len},{j,len},{h,len} ];Riemann=Simplify[Rie1-Rie2+Rie3-Rie4];Print["The Riemann tensor is"];Print[Simplify[Riemann]//MatrixForm];

我們又用了Outer[]算出了 $frac{partial Gamma^alpha_{etadelta}}{partial x^gamma}$ 這個偏導，然後用Transpose[]交換下角標得到 $frac{partial Gamma^alpha_{etagamma}}{partial x^delta}$ 。

然後有個新命令Sum[]，Sum[G[[i,k,l]]G[[l,j,h]],{l, len}]就是 $Gamma^alpha_{gammaepsilon}Gamma^{epsilon}_{etadelta}$ 對 $epsilon$ 求和，但是別的指標是個定值，然後再用Table[]擴充成張量，得到 $Gamma^alpha_{gammaepsilon}Gamma^{epsilon}_{etadelta}$ 所有可能的值。

於是我們就算好了黎曼張量。。。

接下來算下里奇張量 $R_{alphaeta}equiv R^{gamma}_{~alphagammaeta}$ ，我們就用Sum命令和Table命令，把黎曼張量的第一個指標和第三個指標加起來即可。

Ruv=Table[Sum[Riemann[[k,i,k,j]],{k,len}],{i,len},{j,len} ]

然後我們就算完啦～手算是不可能手算的，這輩子不可能手算的～