PRO Kerr, r, theta, phi, Metricg, Affineg, A=A, DetG=DetG, M=M

;  Compute metric coeffs and connection coeffs, in Boyer-Lindquist
;  coordinates for a Kerr black hole.

default,a,0.998
default,M,1.0

;  MathTensor expressions are with a=1,2,3,4, with t=x[4]
Metricg = fltarr(5,5)
AffineG = fltarr(5,5,5)

;  Scratch results
a2 = a^2
r2 = r^2
r4 = r2^2
sth = Sin(theta)
sth2 = sth^2
sth4 = sth2^2
cthsth = Cos(theta)*sth
cthsth3 = cthsth *sth2
cthsth5 = cthsth3*sth2
a2sth2 = a2*sth2
r4sth2 = r4*sth2
sth6 = sth2*sth4
cotth = Cos(theta)/(sth > 1e-10)

Detg = -(sth2*(-a2 - r^2 + a2sth2)^2)

; Metric coeffs
Metricg(4, 4) = (a2 - 2*M*r + r2 - a2sth2)/(-a2 - r2 + a2sth2) 
Metricg(4, 3) = (2*a*M*r*sth2)/(-a2 - r2 + a2sth2)
Metricg(3, 3) = (-(a2*a2sth2) - 2*a2*r2*sth2 - r4sth2 + $
  a^4*sth4 - 2*a2*M*r*sth4 + $
  a2*r2*sth4)/(-a2 - r2 + a2sth2)
Metricg(2, 2) = a2 + r2 - a2sth2
Metricg(1, 1) = (a2 + r2 - a2sth2)/(a2 - 2*M*r + r2)

;  Connection coeffs
AffineG(1, 4, 4) = (a^4*M - 2*a2*M^2*r + 2*M^2*r^3 - M*r^4 - a^4*M*sth2 + $
  2*a2*M^2*r*sth2 - a2*M*r2*sth2)/(-a2 - r2 + a2sth2)^3 
AffineG(1, 4, 3) = (-(a^5*M*sth2) + 2*a^3*M^2*r*sth2 - $
  2*a*M^2*r^3*sth2 + a*M*r4sth2 + a^5*M*sth4 - $
  2*a^3*M^2*r*sth4 + a^3*M*r2*sth4)/(-a2 - r2 + a2sth2)^3 
AffineG(1, 3, 3) = (a^6*r*sth2 - 2*a^4*M*r2*sth2 + $
  3*a^4*r^3*sth2 - 4*a2*M*r4sth2 + 3*a2*r^5*sth2 - $
  2*M*r^6*sth2 + r^7*sth2 + a^6*M*sth4 - $
  2*a^6*r*sth4 - 2*a^4*M^2*r*sth4 + 4*a^4*M*r2*sth4 - $
  4*a^4*r^3*sth4 + 2*a2*M^2*r^3*sth4 + 3*a2*M*r^4*sth4 - $
  2*a2*r^5*sth4 - a^6*M*sth6 + a^6*r*sth6 + $
  2*a^4*M^2*r*sth6 - 3*a^4*M*r2*sth6 + $
  a^4*r^3*sth6)/(-a2 - r2 + a2sth2)^3 
AffineG(1, 2, 1) = (a2*cthsth)/(-a2 - r2 + a2sth2)
AffineG(1, 1, 1) = (-(a2*M) + M*r2 + a2*M*sth2 - $
  a2*r*sth2)/((a2 - 2*M*r + r2)*(-a2 - r2 + a2sth2))
AffineG(2, 4, 4) = (2*a2*M*r*cthsth)/(-a2 - r2 + a2sth2)^3
AffineG(2, 4, 3) = (-2*(a^3*M*r*cthsth + $
  a*M*r^3*cthsth))/(-a2 - r2 + a2sth2)^3
AffineG(2, 3, 3) = (a^6*cthsth + 3*a^4*r2*cthsth + $
  3*a2*r^4*cthsth + r^6*cthsth - $
  2*a^6*cthsth3 + 4*a^4*M*r*cthsth3 - $
  4*a^4*r2*cthsth3 + 4*a2*M*r^3*cthsth3 - $
  2*a2*r^4*cthsth3 + a^6*cthsth5 - $
  2*a^4*M*r*cthsth5 + $
  a^4*r2*cthsth5)/(-a2 - r2 + a2sth2)^3
AffineG(2, 2, 2) = (a2*cthsth)/(-a2 - r2 + a2sth2)
AffineG(2, 2, 1) = r/(a2 + r2 - a2sth2) 
AffineG(2, 1, 1) = -((a2*cthsth)/((a2 - 2*M*r + r2)*(-a2 - r2 + a2sth2)))
AffineG(3, 4, 2) = (-2*a*M*r*cotth)/(-a2 - r2 + a2sth2)^2
AffineG(3, 4, 1) = (-(a^3*M) + a*M*r2 + $
  a^3*M*sth2)/((a2 - 2*M*r + r2)*(-a2 - r2 + a2sth2)^2) 
AffineG(3, 3, 2) = (a^4*cotth + 2*a2*r2*cotth + r^4*cotth - $
  2*a^4*cthsth + 2*a2*M*r*cthsth - $
  2*a2*r2*cthsth + $
  a^4*cthsth3)/(-a2 - r2 + a2sth2)^2 
AffineG(3, 3, 1) = (a^4*r - 2*a2*M*r2 + 2*a2*r^3 - 2*M*r^4 + r^5 + $
  a^4*M*sth2 - 2*a^4*r*sth2 + a2*M*r2*sth2 - $
  2*a2*r^3*sth2 - a^4*M*sth4 + $
  a^4*r*sth4)/((a2 - 2*M*r + r2)*(-a2 - r2 + a2sth2)^2)
AffineG(4, 4, 2) = (-2*a2*M*r*cthsth)/(-a2 - r2 + a2sth2)^2 
AffineG(4, 4, 1) = (-(a^4*M) + M*r^4 + a^4*M*sth2 + $
  a2*M*r2*sth2)/((a2 - 2*M*r + r2)*(-a2 - r2 + a2sth2)^2) 
AffineG(4, 3, 2) = (2*a^3*M*r*cthsth3)/(-a2 - r2 + a2sth2)^2
AffineG(4, 3, 1) = (a^5*M*sth2 - 2*a^3*M*r2*sth2 - $
  3*a*M*r4sth2 - a^5*M*sth4 + $
  a^3*M*r2*sth4)/((a2 - 2*M*r + r2)*(-a2 - r2 + a2sth2)^2)

; Fill out with symmetries
for j=1,4 do begin
  for k=j+1,4 do begin
      Metricg(j,k)=Metricg(k,j)
      AffineG(*,j,k) = AffineG(*,k,j)
  end
end

; Shift from 1,2,3,4 to 0,1,2,3
Metricg = Metricg(1:4,1:4)
Metricg = shift(Metricg,1,0)
Metricg = shift(Metricg,0,1)

Affineg = Affineg(1:4,1:4,1:4)
Affineg = shift(Affineg,1,0,0)
Affineg = shift(Affineg,0,1,0)
Affineg = shift(Affineg,0,0,1)

END