clearAll; 
"Kerr-Newman metrika"
Coords = {t, r, θ, ϕ}; 
ρ = r^2 + a^2Cos[θ]^2; 
Δ = r^2 - 2M r + a^2 + (Q^2); 
tt = (2M r - (Q^2))/ρ - 1; 
rr = ρ/Δ; 
θθ = ρ; 
ϕϕ = (r^2 + a^2 + ((2M r -  (Q^2))a^2Sin[θ]^2)/ρ)Sin[θ]^2; 
tϕ = -((2 M r - (Q^2)) a Sin[θ]^2)/ρ; 

Metric = {{tt, 0, 0, tϕ}, {0, rr, 0, 0}, {0, 0, θθ, 0}, {tϕ, 0, 0, ϕϕ}}; 
d = Length[Coords]; 
InverseMetric=Simplify[Inverse[Metric]] ;
Nonzero[x_]:=Not[StringMatchQ[ToString[x],"0"]]
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Metric[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[\(g\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],Metric[[i]][[j]]]]]]]
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[InverseMetric[[i]][[j]]],Print[StringForm["\!\(\*SuperscriptBox[\(g\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],InverseMetric[[i]][[j]]]]]]]
"Leva stran Einstainovih enačb"
Christoffel=Simplify[1/2*Table[Sum[InverseMetric[[i]][[l]]*(D[Metric[[l]][[j]],Coords[[k]]]+D[Metric[[l]][[k]],Coords[[j]]]-D[Metric[[j]][[k]],Coords[[l]]]),{l,d}],{i,d},{j,d},{k,d}]] ;
For[i=1,i≤d,i++,For[j=1,j≤d,j++,For[k=j,k≤d,k++,If[Nonzero[Christoffel[[i]][[j]][[k]]],Print[StringForm["\!\(\*SubscriptBox[SuperscriptBox[\(Γ\), \(`1`\)], \(`2`\\ `3`\)]\)=`4`",Coords[[i]],Coords[[j]],Coords[[k]],Christoffel[[i]][[j]][[k]]]]]]]]
Riemann=Simplify[Table[D[Christoffel[[l]][[i]][[k]],Coords[[j]]]-
D[Christoffel[[l]][[j]][[k]],Coords[[i]]]+
Sum[ 
Christoffel[[m]][[i]][[k]]*Christoffel[[l]][[m]][[j]]-Christoffel[[m]][[j]][[k]]*Christoffel[[l]][[m]][[i]],{m,1,d}],{i,d},{j,d},{k,d},{l,d}]] ;
For[i=1,i≤d,i++,For[j=1,j≤d,j++,For[k=1,k≤d,k++,For[l=1,l≤d,l++,If[Nonzero[Riemann[[i]][[j]][[k]][[l]]],Print[StringForm["\!\(\*SuperscriptBox[SubscriptBox[\(R\), \(`1`\\ `2`\\ `3`\)], \(`4`\)]\)=`5`",Coords[[i]],Coords[[j]],Coords[[k]],Coords[[l]],Riemann[[i]][[j]][[k]][[l]]]]]]]]]
RicciTensor=FullSimplify[Table[Sum[Riemann[[i]][[k]][[j]][[k]],{k,1,d}],{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[RicciTensor[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[\(R\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],RicciTensor[[i]][[j]]]]]]]
RicciScalar=FullSimplify[Sum[InverseMetric[[i]][[j]]RicciTensor[[i]][[j]],{i,1,d},{j,1,d}]] ;
Print[StringForm["R=`1`",RicciScalar]]
Einstein=FullSimplify[Table[RicciTensor[[i]][[j]]-1/2Metric[[i]][[j]]*RicciScalar,{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Einstein[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[\(G\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],Einstein[[i]][[j]]]]]]]

"Leva stran Einstainovih enačb"
Potential={At,0,0,Aϕ};
At=-Q r/ρ ;
Aϕ=Q r a Sin[θ]^2/ρ ;
Electromag=Simplify[Table[-D[Potential[[j]],Coords[[i]]]+D[Potential[[i]],Coords[[j]]],{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Electromag[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[\(F\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],Electromag[[i]][[j]]]]]]]
Electromag1=FullSimplify[Table[Sum[InverseMetric[[i]][[k]]Electromag[[j]][[k]],{k,1,d}],{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Electromag1[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[SuperscriptBox[\(F\), \(\(`1`\)\(\\ \)\)], \(`2`\)]\)=`3`",Coords[[i]],Coords[[j]],Electromag1[[i]][[j]]]]]]]
Electromag2=FullSimplify[Table[Sum[InverseMetric[[i]][[k]]Electromag1[[j]][[k]],{k,1,d}],{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Electromag2[[i]][[j]]],Print[StringForm["\!\(\*SuperscriptBox[\(F\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],Electromag2[[i]][[j]]]]]]]
Electromag3=FullSimplify[Table[Sum[Electromag[[i]][[k]]Electromag1[[k]][[j]],{k,1,d}],{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Electromag3[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[\(FP\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],Electromag3[[i]][[j]]]]]]]
ScalarF=FullSimplify[Sum[Electromag[[i]][[j]]Electromag2[[i]][[j]],{i,1,d},{j,1,d}]] ;
Print[StringForm["F=`1`",ScalarF]]
Stress=FullSimplify[Table[(Electromag3[[i]][[j]]-1/4Metric[[i]][[j]]*ScalarF)/(4π),{i,d},{j,d}]] ;
For[i=1,i≤d,i++,For[j=i,j≤d,j++,If[Nonzero[Stress[[i]][[j]]],Print[StringForm["\!\(\*SubscriptBox[\(T\), \(`1`\\ `2`\)]\)=`3`",Coords[[i]],Coords[[j]],Stress[[i]][[j]]]]]]]

"Razvoj metrike"
Simplify[Series[(2M r - (Q^2))/ρ - 1,{r,Infinity,3}]]
Simplify[Series[ρ/Δ, {r,Infinity,3}]]
Simplify[Series[ρ, {r,Infinity,3}]]
Simplify[Series[(r^2 + a^2 + ((2M r -  (Q^2))a^2Sin[θ]^2)/ρ)Sin[θ]^2, {r,Infinity,3}]]
Simplify[Series[ -((2 M r - (Q^2)) a Sin[θ]^2)/ρ, {r,Infinity,3}]]

"Razvoj elektromagnetnega polja"
Simplify[Series[((Q((r^2)-(a^2)((Cos[θ])^2)))/((((r^2)+(a^2) ((Cos[θ])^2)))^2)),{r,Infinity,6}]]
Simplify[Series[-((2Q r (a^2)Sin[θ]Cos[θ])/((((r^2)+(a^2) ((Cos[θ])^2)))^2)),{r,Infinity,6}]]
Simplify[Series[-((Q a ((Sin[θ])^2)((r^2)-(a^2)((Cos[θ])^2)))/((((r^2)+(a^2) ((Cos[θ])^2)))^2)),{r,Infinity,6}]]
Simplify[Series[-((2Q r a Sin[θ]Cos[θ])/((((r^2)+(a^2) ((Cos[θ])^2)))^2))((r^2)+(a^2)),{r,Infinity,6}]]