subroutine ricci_gamma(ex, X, Y, Z, & chi, & dxx , gxy , gxz , dyy , gyz , dzz,& Gamx , Gamy , Gamz , & Gamxxx,Gamxxy,Gamxxz,Gamxyy,Gamxyz,Gamxzz,& Gamyxx,Gamyxy,Gamyxz,Gamyyy,Gamyyz,Gamyzz,& Gamzxx,Gamzxy,Gamzxz,Gamzyy,Gamzyz,Gamzzz,& Rxx,Rxy,Rxz,Ryy,Ryz,Rzz,& Symmetry) implicit none !~~~~~~> Input parameters: integer,intent(in ):: ex(1:3), Symmetry real*8, intent(in ):: X(1:ex(1)),Y(1:ex(2)),Z(1:ex(3)) real*8, dimension(ex(1),ex(2),ex(3)),intent(in ) :: chi real*8, dimension(ex(1),ex(2),ex(3)),intent(in ) :: dxx,gxy,gxz,dyy,gyz,dzz real*8, dimension(ex(1),ex(2),ex(3)),intent(in ) :: Gamx,Gamy,Gamz ! when out, physical second kind of connection real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamxxx, Gamxxy, Gamxxz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamxyy, Gamxyz, Gamxzz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamyxx, Gamyxy, Gamyxz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamyyy, Gamyyz, Gamyzz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamzxx, Gamzxy, Gamzxz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamzyy, Gamzyz, Gamzzz ! when out, physical Ricci tensor real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Rxx,Rxy,Rxz,Ryy,Ryz,Rzz !~~~~~~> Other variables: real*8, dimension(ex(1),ex(2),ex(3)) :: gxx,gyy,gzz real*8, dimension(ex(1),ex(2),ex(3)) :: chix,chiy,chiz real*8, dimension(ex(1),ex(2),ex(3)) :: gxxx,gxyx,gxzx,gyyx,gyzx,gzzx real*8, dimension(ex(1),ex(2),ex(3)) :: gxxy,gxyy,gxzy,gyyy,gyzy,gzzy real*8, dimension(ex(1),ex(2),ex(3)) :: gxxz,gxyz,gxzz,gyyz,gyzz,gzzz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamxx,Gamxy,Gamxz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamyx,Gamyy,Gamyz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamzx,Gamzy,Gamzz real*8, dimension(ex(1),ex(2),ex(3)) :: f,fxx,fxy,fxz,fyy,fyz,fzz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamxa,Gamya,Gamza,chin1 real*8, dimension(ex(1),ex(2),ex(3)) :: gupxx,gupxy,gupxz real*8, dimension(ex(1),ex(2),ex(3)) :: gupyy,gupyz,gupzz real*8 :: dX, dY, dZ real*8, parameter :: ONE = 1.D0, TWO = 2.D0, FOUR = 4.D0 real*8, parameter :: HALF = 0.5D0, F2o3 = 2.d0/3.d0, F3o2 = 1.5d0 real*8, parameter :: SYM = 1.D0, ANTI= - 1.D0 dX = X(2) - X(1) dY = Y(2) - Y(1) dZ = Z(2) - Z(1) chin1 = chi + ONE gxx = dxx + ONE gyy = dyy + ONE gzz = dzz + ONE call fderivs(ex,chi,chix,chiy,chiz,X,Y,Z,SYM,SYM,SYM,Symmetry,0) call fderivs(ex,dxx,gxxx,gxxy,gxxz,X,Y,Z,SYM ,SYM ,SYM ,Symmetry,0) call fderivs(ex,gxy,gxyx,gxyy,gxyz,X,Y,Z,ANTI,ANTI,SYM ,Symmetry,0) call fderivs(ex,gxz,gxzx,gxzy,gxzz,X,Y,Z,ANTI,SYM ,ANTI,Symmetry,0) call fderivs(ex,dyy,gyyx,gyyy,gyyz,X,Y,Z,SYM ,SYM ,SYM ,Symmetry,0) call fderivs(ex,gyz,gyzx,gyzy,gyzz,X,Y,Z,SYM ,ANTI,ANTI,Symmetry,0) call fderivs(ex,dzz,gzzx,gzzy,gzzz,X,Y,Z,SYM ,SYM ,SYM ,Symmetry,0) ! invert tilted metric gupzz = gxx * gyy * gzz + gxy * gyz * gxz + gxz * gxy * gyz - & gxz * gyy * gxz - gxy * gxy * gzz - gxx * gyz * gyz gupxx = ( gyy * gzz - gyz * gyz ) / gupzz gupxy = - ( gxy * gzz - gyz * gxz ) / gupzz gupxz = ( gxy * gyz - gyy * gxz ) / gupzz gupyy = ( gxx * gzz - gxz * gxz ) / gupzz gupyz = - ( gxx * gyz - gxy * gxz ) / gupzz gupzz = ( gxx * gyy - gxy * gxy ) / gupzz ! second kind of connection Gamxxx =HALF*( gupxx*gxxx + gupxy*(TWO*gxyx - gxxy ) + gupxz*(TWO*gxzx - gxxz )) Gamyxx =HALF*( gupxy*gxxx + gupyy*(TWO*gxyx - gxxy ) + gupyz*(TWO*gxzx - gxxz )) Gamzxx =HALF*( gupxz*gxxx + gupyz*(TWO*gxyx - gxxy ) + gupzz*(TWO*gxzx - gxxz )) Gamxyy =HALF*( gupxx*(TWO*gxyy - gyyx ) + gupxy*gyyy + gupxz*(TWO*gyzy - gyyz )) Gamyyy =HALF*( gupxy*(TWO*gxyy - gyyx ) + gupyy*gyyy + gupyz*(TWO*gyzy - gyyz )) Gamzyy =HALF*( gupxz*(TWO*gxyy - gyyx ) + gupyz*gyyy + gupzz*(TWO*gyzy - gyyz )) Gamxzz =HALF*( gupxx*(TWO*gxzz - gzzx ) + gupxy*(TWO*gyzz - gzzy ) + gupxz*gzzz) Gamyzz =HALF*( gupxy*(TWO*gxzz - gzzx ) + gupyy*(TWO*gyzz - gzzy ) + gupyz*gzzz) Gamzzz =HALF*( gupxz*(TWO*gxzz - gzzx ) + gupyz*(TWO*gyzz - gzzy ) + gupzz*gzzz) Gamxxy =HALF*( gupxx*gxxy + gupxy*gyyx + gupxz*( gxzy + gyzx - gxyz ) ) Gamyxy =HALF*( gupxy*gxxy + gupyy*gyyx + gupyz*( gxzy + gyzx - gxyz ) ) Gamzxy =HALF*( gupxz*gxxy + gupyz*gyyx + gupzz*( gxzy + gyzx - gxyz ) ) Gamxxz =HALF*( gupxx*gxxz + gupxy*( gxyz + gyzx - gxzy ) + gupxz*gzzx ) Gamyxz =HALF*( gupxy*gxxz + gupyy*( gxyz + gyzx - gxzy ) + gupyz*gzzx ) Gamzxz =HALF*( gupxz*gxxz + gupyz*( gxyz + gyzx - gxzy ) + gupzz*gzzx ) Gamxyz =HALF*( gupxx*( gxyz + gxzy - gyzx ) + gupxy*gyyz + gupxz*gzzy ) Gamyyz =HALF*( gupxy*( gxyz + gxzy - gyzx ) + gupyy*gyyz + gupyz*gzzy ) Gamzyz =HALF*( gupxz*( gxyz + gxzy - gyzx ) + gupyz*gyyz + gupzz*gzzy ) Gamxa = gupxx * Gamxxx + gupyy * Gamxyy + gupzz * Gamxzz + & TWO*( gupxy * Gamxxy + gupxz * Gamxxz + gupyz * Gamxyz ) Gamya = gupxx * Gamyxx + gupyy * Gamyyy + gupzz * Gamyzz + & TWO*( gupxy * Gamyxy + gupxz * Gamyxz + gupyz * Gamyyz ) Gamza = gupxx * Gamzxx + gupyy * Gamzyy + gupzz * Gamzzz + & TWO*( gupxy * Gamzxy + gupxz * Gamzxz + gupyz * Gamzyz ) call fderivs(ex,Gamx,Gamxx,Gamxy,Gamxz,X,Y,Z,ANTI,SYM ,SYM ,Symmetry,0) call fderivs(ex,Gamy,Gamyx,Gamyy,Gamyz,X,Y,Z,SYM ,ANTI,SYM ,Symmetry,0) call fderivs(ex,Gamz,Gamzx,Gamzy,Gamzz,X,Y,Z,SYM ,SYM ,ANTI,Symmetry,0) !first kind of connection stored in gij,k gxxx = gxx * Gamxxx + gxy * Gamyxx + gxz * Gamzxx gxyx = gxx * Gamxxy + gxy * Gamyxy + gxz * Gamzxy gxzx = gxx * Gamxxz + gxy * Gamyxz + gxz * Gamzxz gyyx = gxx * Gamxyy + gxy * Gamyyy + gxz * Gamzyy gyzx = gxx * Gamxyz + gxy * Gamyyz + gxz * Gamzyz gzzx = gxx * Gamxzz + gxy * Gamyzz + gxz * Gamzzz gxxy = gxy * Gamxxx + gyy * Gamyxx + gyz * Gamzxx gxyy = gxy * Gamxxy + gyy * Gamyxy + gyz * Gamzxy gxzy = gxy * Gamxxz + gyy * Gamyxz + gyz * Gamzxz gyyy = gxy * Gamxyy + gyy * Gamyyy + gyz * Gamzyy gyzy = gxy * Gamxyz + gyy * Gamyyz + gyz * Gamzyz gzzy = gxy * Gamxzz + gyy * Gamyzz + gyz * Gamzzz gxxz = gxz * Gamxxx + gyz * Gamyxx + gzz * Gamzxx gxyz = gxz * Gamxxy + gyz * Gamyxy + gzz * Gamzxy gxzz = gxz * Gamxxz + gyz * Gamyxz + gzz * Gamzxz gyyz = gxz * Gamxyy + gyz * Gamyyy + gzz * Gamzyy gyzz = gxz * Gamxyz + gyz * Gamyyz + gzz * Gamzyz gzzz = gxz * Gamxzz + gyz * Gamyzz + gzz * Gamzzz !compute Ricci tensor for tilted metric call fdderivs(ex,dxx,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,SYM ,SYM ,SYM ,Symmetry,0) Rxx = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs(ex,dyy,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,SYM ,SYM ,SYM ,Symmetry,0) Ryy = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs(ex,dzz,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,SYM ,SYM ,SYM ,Symmetry,0) Rzz = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs(ex,gxy,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,ANTI,ANTI,SYM ,Symmetry,0) Rxy = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs(ex,gxz,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,ANTI,SYM ,ANTI,Symmetry,0) Rxz = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs(ex,gyz,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,SYM,ANTI,ANTI ,Symmetry,0) Ryz = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO Rxx = - HALF * Rxx + & gxx * Gamxx+ gxy * Gamyx + gxz * Gamzx + & Gamxa * gxxx + Gamya * gxyx + Gamza * gxzx + & gupxx *( & TWO*(Gamxxx * gxxx + Gamyxx * gxyx + Gamzxx * gxzx) + & Gamxxx * gxxx + Gamyxx * gxxy + Gamzxx * gxxz )+ & gupxy *( & TWO*(Gamxxx * gxyx + Gamyxx * gyyx + Gamzxx * gyzx + & Gamxxy * gxxx + Gamyxy * gxyx + Gamzxy * gxzx) + & Gamxxy * gxxx + Gamyxy * gxxy + Gamzxy * gxxz + & Gamxxx * gxyx + Gamyxx * gxyy + Gamzxx * gxyz )+ & gupxz *( & TWO*(Gamxxx * gxzx + Gamyxx * gyzx + Gamzxx * gzzx + & Gamxxz * gxxx + Gamyxz * gxyx + Gamzxz * gxzx) + & Gamxxz * gxxx + Gamyxz * gxxy + Gamzxz * gxxz + & Gamxxx * gxzx + Gamyxx * gxzy + Gamzxx * gxzz )+ & gupyy *( & TWO*(Gamxxy * gxyx + Gamyxy * gyyx + Gamzxy * gyzx) + & Gamxxy * gxyx + Gamyxy * gxyy + Gamzxy * gxyz )+ & gupyz *( & TWO*(Gamxxy * gxzx + Gamyxy * gyzx + Gamzxy * gzzx + & Gamxxz * gxyx + Gamyxz * gyyx + Gamzxz * gyzx) + & Gamxxz * gxyx + Gamyxz * gxyy + Gamzxz * gxyz + & Gamxxy * gxzx + Gamyxy * gxzy + Gamzxy * gxzz )+ & gupzz *( & TWO*(Gamxxz * gxzx + Gamyxz * gyzx + Gamzxz * gzzx) + & Gamxxz * gxzx + Gamyxz * gxzy + Gamzxz * gxzz ) Ryy = - HALF * Ryy + & gxy * Gamxy+ gyy * Gamyy + gyz * Gamzy + & Gamxa * gxyy + Gamya * gyyy + Gamza * gyzy + & gupxx *( & TWO*(Gamxxy * gxxy + Gamyxy * gxyy + Gamzxy * gxzy) + & Gamxxy * gxyx + Gamyxy * gxyy + Gamzxy * gxyz )+ & gupxy *( & TWO*(Gamxxy * gxyy + Gamyxy * gyyy + Gamzxy * gyzy + & Gamxyy * gxxy + Gamyyy * gxyy + Gamzyy * gxzy) + & Gamxyy * gxyx + Gamyyy * gxyy + Gamzyy * gxyz + & Gamxxy * gyyx + Gamyxy * gyyy + Gamzxy * gyyz )+ & gupxz *( & TWO*(Gamxxy * gxzy + Gamyxy * gyzy + Gamzxy * gzzy + & Gamxyz * gxxy + Gamyyz * gxyy + Gamzyz * gxzy) + & Gamxyz * gxyx + Gamyyz * gxyy + Gamzyz * gxyz + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupyy *( & TWO*(Gamxyy * gxyy + Gamyyy * gyyy + Gamzyy * gyzy) + & Gamxyy * gyyx + Gamyyy * gyyy + Gamzyy * gyyz )+ & gupyz *( & TWO*(Gamxyy * gxzy + Gamyyy * gyzy + Gamzyy * gzzy + & Gamxyz * gxyy + Gamyyz * gyyy + Gamzyz * gyzy) + & Gamxyz * gyyx + Gamyyz * gyyy + Gamzyz * gyyz + & Gamxyy * gyzx + Gamyyy * gyzy + Gamzyy * gyzz )+ & gupzz *( & TWO*(Gamxyz * gxzy + Gamyyz * gyzy + Gamzyz * gzzy) + & Gamxyz * gyzx + Gamyyz * gyzy + Gamzyz * gyzz ) Rzz = - HALF * Rzz + & gxz * Gamxz+ gyz * Gamyz + gzz * Gamzz + & Gamxa * gxzz + Gamya * gyzz + Gamza * gzzz + & gupxx *( & TWO*(Gamxxz * gxxz + Gamyxz * gxyz + Gamzxz * gxzz) + & Gamxxz * gxzx + Gamyxz * gxzy + Gamzxz * gxzz )+ & gupxy *( & TWO*(Gamxxz * gxyz + Gamyxz * gyyz + Gamzxz * gyzz + & Gamxyz * gxxz + Gamyyz * gxyz + Gamzyz * gxzz) + & Gamxyz * gxzx + Gamyyz * gxzy + Gamzyz * gxzz + & Gamxxz * gyzx + Gamyxz * gyzy + Gamzxz * gyzz )+ & gupxz *( & TWO*(Gamxxz * gxzz + Gamyxz * gyzz + Gamzxz * gzzz + & Gamxzz * gxxz + Gamyzz * gxyz + Gamzzz * gxzz) + & Gamxzz * gxzx + Gamyzz * gxzy + Gamzzz * gxzz + & Gamxxz * gzzx + Gamyxz * gzzy + Gamzxz * gzzz )+ & gupyy *( & TWO*(Gamxyz * gxyz + Gamyyz * gyyz + Gamzyz * gyzz) + & Gamxyz * gyzx + Gamyyz * gyzy + Gamzyz * gyzz )+ & gupyz *( & TWO*(Gamxyz * gxzz + Gamyyz * gyzz + Gamzyz * gzzz + & Gamxzz * gxyz + Gamyzz * gyyz + Gamzzz * gyzz) + & Gamxzz * gyzx + Gamyzz * gyzy + Gamzzz * gyzz + & Gamxyz * gzzx + Gamyyz * gzzy + Gamzyz * gzzz )+ & gupzz *( & TWO*(Gamxzz * gxzz + Gamyzz * gyzz + Gamzzz * gzzz) + & Gamxzz * gzzx + Gamyzz * gzzy + Gamzzz * gzzz ) Rxy = HALF*( - Rxy + & gxx * Gamxy + gxy * Gamyy + gxz * Gamzy + & gxy * Gamxx + gyy * Gamyx + gyz * Gamzx + & Gamxa * gxyx + Gamya * gyyx + Gamza * gyzx + & Gamxa * gxxy + Gamya * gxyy + Gamza * gxzy )+ & gupxx *( & Gamxxx * gxxy + Gamyxx * gxyy + Gamzxx * gxzy + & Gamxxy * gxxx + Gamyxy * gxyx + Gamzxy * gxzx + & Gamxxx * gxyx + Gamyxx * gxyy + Gamzxx * gxyz )+ & gupxy *( & Gamxxx * gxyy + Gamyxx * gyyy + Gamzxx * gyzy + & Gamxxy * gxyx + Gamyxy * gyyx + Gamzxy * gyzx + & Gamxxy * gxyx + Gamyxy * gxyy + Gamzxy * gxyz + & Gamxxy * gxxy + Gamyxy * gxyy + Gamzxy * gxzy + & Gamxyy * gxxx + Gamyyy * gxyx + Gamzyy * gxzx + & Gamxxx * gyyx + Gamyxx * gyyy + Gamzxx * gyyz )+ & gupxz *( & Gamxxx * gxzy + Gamyxx * gyzy + Gamzxx * gzzy + & Gamxxy * gxzx + Gamyxy * gyzx + Gamzxy * gzzx + & Gamxxz * gxyx + Gamyxz * gxyy + Gamzxz * gxyz + & Gamxxz * gxxy + Gamyxz * gxyy + Gamzxz * gxzy + & Gamxyz * gxxx + Gamyyz * gxyx + Gamzyz * gxzx + & Gamxxx * gyzx + Gamyxx * gyzy + Gamzxx * gyzz )+ & gupyy *( & Gamxxy * gxyy + Gamyxy * gyyy + Gamzxy * gyzy + & Gamxyy * gxyx + Gamyyy * gyyx + Gamzyy * gyzx + & Gamxxy * gyyx + Gamyxy * gyyy + Gamzxy * gyyz )+ & gupyz *( & Gamxxy * gxzy + Gamyxy * gyzy + Gamzxy * gzzy + & Gamxyy * gxzx + Gamyyy * gyzx + Gamzyy * gzzx + & Gamxxz * gyyx + Gamyxz * gyyy + Gamzxz * gyyz + & Gamxxz * gxyy + Gamyxz * gyyy + Gamzxz * gyzy + & Gamxyz * gxyx + Gamyyz * gyyx + Gamzyz * gyzx + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupzz *( & Gamxxz * gxzy + Gamyxz * gyzy + Gamzxz * gzzy + & Gamxyz * gxzx + Gamyyz * gyzx + Gamzyz * gzzx + & Gamxxz * gyzx + Gamyxz * gyzy + Gamzxz * gyzz ) Rxz = HALF*( - Rxz + & gxx * Gamxz + gxy * Gamyz + gxz * Gamzz + & gxz * Gamxx + gyz * Gamyx + gzz * Gamzx + & Gamxa * gxzx + Gamya * gyzx + Gamza * gzzx + & Gamxa * gxxz + Gamya * gxyz + Gamza * gxzz )+ & gupxx *( & Gamxxx * gxxz + Gamyxx * gxyz + Gamzxx * gxzz + & Gamxxz * gxxx + Gamyxz * gxyx + Gamzxz * gxzx + & Gamxxx * gxzx + Gamyxx * gxzy + Gamzxx * gxzz )+ & gupxy *( & Gamxxx * gxyz + Gamyxx * gyyz + Gamzxx * gyzz + & Gamxxz * gxyx + Gamyxz * gyyx + Gamzxz * gyzx + & Gamxxy * gxzx + Gamyxy * gxzy + Gamzxy * gxzz + & Gamxxy * gxxz + Gamyxy * gxyz + Gamzxy * gxzz + & Gamxyz * gxxx + Gamyyz * gxyx + Gamzyz * gxzx + & Gamxxx * gyzx + Gamyxx * gyzy + Gamzxx * gyzz )+ & gupxz *( & Gamxxx * gxzz + Gamyxx * gyzz + Gamzxx * gzzz + & Gamxxz * gxzx + Gamyxz * gyzx + Gamzxz * gzzx + & Gamxxz * gxzx + Gamyxz * gxzy + Gamzxz * gxzz + & Gamxxz * gxxz + Gamyxz * gxyz + Gamzxz * gxzz + & Gamxzz * gxxx + Gamyzz * gxyx + Gamzzz * gxzx + & Gamxxx * gzzx + Gamyxx * gzzy + Gamzxx * gzzz )+ & gupyy *( & Gamxxy * gxyz + Gamyxy * gyyz + Gamzxy * gyzz + & Gamxyz * gxyx + Gamyyz * gyyx + Gamzyz * gyzx + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupyz *( & Gamxxy * gxzz + Gamyxy * gyzz + Gamzxy * gzzz + & Gamxyz * gxzx + Gamyyz * gyzx + Gamzyz * gzzx + & Gamxxz * gyzx + Gamyxz * gyzy + Gamzxz * gyzz + & Gamxxz * gxyz + Gamyxz * gyyz + Gamzxz * gyzz + & Gamxzz * gxyx + Gamyzz * gyyx + Gamzzz * gyzx + & Gamxxy * gzzx + Gamyxy * gzzy + Gamzxy * gzzz )+ & gupzz *( & Gamxxz * gxzz + Gamyxz * gyzz + Gamzxz * gzzz + & Gamxzz * gxzx + Gamyzz * gyzx + Gamzzz * gzzx + & Gamxxz * gzzx + Gamyxz * gzzy + Gamzxz * gzzz ) Ryz = HALF*( - Ryz + & gxy * Gamxz + gyy * Gamyz + gyz * Gamzz + & gxz * Gamxy + gyz * Gamyy + gzz * Gamzy + & Gamxa * gxzy + Gamya * gyzy + Gamza * gzzy + & Gamxa * gxyz + Gamya * gyyz + Gamza * gyzz )+ & gupxx *( & Gamxxy * gxxz + Gamyxy * gxyz + Gamzxy * gxzz + & Gamxxz * gxxy + Gamyxz * gxyy + Gamzxz * gxzy + & Gamxxy * gxzx + Gamyxy * gxzy + Gamzxy * gxzz )+ & gupxy *( & Gamxxy * gxyz + Gamyxy * gyyz + Gamzxy * gyzz + & Gamxxz * gxyy + Gamyxz * gyyy + Gamzxz * gyzy + & Gamxyy * gxzx + Gamyyy * gxzy + Gamzyy * gxzz + & Gamxyy * gxxz + Gamyyy * gxyz + Gamzyy * gxzz + & Gamxyz * gxxy + Gamyyz * gxyy + Gamzyz * gxzy + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupxz *( & Gamxxy * gxzz + Gamyxy * gyzz + Gamzxy * gzzz + & Gamxxz * gxzy + Gamyxz * gyzy + Gamzxz * gzzy + & Gamxyz * gxzx + Gamyyz * gxzy + Gamzyz * gxzz + & Gamxyz * gxxz + Gamyyz * gxyz + Gamzyz * gxzz + & Gamxzz * gxxy + Gamyzz * gxyy + Gamzzz * gxzy + & Gamxxy * gzzx + Gamyxy * gzzy + Gamzxy * gzzz )+ & gupyy *( & Gamxyy * gxyz + Gamyyy * gyyz + Gamzyy * gyzz + & Gamxyz * gxyy + Gamyyz * gyyy + Gamzyz * gyzy + & Gamxyy * gyzx + Gamyyy * gyzy + Gamzyy * gyzz )+ & gupyz *( & Gamxyy * gxzz + Gamyyy * gyzz + Gamzyy * gzzz + & Gamxyz * gxzy + Gamyyz * gyzy + Gamzyz * gzzy + & Gamxyz * gyzx + Gamyyz * gyzy + Gamzyz * gyzz + & Gamxyz * gxyz + Gamyyz * gyyz + Gamzyz * gyzz + & Gamxzz * gxyy + Gamyzz * gyyy + Gamzzz * gyzy + & Gamxyy * gzzx + Gamyyy * gzzy + Gamzyy * gzzz )+ & gupzz *( & Gamxyz * gxzz + Gamyyz * gyzz + Gamzyz * gzzz + & Gamxzz * gxzy + Gamyzz * gyzy + Gamzzz * gzzy + & Gamxyz * gzzx + Gamyyz * gzzy + Gamzyz * gzzz ) !covariant second derivative of chi respect to tilted metric call fdderivs(ex,chi,fxx,fxy,fxz,fyy,fyz,fzz,X,Y,Z,SYM,SYM,SYM,Symmetry,0) fxx = fxx - Gamxxx * chix - Gamyxx * chiy - Gamzxx * chiz fxy = fxy - Gamxxy * chix - Gamyxy * chiy - Gamzxy * chiz fxz = fxz - Gamxxz * chix - Gamyxz * chiy - Gamzxz * chiz fyy = fyy - Gamxyy * chix - Gamyyy * chiy - Gamzyy * chiz fyz = fyz - Gamxyz * chix - Gamyyz * chiy - Gamzyz * chiz fzz = fzz - Gamxzz * chix - Gamyzz * chiy - Gamzzz * chiz ! Store D^l D_l chi - 3/(2*chi) D^l chi D_l chi in f f = gupxx * ( fxx - F3o2/chin1 * chix * chix ) + & gupyy * ( fyy - F3o2/chin1 * chiy * chiy ) + & gupzz * ( fzz - F3o2/chin1 * chiz * chiz ) + & TWO * gupxy * ( fxy - F3o2/chin1 * chix * chiy ) + & TWO * gupxz * ( fxz - F3o2/chin1 * chix * chiz ) + & TWO * gupyz * ( fyz - F3o2/chin1 * chiy * chiz ) ! Add chi part to Ricci tensor: Rxx = Rxx + (fxx - chix*chix/chin1/TWO + gxx * f)/chin1/TWO Ryy = Ryy + (fyy - chiy*chiy/chin1/TWO + gyy * f)/chin1/TWO Rzz = Rzz + (fzz - chiz*chiz/chin1/TWO + gzz * f)/chin1/TWO Rxy = Rxy + (fxy - chix*chiy/chin1/TWO + gxy * f)/chin1/TWO Rxz = Rxz + (fxz - chix*chiz/chin1/TWO + gxz * f)/chin1/TWO Ryz = Ryz + (fyz - chiy*chiz/chin1/TWO + gyz * f)/chin1/TWO gxxx = (gupxx * chix + gupxy * chiy + gupxz * chiz)/chin1 gxxy = (gupxy * chix + gupyy * chiy + gupyz * chiz)/chin1 gxxz = (gupxz * chix + gupyz * chiy + gupzz * chiz)/chin1 ! now get physical second kind of connection Gamxxx = Gamxxx - ( (chix + chix)/chin1 - gxx * gxxx )*HALF Gamyxx = Gamyxx - ( - gxx * gxxy )*HALF Gamzxx = Gamzxx - ( - gxx * gxxz )*HALF Gamxyy = Gamxyy - ( - gyy * gxxx )*HALF Gamyyy = Gamyyy - ( (chiy + chiy)/chin1 - gyy * gxxy )*HALF Gamzyy = Gamzyy - ( - gyy * gxxz )*HALF Gamxzz = Gamxzz - ( - gzz * gxxx )*HALF Gamyzz = Gamyzz - ( - gzz * gxxy )*HALF Gamzzz = Gamzzz - ( (chiz + chiz)/chin1 - gzz * gxxz )*HALF Gamxxy = Gamxxy - ( chiy /chin1 - gxy * gxxx )*HALF Gamyxy = Gamyxy - ( chix /chin1 - gxy * gxxy )*HALF Gamzxy = Gamzxy - ( - gxy * gxxz )*HALF Gamxxz = Gamxxz - ( chiz /chin1 - gxz * gxxx )*HALF Gamyxz = Gamyxz - ( - gxz * gxxy )*HALF Gamzxz = Gamzxz - ( chix /chin1 - gxz * gxxz )*HALF Gamxyz = Gamxyz - ( - gyz * gxxx )*HALF Gamyyz = Gamyyz - ( chiz /chin1 - gyz * gxxy )*HALF Gamzyz = Gamzyz - ( chiy /chin1 - gyz * gxxz )*HALF return end subroutine ricci_gamma !---------------------------------------------------------------------------- subroutine ricci_gamma_ss(ex,crho,sigma,R,X, Y, Z, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz, & chi, & dxx , gxy , gxz , dyy , gyz , dzz,& Gamx , Gamy , Gamz , & Gamxxx,Gamxxy,Gamxxz,Gamxyy,Gamxyz,Gamxzz,& Gamyxx,Gamyxy,Gamyxz,Gamyyy,Gamyyz,Gamyzz,& Gamzxx,Gamzxy,Gamzxz,Gamzyy,Gamzyz,Gamzzz,& Rxx,Rxy,Rxz,Ryy,Ryz,Rzz,& Symmetry,Lev,sst) implicit none !~~~~~~> Input parameters: integer,intent(in ):: ex(1:3), Symmetry,Lev,sst double precision,intent(in),dimension(ex(1))::crho double precision,intent(in),dimension(ex(2))::sigma double precision,intent(in),dimension(ex(3))::R real*8, intent(in ),dimension(ex(1),ex(2),ex(3)):: X,Y,Z double precision,intent(in),dimension(ex(1),ex(2),ex(3))::drhodx, drhody, drhodz double precision,intent(in),dimension(ex(1),ex(2),ex(3))::dsigmadx,dsigmady,dsigmadz double precision,intent(in),dimension(ex(1),ex(2),ex(3))::dRdx,dRdy,dRdz double precision,intent(in),dimension(ex(1),ex(2),ex(3))::drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz double precision,intent(in),dimension(ex(1),ex(2),ex(3))::dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz double precision,intent(in),dimension(ex(1),ex(2),ex(3))::dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz real*8, dimension(ex(1),ex(2),ex(3)),intent(in ) :: chi real*8, dimension(ex(1),ex(2),ex(3)),intent(in ) :: dxx,gxy,gxz,dyy,gyz,dzz real*8, dimension(ex(1),ex(2),ex(3)),intent(in ) :: Gamx,Gamy,Gamz ! when out, physical second kind of connection real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamxxx, Gamxxy, Gamxxz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamxyy, Gamxyz, Gamxzz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamyxx, Gamyxy, Gamyxz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamyyy, Gamyyz, Gamyzz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamzxx, Gamzxy, Gamzxz real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Gamzyy, Gamzyz, Gamzzz ! when out, physical Ricci tensor real*8, dimension(ex(1),ex(2),ex(3)),intent(out) :: Rxx,Rxy,Rxz,Ryy,Ryz,Rzz !~~~~~~> Other variables: real*8, dimension(ex(1),ex(2),ex(3)) :: gxx,gyy,gzz real*8, dimension(ex(1),ex(2),ex(3)) :: chix,chiy,chiz real*8, dimension(ex(1),ex(2),ex(3)) :: gxxx,gxyx,gxzx,gyyx,gyzx,gzzx real*8, dimension(ex(1),ex(2),ex(3)) :: gxxy,gxyy,gxzy,gyyy,gyzy,gzzy real*8, dimension(ex(1),ex(2),ex(3)) :: gxxz,gxyz,gxzz,gyyz,gyzz,gzzz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamxx,Gamxy,Gamxz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamyx,Gamyy,Gamyz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamzx,Gamzy,Gamzz real*8, dimension(ex(1),ex(2),ex(3)) :: f,fxx,fxy,fxz,fyy,fyz,fzz real*8, dimension(ex(1),ex(2),ex(3)) :: Gamxa,Gamya,Gamza,chin1 real*8, dimension(ex(1),ex(2),ex(3)) :: gupxx,gupxy,gupxz real*8, dimension(ex(1),ex(2),ex(3)) :: gupyy,gupyz,gupzz real*8, parameter :: ONE = 1.D0, TWO = 2.D0, FOUR = 4.D0 real*8, parameter :: HALF = 0.5D0, F2o3 = 2.d0/3.d0, F3o2 = 1.5d0 real*8, parameter :: SYM = 1.D0, ANTI= - 1.D0 chin1 = chi + ONE gxx = dxx + ONE gyy = dyy + ONE gzz = dzz + ONE call fderivs_shc(ex,chi,chix,chiy,chiz,crho,sigma,R, SYM, SYM,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,dxx,gxxx,gxxy,gxxz,crho,sigma,R, SYM, SYM,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,gxy,gxyx,gxyy,gxyz,crho,sigma,R,ANTI,ANTI,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,gxz,gxzx,gxzy,gxzz,crho,sigma,R,ANTI,SYM ,ANTI,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,dyy,gyyx,gyyy,gyyz,crho,sigma,R, SYM, SYM,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,gyz,gyzx,gyzy,gyzz,crho,sigma,R,SYM ,ANTI,ANTI,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,dzz,gzzx,gzzy,gzzz,crho,sigma,R, SYM, SYM,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) ! invert tilted metric gupzz = gxx * gyy * gzz + gxy * gyz * gxz + gxz * gxy * gyz - & gxz * gyy * gxz - gxy * gxy * gzz - gxx * gyz * gyz gupxx = ( gyy * gzz - gyz * gyz ) / gupzz gupxy = - ( gxy * gzz - gyz * gxz ) / gupzz gupxz = ( gxy * gyz - gyy * gxz ) / gupzz gupyy = ( gxx * gzz - gxz * gxz ) / gupzz gupyz = - ( gxx * gyz - gxy * gxz ) / gupzz gupzz = ( gxx * gyy - gxy * gxy ) / gupzz ! second kind of connection Gamxxx =HALF*( gupxx*gxxx + gupxy*(TWO*gxyx - gxxy ) + gupxz*(TWO*gxzx - gxxz )) Gamyxx =HALF*( gupxy*gxxx + gupyy*(TWO*gxyx - gxxy ) + gupyz*(TWO*gxzx - gxxz )) Gamzxx =HALF*( gupxz*gxxx + gupyz*(TWO*gxyx - gxxy ) + gupzz*(TWO*gxzx - gxxz )) Gamxyy =HALF*( gupxx*(TWO*gxyy - gyyx ) + gupxy*gyyy + gupxz*(TWO*gyzy - gyyz )) Gamyyy =HALF*( gupxy*(TWO*gxyy - gyyx ) + gupyy*gyyy + gupyz*(TWO*gyzy - gyyz )) Gamzyy =HALF*( gupxz*(TWO*gxyy - gyyx ) + gupyz*gyyy + gupzz*(TWO*gyzy - gyyz )) Gamxzz =HALF*( gupxx*(TWO*gxzz - gzzx ) + gupxy*(TWO*gyzz - gzzy ) + gupxz*gzzz) Gamyzz =HALF*( gupxy*(TWO*gxzz - gzzx ) + gupyy*(TWO*gyzz - gzzy ) + gupyz*gzzz) Gamzzz =HALF*( gupxz*(TWO*gxzz - gzzx ) + gupyz*(TWO*gyzz - gzzy ) + gupzz*gzzz) Gamxxy =HALF*( gupxx*gxxy + gupxy*gyyx + gupxz*( gxzy + gyzx - gxyz ) ) Gamyxy =HALF*( gupxy*gxxy + gupyy*gyyx + gupyz*( gxzy + gyzx - gxyz ) ) Gamzxy =HALF*( gupxz*gxxy + gupyz*gyyx + gupzz*( gxzy + gyzx - gxyz ) ) Gamxxz =HALF*( gupxx*gxxz + gupxy*( gxyz + gyzx - gxzy ) + gupxz*gzzx ) Gamyxz =HALF*( gupxy*gxxz + gupyy*( gxyz + gyzx - gxzy ) + gupyz*gzzx ) Gamzxz =HALF*( gupxz*gxxz + gupyz*( gxyz + gyzx - gxzy ) + gupzz*gzzx ) Gamxyz =HALF*( gupxx*( gxyz + gxzy - gyzx ) + gupxy*gyyz + gupxz*gzzy ) Gamyyz =HALF*( gupxy*( gxyz + gxzy - gyzx ) + gupyy*gyyz + gupyz*gzzy ) Gamzyz =HALF*( gupxz*( gxyz + gxzy - gyzx ) + gupyz*gyyz + gupzz*gzzy ) Gamxa = gupxx * Gamxxx + gupyy * Gamxyy + gupzz * Gamxzz + & TWO*( gupxy * Gamxxy + gupxz * Gamxxz + gupyz * Gamxyz ) Gamya = gupxx * Gamyxx + gupyy * Gamyyy + gupzz * Gamyzz + & TWO*( gupxy * Gamyxy + gupxz * Gamyxz + gupyz * Gamyyz ) Gamza = gupxx * Gamzxx + gupyy * Gamzyy + gupzz * Gamzzz + & TWO*( gupxy * Gamzxy + gupxz * Gamzxz + gupyz * Gamzyz ) call fderivs_shc(ex,Gamx,Gamxx,Gamxy,Gamxz,crho,sigma,R,ANTI,SYM ,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,Gamy,Gamyx,Gamyy,Gamyz,crho,sigma,R,SYM ,ANTI,SYM,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) call fderivs_shc(ex,Gamz,Gamzx,Gamzy,Gamzz,crho,sigma,R,SYM ,SYM ,ANTI,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz) !first kind of connection stored in gij,k gxxx = gxx * Gamxxx + gxy * Gamyxx + gxz * Gamzxx gxyx = gxx * Gamxxy + gxy * Gamyxy + gxz * Gamzxy gxzx = gxx * Gamxxz + gxy * Gamyxz + gxz * Gamzxz gyyx = gxx * Gamxyy + gxy * Gamyyy + gxz * Gamzyy gyzx = gxx * Gamxyz + gxy * Gamyyz + gxz * Gamzyz gzzx = gxx * Gamxzz + gxy * Gamyzz + gxz * Gamzzz gxxy = gxy * Gamxxx + gyy * Gamyxx + gyz * Gamzxx gxyy = gxy * Gamxxy + gyy * Gamyxy + gyz * Gamzxy gxzy = gxy * Gamxxz + gyy * Gamyxz + gyz * Gamzxz gyyy = gxy * Gamxyy + gyy * Gamyyy + gyz * Gamzyy gyzy = gxy * Gamxyz + gyy * Gamyyz + gyz * Gamzyz gzzy = gxy * Gamxzz + gyy * Gamyzz + gyz * Gamzzz gxxz = gxz * Gamxxx + gyz * Gamyxx + gzz * Gamzxx gxyz = gxz * Gamxxy + gyz * Gamyxy + gzz * Gamzxy gxzz = gxz * Gamxxz + gyz * Gamyxz + gzz * Gamzxz gyyz = gxz * Gamxyy + gyz * Gamyyy + gzz * Gamzyy gyzz = gxz * Gamxyz + gyz * Gamyyz + gzz * Gamzyz gzzz = gxz * Gamxzz + gyz * Gamyzz + gzz * Gamzzz !compute Ricci tensor for tilted metric call fdderivs_shc(ex,dxx,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R, SYM, SYM,SYM ,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) Rxx = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs_shc(ex,dyy,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R, SYM, SYM,SYM ,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) Ryy = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs_shc(ex,dzz,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R, SYM, SYM,SYM ,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) Rzz = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs_shc(ex,gxy,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R,ANTI,ANTI,SYM ,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) Rxy = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs_shc(ex,gxz,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R,ANTI,SYM ,ANTI,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) Rxz = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO call fdderivs_shc(ex,gyz,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R,SYM ,ANTI,ANTI,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) Ryz = gupxx * fxx + gupyy * fyy + gupzz * fzz + & ( gupxy * fxy + gupxz * fxz + gupyz * fyz ) * TWO Rxx = - HALF * Rxx + & gxx * Gamxx+ gxy * Gamyx + gxz * Gamzx + & Gamxa * gxxx + Gamya * gxyx + Gamza * gxzx + & gupxx *( & TWO*(Gamxxx * gxxx + Gamyxx * gxyx + Gamzxx * gxzx) + & Gamxxx * gxxx + Gamyxx * gxxy + Gamzxx * gxxz )+ & gupxy *( & TWO*(Gamxxx * gxyx + Gamyxx * gyyx + Gamzxx * gyzx + & Gamxxy * gxxx + Gamyxy * gxyx + Gamzxy * gxzx) + & Gamxxy * gxxx + Gamyxy * gxxy + Gamzxy * gxxz + & Gamxxx * gxyx + Gamyxx * gxyy + Gamzxx * gxyz )+ & gupxz *( & TWO*(Gamxxx * gxzx + Gamyxx * gyzx + Gamzxx * gzzx + & Gamxxz * gxxx + Gamyxz * gxyx + Gamzxz * gxzx) + & Gamxxz * gxxx + Gamyxz * gxxy + Gamzxz * gxxz + & Gamxxx * gxzx + Gamyxx * gxzy + Gamzxx * gxzz )+ & gupyy *( & TWO*(Gamxxy * gxyx + Gamyxy * gyyx + Gamzxy * gyzx) + & Gamxxy * gxyx + Gamyxy * gxyy + Gamzxy * gxyz )+ & gupyz *( & TWO*(Gamxxy * gxzx + Gamyxy * gyzx + Gamzxy * gzzx + & Gamxxz * gxyx + Gamyxz * gyyx + Gamzxz * gyzx) + & Gamxxz * gxyx + Gamyxz * gxyy + Gamzxz * gxyz + & Gamxxy * gxzx + Gamyxy * gxzy + Gamzxy * gxzz )+ & gupzz *( & TWO*(Gamxxz * gxzx + Gamyxz * gyzx + Gamzxz * gzzx) + & Gamxxz * gxzx + Gamyxz * gxzy + Gamzxz * gxzz ) Ryy = - HALF * Ryy + & gxy * Gamxy+ gyy * Gamyy + gyz * Gamzy + & Gamxa * gxyy + Gamya * gyyy + Gamza * gyzy + & gupxx *( & TWO*(Gamxxy * gxxy + Gamyxy * gxyy + Gamzxy * gxzy) + & Gamxxy * gxyx + Gamyxy * gxyy + Gamzxy * gxyz )+ & gupxy *( & TWO*(Gamxxy * gxyy + Gamyxy * gyyy + Gamzxy * gyzy + & Gamxyy * gxxy + Gamyyy * gxyy + Gamzyy * gxzy) + & Gamxyy * gxyx + Gamyyy * gxyy + Gamzyy * gxyz + & Gamxxy * gyyx + Gamyxy * gyyy + Gamzxy * gyyz )+ & gupxz *( & TWO*(Gamxxy * gxzy + Gamyxy * gyzy + Gamzxy * gzzy + & Gamxyz * gxxy + Gamyyz * gxyy + Gamzyz * gxzy) + & Gamxyz * gxyx + Gamyyz * gxyy + Gamzyz * gxyz + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupyy *( & TWO*(Gamxyy * gxyy + Gamyyy * gyyy + Gamzyy * gyzy) + & Gamxyy * gyyx + Gamyyy * gyyy + Gamzyy * gyyz )+ & gupyz *( & TWO*(Gamxyy * gxzy + Gamyyy * gyzy + Gamzyy * gzzy + & Gamxyz * gxyy + Gamyyz * gyyy + Gamzyz * gyzy) + & Gamxyz * gyyx + Gamyyz * gyyy + Gamzyz * gyyz + & Gamxyy * gyzx + Gamyyy * gyzy + Gamzyy * gyzz )+ & gupzz *( & TWO*(Gamxyz * gxzy + Gamyyz * gyzy + Gamzyz * gzzy) + & Gamxyz * gyzx + Gamyyz * gyzy + Gamzyz * gyzz ) Rzz = - HALF * Rzz + & gxz * Gamxz+ gyz * Gamyz + gzz * Gamzz + & Gamxa * gxzz + Gamya * gyzz + Gamza * gzzz + & gupxx *( & TWO*(Gamxxz * gxxz + Gamyxz * gxyz + Gamzxz * gxzz) + & Gamxxz * gxzx + Gamyxz * gxzy + Gamzxz * gxzz )+ & gupxy *( & TWO*(Gamxxz * gxyz + Gamyxz * gyyz + Gamzxz * gyzz + & Gamxyz * gxxz + Gamyyz * gxyz + Gamzyz * gxzz) + & Gamxyz * gxzx + Gamyyz * gxzy + Gamzyz * gxzz + & Gamxxz * gyzx + Gamyxz * gyzy + Gamzxz * gyzz )+ & gupxz *( & TWO*(Gamxxz * gxzz + Gamyxz * gyzz + Gamzxz * gzzz + & Gamxzz * gxxz + Gamyzz * gxyz + Gamzzz * gxzz) + & Gamxzz * gxzx + Gamyzz * gxzy + Gamzzz * gxzz + & Gamxxz * gzzx + Gamyxz * gzzy + Gamzxz * gzzz )+ & gupyy *( & TWO*(Gamxyz * gxyz + Gamyyz * gyyz + Gamzyz * gyzz) + & Gamxyz * gyzx + Gamyyz * gyzy + Gamzyz * gyzz )+ & gupyz *( & TWO*(Gamxyz * gxzz + Gamyyz * gyzz + Gamzyz * gzzz + & Gamxzz * gxyz + Gamyzz * gyyz + Gamzzz * gyzz) + & Gamxzz * gyzx + Gamyzz * gyzy + Gamzzz * gyzz + & Gamxyz * gzzx + Gamyyz * gzzy + Gamzyz * gzzz )+ & gupzz *( & TWO*(Gamxzz * gxzz + Gamyzz * gyzz + Gamzzz * gzzz) + & Gamxzz * gzzx + Gamyzz * gzzy + Gamzzz * gzzz ) Rxy = HALF*( - Rxy + & gxx * Gamxy + gxy * Gamyy + gxz * Gamzy + & gxy * Gamxx + gyy * Gamyx + gyz * Gamzx + & Gamxa * gxyx + Gamya * gyyx + Gamza * gyzx + & Gamxa * gxxy + Gamya * gxyy + Gamza * gxzy )+ & gupxx *( & Gamxxx * gxxy + Gamyxx * gxyy + Gamzxx * gxzy + & Gamxxy * gxxx + Gamyxy * gxyx + Gamzxy * gxzx + & Gamxxx * gxyx + Gamyxx * gxyy + Gamzxx * gxyz )+ & gupxy *( & Gamxxx * gxyy + Gamyxx * gyyy + Gamzxx * gyzy + & Gamxxy * gxyx + Gamyxy * gyyx + Gamzxy * gyzx + & Gamxxy * gxyx + Gamyxy * gxyy + Gamzxy * gxyz + & Gamxxy * gxxy + Gamyxy * gxyy + Gamzxy * gxzy + & Gamxyy * gxxx + Gamyyy * gxyx + Gamzyy * gxzx + & Gamxxx * gyyx + Gamyxx * gyyy + Gamzxx * gyyz )+ & gupxz *( & Gamxxx * gxzy + Gamyxx * gyzy + Gamzxx * gzzy + & Gamxxy * gxzx + Gamyxy * gyzx + Gamzxy * gzzx + & Gamxxz * gxyx + Gamyxz * gxyy + Gamzxz * gxyz + & Gamxxz * gxxy + Gamyxz * gxyy + Gamzxz * gxzy + & Gamxyz * gxxx + Gamyyz * gxyx + Gamzyz * gxzx + & Gamxxx * gyzx + Gamyxx * gyzy + Gamzxx * gyzz )+ & gupyy *( & Gamxxy * gxyy + Gamyxy * gyyy + Gamzxy * gyzy + & Gamxyy * gxyx + Gamyyy * gyyx + Gamzyy * gyzx + & Gamxxy * gyyx + Gamyxy * gyyy + Gamzxy * gyyz )+ & gupyz *( & Gamxxy * gxzy + Gamyxy * gyzy + Gamzxy * gzzy + & Gamxyy * gxzx + Gamyyy * gyzx + Gamzyy * gzzx + & Gamxxz * gyyx + Gamyxz * gyyy + Gamzxz * gyyz + & Gamxxz * gxyy + Gamyxz * gyyy + Gamzxz * gyzy + & Gamxyz * gxyx + Gamyyz * gyyx + Gamzyz * gyzx + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupzz *( & Gamxxz * gxzy + Gamyxz * gyzy + Gamzxz * gzzy + & Gamxyz * gxzx + Gamyyz * gyzx + Gamzyz * gzzx + & Gamxxz * gyzx + Gamyxz * gyzy + Gamzxz * gyzz ) Rxz = HALF*( - Rxz + & gxx * Gamxz + gxy * Gamyz + gxz * Gamzz + & gxz * Gamxx + gyz * Gamyx + gzz * Gamzx + & Gamxa * gxzx + Gamya * gyzx + Gamza * gzzx + & Gamxa * gxxz + Gamya * gxyz + Gamza * gxzz )+ & gupxx *( & Gamxxx * gxxz + Gamyxx * gxyz + Gamzxx * gxzz + & Gamxxz * gxxx + Gamyxz * gxyx + Gamzxz * gxzx + & Gamxxx * gxzx + Gamyxx * gxzy + Gamzxx * gxzz )+ & gupxy *( & Gamxxx * gxyz + Gamyxx * gyyz + Gamzxx * gyzz + & Gamxxz * gxyx + Gamyxz * gyyx + Gamzxz * gyzx + & Gamxxy * gxzx + Gamyxy * gxzy + Gamzxy * gxzz + & Gamxxy * gxxz + Gamyxy * gxyz + Gamzxy * gxzz + & Gamxyz * gxxx + Gamyyz * gxyx + Gamzyz * gxzx + & Gamxxx * gyzx + Gamyxx * gyzy + Gamzxx * gyzz )+ & gupxz *( & Gamxxx * gxzz + Gamyxx * gyzz + Gamzxx * gzzz + & Gamxxz * gxzx + Gamyxz * gyzx + Gamzxz * gzzx + & Gamxxz * gxzx + Gamyxz * gxzy + Gamzxz * gxzz + & Gamxxz * gxxz + Gamyxz * gxyz + Gamzxz * gxzz + & Gamxzz * gxxx + Gamyzz * gxyx + Gamzzz * gxzx + & Gamxxx * gzzx + Gamyxx * gzzy + Gamzxx * gzzz )+ & gupyy *( & Gamxxy * gxyz + Gamyxy * gyyz + Gamzxy * gyzz + & Gamxyz * gxyx + Gamyyz * gyyx + Gamzyz * gyzx + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupyz *( & Gamxxy * gxzz + Gamyxy * gyzz + Gamzxy * gzzz + & Gamxyz * gxzx + Gamyyz * gyzx + Gamzyz * gzzx + & Gamxxz * gyzx + Gamyxz * gyzy + Gamzxz * gyzz + & Gamxxz * gxyz + Gamyxz * gyyz + Gamzxz * gyzz + & Gamxzz * gxyx + Gamyzz * gyyx + Gamzzz * gyzx + & Gamxxy * gzzx + Gamyxy * gzzy + Gamzxy * gzzz )+ & gupzz *( & Gamxxz * gxzz + Gamyxz * gyzz + Gamzxz * gzzz + & Gamxzz * gxzx + Gamyzz * gyzx + Gamzzz * gzzx + & Gamxxz * gzzx + Gamyxz * gzzy + Gamzxz * gzzz ) Ryz = HALF*( - Ryz + & gxy * Gamxz + gyy * Gamyz + gyz * Gamzz + & gxz * Gamxy + gyz * Gamyy + gzz * Gamzy + & Gamxa * gxzy + Gamya * gyzy + Gamza * gzzy + & Gamxa * gxyz + Gamya * gyyz + Gamza * gyzz )+ & gupxx *( & Gamxxy * gxxz + Gamyxy * gxyz + Gamzxy * gxzz + & Gamxxz * gxxy + Gamyxz * gxyy + Gamzxz * gxzy + & Gamxxy * gxzx + Gamyxy * gxzy + Gamzxy * gxzz )+ & gupxy *( & Gamxxy * gxyz + Gamyxy * gyyz + Gamzxy * gyzz + & Gamxxz * gxyy + Gamyxz * gyyy + Gamzxz * gyzy + & Gamxyy * gxzx + Gamyyy * gxzy + Gamzyy * gxzz + & Gamxyy * gxxz + Gamyyy * gxyz + Gamzyy * gxzz + & Gamxyz * gxxy + Gamyyz * gxyy + Gamzyz * gxzy + & Gamxxy * gyzx + Gamyxy * gyzy + Gamzxy * gyzz )+ & gupxz *( & Gamxxy * gxzz + Gamyxy * gyzz + Gamzxy * gzzz + & Gamxxz * gxzy + Gamyxz * gyzy + Gamzxz * gzzy + & Gamxyz * gxzx + Gamyyz * gxzy + Gamzyz * gxzz + & Gamxyz * gxxz + Gamyyz * gxyz + Gamzyz * gxzz + & Gamxzz * gxxy + Gamyzz * gxyy + Gamzzz * gxzy + & Gamxxy * gzzx + Gamyxy * gzzy + Gamzxy * gzzz )+ & gupyy *( & Gamxyy * gxyz + Gamyyy * gyyz + Gamzyy * gyzz + & Gamxyz * gxyy + Gamyyz * gyyy + Gamzyz * gyzy + & Gamxyy * gyzx + Gamyyy * gyzy + Gamzyy * gyzz )+ & gupyz *( & Gamxyy * gxzz + Gamyyy * gyzz + Gamzyy * gzzz + & Gamxyz * gxzy + Gamyyz * gyzy + Gamzyz * gzzy + & Gamxyz * gyzx + Gamyyz * gyzy + Gamzyz * gyzz + & Gamxyz * gxyz + Gamyyz * gyyz + Gamzyz * gyzz + & Gamxzz * gxyy + Gamyzz * gyyy + Gamzzz * gyzy + & Gamxyy * gzzx + Gamyyy * gzzy + Gamzyy * gzzz )+ & gupzz *( & Gamxyz * gxzz + Gamyyz * gyzz + Gamzyz * gzzz + & Gamxzz * gxzy + Gamyzz * gyzy + Gamzzz * gzzy + & Gamxyz * gzzx + Gamyyz * gzzy + Gamzyz * gzzz ) !covariant second derivative of chi respect to tilted metric call fdderivs_shc(ex,chi,fxx,fxy,fxz,fyy,fyz,fzz,crho,sigma,R,SYM ,SYM ,SYM ,Symmetry,Lev,sst, & drhodx, drhody, drhodz, & dsigmadx,dsigmady,dsigmadz, & dRdx,dRdy,dRdz, & drhodxx,drhodxy,drhodxz,drhodyy,drhodyz,drhodzz, & dsigmadxx,dsigmadxy,dsigmadxz,dsigmadyy,dsigmadyz,dsigmadzz, & dRdxx,dRdxy,dRdxz,dRdyy,dRdyz,dRdzz) fxx = fxx - Gamxxx * chix - Gamyxx * chiy - Gamzxx * chiz fxy = fxy - Gamxxy * chix - Gamyxy * chiy - Gamzxy * chiz fxz = fxz - Gamxxz * chix - Gamyxz * chiy - Gamzxz * chiz fyy = fyy - Gamxyy * chix - Gamyyy * chiy - Gamzyy * chiz fyz = fyz - Gamxyz * chix - Gamyyz * chiy - Gamzyz * chiz fzz = fzz - Gamxzz * chix - Gamyzz * chiy - Gamzzz * chiz ! Store D^l D_l chi - 3/(2*chi) D^l chi D_l chi in f f = gupxx * ( fxx - F3o2/chin1 * chix * chix ) + & gupyy * ( fyy - F3o2/chin1 * chiy * chiy ) + & gupzz * ( fzz - F3o2/chin1 * chiz * chiz ) + & TWO * gupxy * ( fxy - F3o2/chin1 * chix * chiy ) + & TWO * gupxz * ( fxz - F3o2/chin1 * chix * chiz ) + & TWO * gupyz * ( fyz - F3o2/chin1 * chiy * chiz ) ! Add chi part to Ricci tensor: Rxx = Rxx + (fxx - chix*chix/chin1/TWO + gxx * f)/chin1/TWO Ryy = Ryy + (fyy - chiy*chiy/chin1/TWO + gyy * f)/chin1/TWO Rzz = Rzz + (fzz - chiz*chiz/chin1/TWO + gzz * f)/chin1/TWO Rxy = Rxy + (fxy - chix*chiy/chin1/TWO + gxy * f)/chin1/TWO Rxz = Rxz + (fxz - chix*chiz/chin1/TWO + gxz * f)/chin1/TWO Ryz = Ryz + (fyz - chiy*chiz/chin1/TWO + gyz * f)/chin1/TWO gxxx = (gupxx * chix + gupxy * chiy + gupxz * chiz)/chin1 gxxy = (gupxy * chix + gupyy * chiy + gupyz * chiz)/chin1 gxxz = (gupxz * chix + gupyz * chiy + gupzz * chiz)/chin1 ! now get physical second kind of connection Gamxxx = Gamxxx - ( (chix + chix)/chin1 - gxx * gxxx )*HALF Gamyxx = Gamyxx - ( - gxx * gxxy )*HALF Gamzxx = Gamzxx - ( - gxx * gxxz )*HALF Gamxyy = Gamxyy - ( - gyy * gxxx )*HALF Gamyyy = Gamyyy - ( (chiy + chiy)/chin1 - gyy * gxxy )*HALF Gamzyy = Gamzyy - ( - gyy * gxxz )*HALF Gamxzz = Gamxzz - ( - gzz * gxxx )*HALF Gamyzz = Gamyzz - ( - gzz * gxxy )*HALF Gamzzz = Gamzzz - ( (chiz + chiz)/chin1 - gzz * gxxz )*HALF Gamxxy = Gamxxy - ( chiy /chin1 - gxy * gxxx )*HALF Gamyxy = Gamyxy - ( chix /chin1 - gxy * gxxy )*HALF Gamzxy = Gamzxy - ( - gxy * gxxz )*HALF Gamxxz = Gamxxz - ( chiz /chin1 - gxz * gxxx )*HALF Gamyxz = Gamyxz - ( - gxz * gxxy )*HALF Gamzxz = Gamzxz - ( chix /chin1 - gxz * gxxz )*HALF Gamxyz = Gamxyz - ( - gyz * gxxx )*HALF Gamyyz = Gamyyz - ( chiz /chin1 - gyz * gxxy )*HALF Gamzyz = Gamzyz - ( chiy /chin1 - gyz * gxxz )*HALF return end subroutine ricci_gamma_ss