0

I'm trying to calculate the Kretschmann scalar in mathematica, it is given by:

$c = R^{abcd} R_{abcd}$

Where $R^{abcd}$ is the Riemann tensor. I'm following this MSE post so I modified it to get the corresponding $R^{abcd}$ in the following manner:

RiemannTensorarriba[g_, xx_] := Block[{n, Rie, ig, res}, n = 4; 
Rie = RiemannTensor[g, xx]; 
ig = InverseMetric[g];                       
res = Table[
Sum[ig[[i, a]]*ig[[j, b]]*ig[[k, c]]*Rie[l, a, b, c], {a, 1, 
  n}, {b, 1, n}, {c, 1, n}], {i, 1, n}, {j, 1, n}, {k, 1, n}, {l, 
 1, n}];
Simplify[res]]

And I'm raising the first index in the Riemann Tensor of the MSE:

Riemannabajo[g_, xx_] := Block[{n, rie, res}, n = 4;
rie = RiemannTensor[g, xx];
res = Table[
Sum[g[[i, j]]*rie[[j, a, b, c]], {j, 1, n}], {i, 1, n}, {a, 1, 
 n}, {b, 1, n}, {c, 1, n}];
Simplify[res]]

Now if I have the both tensors that I need, so to contract the indices I'm trying both:

Riemannabajo[g,xx].RiemannTensorArriba[g,xx]

And programing it via a summation like:

Kreztchman[g_, xx_] := Block[{n, Rie, Riea, res}, n = 4;
Rie = Riemannabajo[g, xx];
Riea = RiemannTensorarriba[g, xx];
res = Sum[
Riea[[i, k, l, m]]*Rie[[i, k, l, m]], {i, 1, n}, {k, 1, n}, {l, 1,
  n}, {m, 1, n}];
Simplify[res]]

But both answers don't give me a scalar, instead I get some kind of tensor. Where I'm missing the point?

Edit: In order to be consistent, Im using the solution of the MSE as proposed by @Arte, ie:

InverseMetric[g_]:=Simplify[Inverse[g]]

ChristoffelSymbol[g_,xx_]:=Block[{n,ig,res},n=4;ig=InverseMetric[g];
res=Table[(1/2)*Sum[ig[[i,s]]*(-
D[g[[j,k]],xx[[s]]]+D[g[[j,s]],xx[[k]]]+D[g[[s,k]],xx[[j]]]),{s,1,n}],{i,1,n},{j,1,n},{k,1,n}];
Simplify[res]]
RiemannTensor[g_,xx_]:=Block[{n,Chr,res},n=4;Chr=ChristoffelSymbol[g,xx];
res=Table[D[Chr[[i,k,m]],xx[[l]]]-D[Chr[[i,k,l]],xx[[m]]]+Sum[Chr[[i,s,l]]*Chr[[s,k,m]],{s,1,n}]-Sum[Chr[[i,s,m]]*Chr[[s,k,l]],{s,1,n}],{i,1,n},{k,1,n},{l,1,n},{m,1,n}];
Simplify[res]]

RicciTensor[g_,xx_]:=Block[{Rie,res,n},n=4;Rie=RiemannTensor[g,xx];
    res=Table[Sum[Rie[[s,i,s,j]],{s,1,n}],{i,1,n},{j,1,n}];
    Simplify[res]]
RicciScalar[g_,xx_]:=Block[{Ricc,ig,res,n},n=4;Ricc=RicciTensor[g,xx];ig=InverseMetric[g];
res=Sum[ig[[s,i]] Ricc[[s,i]],{s,1,n},{i,1,n}];
Simplify[res]]

RiemannTensorarriba[g_,xx_]:=Block[{n,Rie,ig,res},n=4;
Rie=RiemannTensor[g,xx];
ig=InverseMetric[g];
res=Table[Sum[ig[[i,a]]*ig[[j,b]]*ig[[k,c]]*Rie[l,a,b,c],{a,1,n},{b,1,n},{c,1,n}],{i,1,n},{j,1,n},{k,1,n},{l,1,n}];
Simplify[res]]

Riemannabajo[g_,xx_]:=Block[{n,rie,res},n=4;
rie=RiemannTensor[g,xx];
res=Table[Sum[g[[i,j]]*rie[[j,a,b,c]],{j,1,n}],{i,1,n},{a,1,n},{b,1,n},{c,1,n}];
Simplify[res]]

And I'm trying to calculete the Kretschmann invariant to the following metric:

xx = {t, x, \[Theta], \[Phi]};
g = {{-(1 - x^2/a^2), 0, 0, 0}, {0, 1/(1 - x^2/a^2), 0, 0}, {0, 0, 
x^2, 0}, {0, 0, 0, x^2 Sin[\[Theta]]^2}};
Jasimud
  • 119
  • 3
  • Even after copying the definitions from the linked question, your code does not execute successfully, instead generating a slew of errors. Please make sure that you a) include all relevant definitions and parameter values, and b) make sure that your code will execute correctly in a fresh Mathematica session. We won't be able to help otherwise. – MarcoB Jun 08 '16 at 00:01
  • I tought it to be complete. Anyway there is the edit with all the code in order to work – Jasimud Jun 08 '16 at 00:18
  • Jasimud, your RiemannTensors are rank-4 tensors. According to the docs, Dot contracts the last index in its first arguments with the first index in its second argument; for your tensors, this will result in a tensor of rank 6, and not a scalar, so this can't possibly work. Perhaps you want to look at Inner or TensorContract. – MarcoB Jun 08 '16 at 00:38
  • Since it's a rank-4 tensor in order tu apply Inner, or TensorContract I would need to put the list of indices I need to contract, for example only for the first pair of indices it would be {1,1,a,b,c},{2,2,a,b,c}{3,3,a,b,c}{4,4,a,b,c}, and then I need to expand all other 3 indeces, there are lot's of contractions, there must be a simpler way than just counting them by hand – Jasimud Jun 08 '16 at 01:21
  • 1
    The capitalization of variable names is inconsistent (e.g., RiemannTensorarriba vs RiemannTensorArriba). Pleass make sure that you don't have stray definitions by restarting the kernel. It might be helpful to English speakers to replace arriba and abajo by Up and Down. Your approach with Sum[...,{i,1,n},{j,1,n},{k,1,n},{l,1,n}] seems reasonable, but it seems that Riemannabajo[[i,j,k,l]] should be Riemannabajo[g,xx][[i,j,k,l]] if I read the code correctly. The order of indices in RiemannTensorarriba seems wrong (should be l, i, j, k). – Bruno Le Floch Jun 08 '16 at 13:40
  • @BrunoLeFloch I edited the post with the full code I was using for the scalar, it still gives me some higher rank tensor, and I don't seem to grasp why or why the order of the contraction for RiemannTensorarriba should be l i j k – Jasimud Jun 08 '16 at 14:45
  • 1
    Rie[l,a,b,c] should have double brackets. Perhaps that mistake occurs elsewhere: you can check the dimension of various parts by Dimensions[Riemannabajo[g,xx]] etc. As for order of indices, if I read Sum[ig[[i, a]]*ig[[j, b]]*ig[[k, c]]*Rie[l, a, b, c], {a, 1, n}, {b, 1, n}, {c, 1, n}] correctly it computes $g^{ia}g^{jb}g^{kc}R^l{}{abc} = R^{lijk}$ not $R^{ijkl}$. But maybe Rie[[l, a, b, c]] means $R{abc}{}^l$. – Bruno Le Floch Jun 08 '16 at 17:49

1 Answers1

2

There might be some extra bits feel free to remove

   ClearAll["Global`*"];

n = 4;
coord = {t, x, y, z};

(*For raising/lowering latin indices like a, b, ...*)
η = {{-1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, 1, 0}, {0, 0, 0, 1}};
η // MatrixForm;
inveta = Inverse[η];
inveta // MatrixForm;

(*"e-Metric";*)
e = {{\[ScriptCapitalN][t], 0, 0, 0}, {0, \[ScriptA][t], 0, 0}, {0, 
    0, \[ScriptA][t], 0}, {0, 0, 0, \[ScriptA][t]}};
\[ScriptCapitalN][t] = -1;
e // MatrixForm;
dete = Det[e];
(*Inverse e-Metric*)
inve = Inverse[e];
inve // MatrixForm;
detinve = Det[inve];

g := g = Simplify[Table[
    ParallelSum[ η[[a, b]]*e[[a, μ]]*e[[b, ν]]
     , {a, 1, n}, {b, 1, n}]
    , {μ, 1, n}, {ν, 1, n}]]
g // MatrixForm;
(*g is used to LOWER indices for Greek indices μ, ν*)
invg = Inverse[g];
invg // MatrixForm;
(*invg is used to RAISE indices for Greek indices μ, ν*)

(* In the form Γ^{x}_{xx}*)
affine := affine = Simplify[Table[(1/2)*Sum[(invg[[i, s]])*
        (D[g[[s, j]], coord[[k]] ] +
          D[g[[s, k]], coord[[j]] ] - D[g[[j, k]], coord[[s]] ]), {s, 
        1, n}],
     {i, 1, n}, {j, 1, n}, {k, 1, n}] ];
listaffine := 
  Table[If[UnsameQ[affine[[i, j, k]], 
     0], {ToString[Γ[i, j, k]], 
     affine[[i, j, k]]}] , {i, 1, n}, {j, 1, n}, {k, 1, j}];
TableForm[Partition[DeleteCases[Flatten[listaffine], Null], 2], 
  TableSpacing -> {2, 2}];

(* In the form R^{x}_{xxx} *)
riemann := riemann = Simplify[Table[
     D[affine[[i, j, l]], coord[[k]] ] - 
      D[affine[[i, j, k]], coord[[l]] ] +
      Sum[
       affine[[s, j, l]] affine[[i, k, s]] - 
        affine[[s, j, k]] affine[[i, l, s]],
       {s, 1, n}],
     {i, 1, n}, {j, 1, n}, {k, 1, n}, {l, 1, n}] ];
listriemann := 
 Table[If[UnsameQ[riemann[[i, j, k, l]], 0], {ToString[R[i, j, k, l]],
     riemann[[i, j, k, l]]}] , {i, 1, n}, {j, 1, n}, {k, 1, n}, {l, 1,
    k - 1}]
TableForm[Partition[DeleteCases[Flatten[listriemann], Null], 2], 
  TableSpacing -> {2, 2}];

(* In the form R_{xxxx} *)
riemann1 := riemann1 = Simplify[Table[
    Sum[
     g[[μ, μ1]]*riemann[[μ1, ν, ρ, σ]]
     , {μ1, 1, n}]
    , {μ, 1, n}, {ν, 1, n}, {ρ, 1, n}, {σ, 1, n}]]
listriemann1 := 
  Table[If[UnsameQ[riemann1[[μ, ν, ρ, σ]], 
     0], {ToString[R1[μ, ν, ρ, σ]], 
     riemann1[[μ, ν, ρ, σ]]}] , {μ, 1, 
    n}, {ν, 1, n}, {ρ, 1, n}, {σ, 1, ρ - 1}];
TableForm[Partition[DeleteCases[Flatten[listriemann1], Null], 2], 
  TableSpacing -> {2, 2}];

(* In the form R^{xxxx} *)
riemann2 := riemann2 = Simplify[Table[
    Sum[
     Sum[
      Sum[
       invg[[ν1, ν]]*invg[[ρ1, ρ]]*
        invg[[σ1, σ]]*
        riemann[[μ, ν1, ρ1, σ1]]
       , {ν1, 1, n}]
      , {ρ1, 1, n}]
     , {σ1, 1, n}]
    , {μ, 1, n}, {ν, 1, n}, {ρ, 1, n}, {σ, 1, n}]]
listriemann2 := 
  Table[If[UnsameQ[riemann2[[μ, ν, ρ, σ]], 
     0], {ToString[R2[μ, ν, ρ, σ]], 
     riemann2[[μ, ν, ρ, σ]]}] , {μ, 1, 
    n}, {ν, 1, n}, {ρ, 1, n}, {σ, 1, ρ - 1}];
TableForm[Partition[DeleteCases[Flatten[listriemann2], Null], 2], 
  TableSpacing -> {2, 2}];

KretschmannScalar = Simplify[
  Sum[riemann1[[a, b, c, d]]*riemann2[[a, b, c, d]], {a, 1, n}, {b, 1,
     n}, {c, 1, n}, {d, 1, n}]]

With that input metric I got this answer:

$\frac{12 \left(\mathit{a}(t)^2 \mathit{a}''(t)^2+\mathit{a}'(t)^4\right)}{\mathit{a}(t)^4}$

Henrik Schumacher
  • 106,770
  • 7
  • 179
  • 309
Mark Pace
  • 316
  • 1
  • 12