From f200066992db3a32d01607df031445b0ce3a05c3 Mon Sep 17 00:00:00 2001 From: schnetter Date: Tue, 26 Apr 2005 14:48:24 +0000 Subject: Add routines for derivatives of the Ricci tensor git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinUtils/TGRtensor/trunk@23 b716e942-a2de-43ad-8f52-f3dfe468e4e7 --- src/ricci.F90 | 72 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 69 insertions(+), 3 deletions(-) diff --git a/src/ricci.F90 b/src/ricci.F90 index ebe01fc..f14eb49 100644 --- a/src/ricci.F90 +++ b/src/ricci.F90 @@ -7,12 +7,14 @@ module ricci private public calc_connections public calc_connectionderivs + public calc_connectionderivs2 public calc_ricci + public calc_ricciderivs contains - subroutine calc_connections (gu, dg, gamma) - CCTK_REAL, intent(in) :: gu(3,3), dg(3,3,3) + subroutine calc_connections (gu, dgg, gamma) + CCTK_REAL, intent(in) :: gu(3,3), dgg(3,3,3) CCTK_REAL, intent(out) :: gamma(3,3,3) integer :: i,j,k,l ! Gamma^i_jk = 1/2 g^il (g_lj,k + g_lk,j - g_jk,l) @@ -22,7 +24,7 @@ contains gamma(i,j,k) = 0 do l=1,3 gamma(i,j,k) = gamma(i,j,k) + 0.5d0 * gu(i,l) & - * (dg(l,j,k) + dg(l,k,j) - dg(j,k,l)) + * (dgg(l,j,k) + dgg(l,k,j) - dgg(j,k,l)) end do end do end do @@ -51,6 +53,32 @@ contains end do end subroutine calc_connectionderivs + subroutine calc_connectionderivs2 (gu, dgg, dgu, ddgg, ddgu, dddgg, ddgamma) + CCTK_REAL, intent(in) :: gu(3,3), dgg(3,3,3), dgu(3,3,3), ddgg(3,3,3,3), ddgu(3,3,3,3), dddgg(3,3,3,3,3) + CCTK_REAL, intent(out) :: ddgamma(3,3,3,3,3) + integer :: i,j,k,l,m,n + ! Gamma^i_jk,l = 1/2 g^im,l (g_mj,k + g_mk,j - g_jk,m) + ! + 1/2 g^im (g_mj,kl + g_mk,jl - g_jk,ml) + do i=1,3 + do j=1,3 + do k=1,3 + do l=1,3 + do m=1,3 + ddgamma(i,j,k,l,m) = 0 + do n=1,3 + ddgamma(i,j,k,l,m) = ddgamma(i,j,k,l,m) & + + 0.5d0 * ddgu(i,n,l,m) * (dgg(n,j,k) + dgg(n,k,j) - dgg(j,k,n)) & + + 0.5d0 * dgu(i,n,l) * (ddgg(n,j,k,m) + ddgg(n,k,j,m) - ddgg(j,k,n,m)) & + + 0.5d0 * dgu(i,n,m) * (ddgg(n,j,k,l) + ddgg(n,k,j,l) - ddgg(j,k,n,l)) & + + 0.5d0 * gu(i,n) * (dddgg(n,j,k,l,m) + dddgg(n,k,j,l,m) - dddgg(j,k,n,l,m)) + end do + end do + end do + end do + end do + end do + end subroutine calc_connectionderivs2 + subroutine calc_ricci (gamma, dgamma, ri) CCTK_REAL, intent(in) :: gamma(3,3,3), dgamma(3,3,3,3) CCTK_REAL, intent(out) :: ri(3,3) @@ -82,4 +110,42 @@ contains if (sum > 1.0e-12 * cnt) call CCTK_WARN (0, "Ricci tensor is not symmetric") end subroutine calc_ricci + subroutine calc_ricciderivs (gamma, dgamma, ddgamma, dri) + CCTK_REAL, intent(in) :: gamma(3,3,3), dgamma(3,3,3,3), ddgamma(3,3,3,3,3) + CCTK_REAL, intent(out) :: dri(3,3,3) + CCTK_REAL :: sum, cnt + integer :: i,j,k,l,m + ! R_ij = Gamma^k_ij,k - Gamma^k_ik,j + ! + Gamma^k_lk Gamma^l_ij - Gamma^k_lj Gamma^l_ki + do i=1,3 + do j=1,3 + do k=1,3 + dri(i,j,k) = 0 + do l=1,3 + dri(i,j,k) = dri(i,j,k) + ddgamma(l,i,j,l,k) & + & - ddgamma(l,i,l,j,k) + do m=1,3 + dri(i,j,k) = dri(i,j,k) + dgamma(l,m,l,k) * gamma(m,i,j) & + & + gamma(l,m,l) * dgamma(m,i,j,k) & + & - dgamma(l,m,j,k) * gamma(m,l,i) & + & - gamma(l,m,j) * dgamma(m,l,i,k) + end do + end do + end do + end do + end do + ! check symmetries + sum = 0 + cnt = 0 + do i=1,3 + do j=1,3 + do k=1,3 + sum = sum + (dri(i,j,k) - dri(j,i,k))**2 + cnt = cnt + dri(i,j,k)**2 + end do + end do + end do + if (sum > 1.0e-12 * cnt) call CCTK_WARN (0, "Ricci tensor is not symmetric") + end subroutine calc_ricciderivs + end module ricci -- cgit v1.2.3