[svn:einsteintoolkit] NullSHRExtract/branches/tapir/ (Rev. 20) - Commits

7 Nov 2013

User: szilagyi
Date: 2013/11/07 11:52 AM
Modified:
 /branches/tapir/
  interface.ccl, schedule.ccl
 /branches/tapir/src/
  NullSHRE_ExtractWT.F90, NullSHRE_InitialWT.F90, NullSHRE_Pointing.F90, NullSHRE_modInverse.F90, NullSHRE_modNullDir.F90, NullSHRE_modRadius.F90, NullSHRE_modVars.F90
 /branches/tapir/src/include/
  NullSHRE_PointingNull.h
Log:
 work on an alterantive way of computing the bondi r
File Changes:
Directory: /branches/tapir/src/include/
=======================================
File [modified]: NullSHRE_PointingNull.h
Delta lines: +5 -2
===================================================================

--- branches/tapir/src/include/NullSHRE_PointingNull.h	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/include/NullSHRE_PointingNull.h	2013-11-07 17:52:13 UTC (rev 20)
@@ -77,9 +77,12 @@
     call my_nullify2(1,4,1,4, g1_n)
     call my_nullify2(1,4,1,4, g1_s)
-    call my_nullify2(1,4,1,4, gup_n)
-    call my_nullify2(1,4,1,4, gup_s)
+    call my_nullify2(1,3,1,3, g3up_n)
+    call my_nullify2(1,3,1,3, g3up_s)
+    call my_nullify2(1,4,1,4, g4up_n)
+    call my_nullify2(1,4,1,4, g4up_s)
+
     call my_nullify2(1,4,1,4, j0_n)
     call my_nullify2(1,4,1,4, j0_s)
Directory: /branches/tapir/src/
===============================
File [modified]: NullSHRE_ExtractWT.F90
Delta lines: +10 -10
===================================================================
--- branches/tapir/src/NullSHRE_ExtractWT.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_ExtractWT.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -68,8 +68,8 @@
   !---------------------------------------------------------------------------
   ! upstairs 3-metric
   !---------------------------------------------------------------------------
-  call wt_gup3 (g_n, detg_n, gup_n)
-  call wt_gup3 (g_s, detg_s, gup_s)
+  call wt_gup3 (g_n(1:3,1:3), detg_n, g3up_n)
+  call wt_gup3 (g_s(1:3,1:3), detg_s, g3up_s)
!---------------------------------------------------------------------------
   ! normals, null direction, geodesic expansion
@@ -78,8 +78,8 @@
   call wt_na (alpha_n, beta_n, na_n)
   call wt_na (alpha_s, beta_s, na_s)
-  call wt_sa (gup_n, sigma_n, sigma2_n, sa_n)
-  call wt_sa (gup_s, sigma_s, sigma2_s, sa_s)
+  call wt_sa (g3up_n, sigma_n, sigma2_n, sa_n)
+  call wt_sa (g3up_s, sigma_s, sigma2_s, sa_s)
elld_min = elld_min_coef * cctk_delta_time**elld_min_pow
@@ -102,11 +102,11 @@
   call wt_dna (alpha_n, beta_n, dalpha_n, dbeta_n, dna_n)
   call wt_dna (alpha_s, beta_s, dalpha_s, dbeta_s, dna_s)
-  call wt_dsa (sigma2_n, dsigma_n, sa_n, gup_n, dg_n, dsa_n)
-  call wt_dsa (sigma2_s, dsigma_s, sa_s, gup_s, dg_s, dsa_s)
+  call wt_dsa (sigma2_n, dsigma_n, sa_n, g3up_n, dg_n, dsa_n)
+  call wt_dsa (sigma2_s, dsigma_s, sa_s, g3up_s, dg_s, dsa_s)
-  call wt_gup4 (alpha_n, na_n, gup_n)
-  call wt_gup4 (alpha_s, na_s, gup_s)
+  call wt_gup4 (alpha_n, na_n, g3up_n, g4up_n)
+  call wt_gup4 (alpha_s, na_s, g3up_s, g4up_s)
!-------------------------------------------------------------------------
   ! derivatives with respect to the affine parameter lambda
@@ -115,9 +115,9 @@
   call wt_g1 (dg_n, ell_n, j0inv_n, g1_n)
   call wt_g1 (dg_s, ell_s, j0inv_s, g1_s)
-  call wt_j1 (g_n, gup_n, dg_n, dalpha_n, beta_n, dbeta_n, na_n,&
+  call wt_j1 (g_n, g4up_n, dg_n, dalpha_n, beta_n, dbeta_n, na_n,&
        sa_n, ell_n, dsa_n, dna_n, j0inv_n, elld_n, delld_n, j1_n)
-  call wt_j1 (g_s, gup_s, dg_s, dalpha_s, beta_s, dbeta_s, na_s,&
+  call wt_j1 (g_s, g4up_s, dg_s, dalpha_s, beta_s, dbeta_s, na_s,&
        sa_s, ell_s, dsa_s, dna_s, j0inv_s, elld_s, delld_s, j1_s)
!---------------------------------------------------------------------------
File [modified]: NullSHRE_InitialWT.F90
Delta lines: +10 -10
===================================================================
--- branches/tapir/src/NullSHRE_InitialWT.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_InitialWT.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -68,8 +68,8 @@
   !---------------------------------------------------------------------------
   ! upstairs 3-metric
   !---------------------------------------------------------------------------
-  call wt_gup3 (g_n, detg_n, gup_n)
-  call wt_gup3 (g_s, detg_s, gup_s)
+  call wt_gup3 (g_n(1:3,1:3), detg_n, g3up_n)
+  call wt_gup3 (g_s(1:3,1:3), detg_s, g3up_s)
!---------------------------------------------------------------------------
   ! normals, null direction, geodesic expansion
@@ -78,8 +78,8 @@
   call wt_na (alpha_n, beta_n, na_n)
   call wt_na (alpha_s, beta_s, na_s)
-  call wt_sa (gup_n, sigma_n, sigma2_n, sa_n)
-  call wt_sa (gup_s, sigma_s, sigma2_s, sa_s)
+  call wt_sa (g3up_n, sigma_n, sigma2_n, sa_n)
+  call wt_sa (g3up_s, sigma_s, sigma2_s, sa_s)
elld_min = elld_min_coef * cctk_delta_time**elld_min_pow
@@ -102,11 +102,11 @@
   call wt_dna (alpha_n, beta_n, dalpha_n, dbeta_n, dna_n)
   call wt_dna (alpha_s, beta_s, dalpha_s, dbeta_s, dna_s)
-  call wt_dsa (sigma2_n, dsigma_n, sa_n, gup_n, dg_n, dsa_n)
-  call wt_dsa (sigma2_s, dsigma_s, sa_s, gup_s, dg_s, dsa_s)
+  call wt_dsa (sigma2_n, dsigma_n, sa_n, g3up_n, dg_n, dsa_n)
+  call wt_dsa (sigma2_s, dsigma_s, sa_s, g3up_s, dg_s, dsa_s)
-  call wt_gup4 (alpha_n, na_n, gup_n)
-  call wt_gup4 (alpha_s, na_s, gup_s)
+  call wt_gup4 (alpha_n, na_n, g3up_n, g4up_n)
+  call wt_gup4 (alpha_s, na_s, g3up_s, g4up_s)
!-------------------------------------------------------------------------
   ! derivatives with respect to the affine parameter lambda
@@ -115,9 +115,9 @@
   call wt_g1 (dg_n, ell_n, j0inv_n, g1_n)
   call wt_g1 (dg_s, ell_s, j0inv_s, g1_s)
-  call wt_j1 (g_n, gup_n, dg_n, dalpha_n, beta_n, dbeta_n, na_n,&
+  call wt_j1 (g_n, g4up_n, dg_n, dalpha_n, beta_n, dbeta_n, na_n,&
        sa_n, ell_n, dsa_n, dna_n, j0inv_n, elld_n, delld_n, j1_n)
-  call wt_j1 (g_s, gup_s, dg_s, dalpha_s, beta_s, dbeta_s, na_s,&
+  call wt_j1 (g_s, g4up_s, dg_s, dalpha_s, beta_s, dbeta_s, na_s,&
        sa_s, ell_s, dsa_s, dna_s, j0inv_s, elld_s, delld_s, j1_s)
!---------------------------------------------------------------------------
File [modified]: NullSHRE_Pointing.F90
Delta lines: +3 -2
===================================================================
--- branches/tapir/src/NullSHRE_Pointing.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_Pointing.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -18,7 +18,7 @@
DECLARE_CCTK_ARGUMENTS
   DECLARE_CCTK_PARAMETERS
-  CCTK_INT :: i, j
+! CCTK_INT :: i, j
logical, save :: first_call = .true.
   if (first_call) call nullify_all
@@ -81,7 +81,8 @@
call point_v2(null_lsh, 12, 1, 3, 2, 3, dsigma_n, dsigma_s, WT_dsigma_vect)
-  call point_v2sym(null_lsh, 20, 1, 4, gup_n, gup_s, WT_gup_vect)
+  call point_v2sym(null_lsh, 12, 1, 3, g3up_n, g3up_s, WT_g3up_vect)
+  call point_v2sym(null_lsh, 20, 1, 4, g4up_n, g4up_s, WT_g4up_vect)
   call point_v2sym(null_lsh, 20, 1, 4, g1_n, g1_s, WT_g1_vect)
call point_v2(null_lsh, 24, 1, 4, 2, 4, dsa_n, dsa_s, WT_DSnorm_vect)
File [modified]: NullSHRE_modInverse.F90
Delta lines: +20 -19
===================================================================
--- branches/tapir/src/NullSHRE_modInverse.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_modInverse.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -5,51 +5,52 @@
contains
-  subroutine wt_gup3(g, detg, gup)
+  subroutine wt_gup3(g3, detg, g3up)
use NullSHRE_modGFdef
     implicit none
-    type (gf2d), dimension (4,4), intent (in)    :: g
+    type (gf2d), dimension (3,3), intent (in)    :: g3
     type (gf2d),                  intent (inout) :: detg
-    type (gf2d), dimension (4,4), intent (inout) :: gup
+    type (gf2d), dimension (3,3), intent (inout) :: g3up
-    detg%d = g(1,1)%d * g(2,2)%d * g(3,3)%d &
-           - g(1,1)%d * g(2,3)%d * g(2,3)%d &
-           - g(1,2)%d * g(1,2)%d * g(3,3)%d &
-           + g(1,2)%d * g(1,3)%d * g(2,3)%d &
-           + g(1,3)%d * g(1,2)%d * g(2,3)%d &
-           - g(1,3)%d * g(1,3)%d * g(2,2)%d
+    detg%d = g3(1,1)%d * g3(2,2)%d * g3(3,3)%d &
+           - g3(1,1)%d * g3(2,3)%d * g3(2,3)%d &
+           - g3(1,2)%d * g3(1,2)%d * g3(3,3)%d &
+           + g3(1,2)%d * g3(1,3)%d * g3(2,3)%d &
+           + g3(1,3)%d * g3(1,2)%d * g3(2,3)%d &
+           - g3(1,3)%d * g3(1,3)%d * g3(2,2)%d
-    gup(1,1)%d = (  g(2,2)%d * g(3,3)%d - g(2,3)%d * g(2,3)%d ) / detg%d
-    gup(1,2)%d = (- g(1,2)%d * g(3,3)%d + g(1,3)%d * g(2,3)%d ) / detg%d
-    gup(1,3)%d = (  g(1,2)%d * g(2,3)%d - g(1,3)%d * g(2,2)%d ) / detg%d
-    gup(2,2)%d = (  g(1,1)%d * g(3,3)%d - g(1,3)%d * g(1,3)%d ) / detg%d
-    gup(2,3)%d = (- g(1,1)%d * g(2,3)%d + g(1,3)%d * g(1,2)%d ) / detg%d
-    gup(3,3)%d = (  g(1,1)%d * g(2,2)%d - g(1,2)%d * g(1,2)%d ) / detg%d
+    g3up(1,1)%d = (  g3(2,2)%d * g3(3,3)%d - g3(2,3)%d * g3(2,3)%d ) / detg%d
+    g3up(1,2)%d = (- g3(1,2)%d * g3(3,3)%d + g3(1,3)%d * g3(2,3)%d ) / detg%d
+    g3up(1,3)%d = (  g3(1,2)%d * g3(2,3)%d - g3(1,3)%d * g3(2,2)%d ) / detg%d
+    g3up(2,2)%d = (  g3(1,1)%d * g3(3,3)%d - g3(1,3)%d * g3(1,3)%d ) / detg%d
+    g3up(2,3)%d = (- g3(1,1)%d * g3(2,3)%d + g3(1,3)%d * g3(1,2)%d ) / detg%d
+    g3up(3,3)%d = (  g3(1,1)%d * g3(2,2)%d - g3(1,2)%d * g3(1,2)%d ) / detg%d
end subroutine wt_gup3
-  subroutine wt_gup4(alpha, na, gup)
+  subroutine wt_gup4(alpha, na, g3up, g4up)
use NullSHRE_modGFdef
     implicit none
type (gf2d),                  intent (in)    :: alpha
     type (gf2d), dimension (4),   intent (in)    :: na
-    type (gf2d), dimension (4,4), intent (inout) :: gup
+    type (gf2d), dimension (3,3), intent (in)    :: g3up
+    type (gf2d), dimension (4,4), intent (inout) :: g4up
integer :: i, j
do i = 1, 3
        do j = i, 3
-          gup(i,j)%d = gup(i,j)%d - na(i)%d * na(j)%d
+          g4up(i,j)%d = g3up(i,j)%d - na(i)%d * na(j)%d
        end do
     end do
do i = 1, 4
-       gup(i,4)%d = - na(i)%d / alpha%d
+       g4up(i,4)%d = - na(i)%d / alpha%d
     end do
end subroutine wt_gup4
File [modified]: NullSHRE_modNullDir.F90
Delta lines: +8 -8
===================================================================
--- branches/tapir/src/NullSHRE_modNullDir.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_modNullDir.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -26,12 +26,12 @@
   end subroutine wt_na
-  subroutine wt_sa(gup, sigma, sigma2, sa) ! eq. (25)
+  subroutine wt_sa(g3up, sigma, sigma2, sa) ! eq. (25)
use NullSHRE_modGFdef
     implicit none
-    type (gf2d), dimension (4,4), intent (in)    :: gup
+    type (gf2d), dimension (3,3), intent (in)    :: g3up
     type (gf2d), dimension (3),   intent (in)    :: sigma
     type (gf2d),                  intent (inout) :: sigma2
     type (gf2d), dimension (4),   intent (inout) :: sa
@@ -43,7 +43,7 @@
     do i = 1, 3
        sa(i)%d = 0.d0
        do j = 1, 3
-          sa(i)%d = sa(i)%d + gup(i,j)%d * sigma(j)%d
+          sa(i)%d = sa(i)%d + g3up(i,j)%d * sigma(j)%d
        end do
     end do
@@ -176,7 +176,7 @@
   end subroutine wt_dna
-  subroutine wt_dsa (sigma2, dsigma, sa, gup, dg, dsa)  ! eq (41..43)
+  subroutine wt_dsa (sigma2, dsigma, sa, g3up, dg, dsa)  ! eq (41..43)
use NullSHRE_modGFdef
     implicit none
@@ -184,7 +184,7 @@
     type (gf2d),                    intent (in)    :: sigma2
     type (gf2d), dimension (3,2:3), intent (in)    :: dsigma
     type (gf2d), dimension (4),     intent (in)    :: sa
-    type (gf2d), dimension (4,4),   intent (in)    :: gup
+    type (gf2d), dimension (3,3),   intent (in)    :: g3up
     type (gf2d), dimension (4,4,4), intent (in)    :: dg
     type (gf2d), dimension (4,2:4), intent (inout) :: dsa
@@ -197,12 +197,12 @@
           dsa(i,a)%d = 0.d0
           do k = 1, 3
              do l = 1, 3
-                dsa(i,a)%d = dsa(i,a)%d + (- gup(i,k)%d + 0.5d0 *&
+                dsa(i,a)%d = dsa(i,a)%d + (- g3up(i,k)%d + 0.5d0 *&
                         & sa(i)%d * sa(k)%d) * sa(l)%d * dg(k,l,a)%d
              end do
           end do
           do k = 1, 3
-             dsa(i,a)%d = dsa(i,a)%d + (gup(i,k)%d - sa(i)%d * sa(k)&
+             dsa(i,a)%d = dsa(i,a)%d + (g3up(i,k)%d - sa(i)%d * sa(k)&
                   & %d) * dsigma(k,a)%d / sqrt(sigma2%d)
           end do
        end do
@@ -215,7 +215,7 @@
        do k = 1, 3
           do l = 1, 3
              dsa(i,4)%d = dsa(i,4)%d + dg(k,l,4)%d * sa(l)%d * (&
-                  & -gup(i,k)%d + 0.5d0 * sa(i)%d * sa(k)%d)
+                  & -g3up(i,k)%d + 0.5d0 * sa(i)%d * sa(k)%d)
           end do
        end do
     end do
File [modified]: NullSHRE_modRadius.F90
Delta lines: +24 -0
===================================================================
--- branches/tapir/src/NullSHRE_modRadius.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_modRadius.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -5,6 +5,30 @@
contains
+  subroutine wt_r_long(cr, temp, detg, g3up, j0inv, r0)
+    use NullSHRE_modGFdef
+    implicit none
+    CCTK_REAL,                     intent (in)    :: cr
+    type (gf2d),                   intent (in)    :: detg
+    type (gf2d),                   intent (inout) :: temp
+    type (gf2d), dimension (3,3),  intent (in)    :: g3up
+    type (gf2d), dimension (4,4),  intent (in)    :: j0inv
+    type (gf2d),                   intent (inout) :: r0
+
+    integer :: i, j
+    ! temp stands for g^rr = j0inv^r_i j0inv^r_j g^ij
+
+    temp%d = 0
+    do i = 1, 3
+       do j = 1, 3
+       temp%d = temp%d + j0inv(1,i)%d * j0inv(1,j)%d * g3up(i,j)%d
+       end do
+    end do
+
+    r0%d = cr * sqrt(sqrt(detg%d*temp%d))
+
+  end subroutine wt_r_long
+
   subroutine wt_r(nn1, nn2, pp, eta0, r0)
 ! the bondi radius eq. (53) 
     use NullSHRE_modGFdef
File [modified]: NullSHRE_modVars.F90
Delta lines: +2 -1
===================================================================
--- branches/tapir/src/NullSHRE_modVars.F90	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/src/NullSHRE_modVars.F90	2013-11-07 17:52:13 UTC (rev 20)
@@ -23,7 +23,8 @@
type (gf2d), save, dimension (3)         :: sigma_n, sigma_s
    type (gf2d), save, dimension (4,4)       :: g1_n, g1_s
-   type (gf2d), save, dimension (4,4)       :: gup_n, gup_s
+   type (gf2d), save, dimension (4,4)       :: g4up_n, g4up_s
+   type (gf2d), save, dimension (3,3)       :: g3up_n, g3up_s
    type (gf2d), save, dimension (4,4)       :: j0_n, j0_s
    type (gf2d), save, dimension (4,4)       :: j0inv_n, j0inv_s
    type (gf2d), save, dimension (4,4)       :: j1_n, j1_s
Directory: /branches/tapir/
===========================
File [modified]: interface.ccl
Delta lines: +7 -1
===================================================================
--- branches/tapir/interface.ccl	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/interface.ccl	2013-11-07 17:52:13 UTC (rev 20)
@@ -244,12 +244,18 @@
    GHOSTSIZE=NullGrid::N_ang_ghost_pts,NullGrid::N_ang_ghost_pts \
  "lambda deriv. of covariant metric"
-CCTK_REAL WT_gup_vect[20] TYPE=ARRAY dim=2 size=\
+CCTK_REAL WT_g4up_vect[20] TYPE=ARRAY dim=2 size=\
 (NullGrid::N_ang_pts_inside_eq+2*(NullGrid::N_ang_ev_outside_eq+NullGrid::N_ang_stencil_size)),\
 (NullGrid::N_ang_pts_inside_eq+2*(NullGrid::N_ang_ev_outside_eq+NullGrid::N_ang_stencil_size)) \
    GHOSTSIZE=NullGrid::N_ang_ghost_pts,NullGrid::N_ang_ghost_pts \
 "contravariant metric"
+CCTK_REAL WT_g3up_vect[12] TYPE=ARRAY dim=2 size=\
+(NullGrid::N_ang_pts_inside_eq+2*(NullGrid::N_ang_ev_outside_eq+NullGrid::N_ang_stencil_size)),\
+(NullGrid::N_ang_pts_inside_eq+2*(NullGrid::N_ang_ev_outside_eq+NullGrid::N_ang_stencil_size)) \
+   GHOSTSIZE=NullGrid::N_ang_ghost_pts,NullGrid::N_ang_ghost_pts \
+"contravariant metric"
+
 CCTK_REAL WT_jac[32] TYPE=ARRAY dim=2 size=\
 (NullGrid::N_ang_pts_inside_eq+2*(NullGrid::N_ang_ev_outside_eq+NullGrid::N_ang_stencil_size)),\
 (NullGrid::N_ang_pts_inside_eq+2*(NullGrid::N_ang_ev_outside_eq+NullGrid::N_ang_stencil_size)) \
File [modified]: schedule.ccl
Delta lines: +2 -2
===================================================================
--- branches/tapir/schedule.ccl	2013-11-07 01:02:36 UTC (rev 19)
+++ branches/tapir/schedule.ccl	2013-11-07 17:52:13 UTC (rev 20)
@@ -17,7 +17,7 @@
STORAGE: WT_CmpxCoord
   STORAGE: WT_scalars, WT_3vectors, WT_4vectors, WT_d4vectors
-  STORAGE: WT_dsigma_vect, WT_g1_vect, WT_gup_vect
+  STORAGE: WT_dsigma_vect, WT_g1_vect, WT_g4up_vect, WT_g3up_vect
   STORAGE: WT_jac, WT_affine_metric, WT_deta0_vect
   STORAGE: x_wt[2], j_wt[2], j_l, beta_wt[2], beta_l
   STORAGE: q_wt[2], u_wt[2], u_l, u_x_wt, w_wt[2], w_l
@@ -67,7 +67,7 @@
    {
      STORAGE: WT_CmpxCoord
      STORAGE: WT_scalars, WT_3vectors, WT_4vectors, WT_d4vectors
-     STORAGE: WT_dsigma_vect, WT_g1_vect, WT_gup_vect
+     STORAGE: WT_dsigma_vect, WT_g1_vect, WT_g4up_vect, WT_g3up_vect
      STORAGE: WT_jac, WT_affine_metric, WT_deta0_vect
      STORAGE: x_wt[2], j_wt[2], j_l, beta_wt[2], beta_l
      STORAGE: u_wt[2], u_l, u_x_wt, w_wt[2], w_l, q_wt[2]

    