second iteration of constraint transport

  * fix some indices
  * move poison loop to proper spot
  * compute conserved divergence in track_divB
  * add Maple worksheet to check constraint transport indices


git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinEvolve/GRHydro/trunk@298 c83d129a-5a75-4d5a-9c4d-ed3a5855bf45

1 files changed, 97 insertions, 3 deletions

diff --git a/src/GRHydro_SourceM.F90 b/src/GRHydro_SourceM.F90
index 1456500..cdf9ef0 100644
--- a/src/GRHydro_SourceM.F90
+++ b/src/GRHydro_SourceM.F90
@@ -113,6 +113,9 @@ subroutine SourceTermsM(CCTK_ARGUMENTS)
   Bconsrhs=0.d0
 
   if (evolve_tracer .ne. 0) then
+     Bconsrhsx(i,j,k) = 0.d0
+     Bconsrhsy(i,j,k) = 0.d0
+     Bconsrhsz(i,j,k) = 0.d0
     cons_tracerrhs = 0.d0
   end if
 
@@ -130,9 +133,6 @@ subroutine SourceTermsM(CCTK_ARGUMENTS)
 
   if (transport_constraints .ne. 0) then
      Evec = 0.d0
-     Bconsrhsx(i,j,k) = 0.d0
-     Bconsrhsy(i,j,k) = 0.d0
-     Bconsrhsz(i,j,k) = 0.d0
   endif
   
 
@@ -495,6 +495,100 @@ subroutine SourceTermsM(CCTK_ARGUMENTS)
 
   deallocate(force_spatial_second_order)
 
+#if(0) /* poison edges of domain */
+  if(last_iteration_seen .ne. cctk_iteration .or. reflevel .ne. grhydro_reflevel) then
+    last_iteration_seen = cctk_iteration
+    reflevel = grhydro_reflevel
+    mol_substep = 0
+  else
+    mol_substep = mol_substep + 1
+  end if
+  do k = 1, GRHydro_stencil*mol_substep
+    do j = 1, cctk_lsh(2)
+      do i = 1, cctk_lsh(1)
+        dens(i,j,k) = -1d100
+        Scon(i,j,k,1) = -1d100
+        Scon(i,j,k,2) = -1d100
+        Scon(i,j,k,3) = -1d100
+        tau(i,j,k) = -1d100
+        Bcons(i,j,k,1) = -1d100
+        Bcons(i,j,k,2) = -1d100
+        Bcons(i,j,k,3) = -1d100
+      end do
+    end do
+  end do
+  do k = cctk_lsh(3)-GRHydro_stencil*mol_substep+1, cctk_lsh(3)
+    do j = 1, cctk_lsh(2)
+      do i = 1, cctk_lsh(1)
+        dens(i,j,k) = -1d100
+        Scon(i,j,k,1) = -1d100
+        Scon(i,j,k,2) = -1d100
+        Scon(i,j,k,3) = -1d100
+        tau(i,j,k) = -1d100
+        Bcons(i,j,k,1) = -1d100
+        Bcons(i,j,k,2) = -1d100
+        Bcons(i,j,k,3) = -1d100
+      end do
+    end do
+  end do
+  do i = 1, GRHydro_stencil*mol_substep
+    do k = 1, cctk_lsh(3)
+      do j = 1, cctk_lsh(2)
+        dens(i,j,k) = -1d100
+        Scon(i,j,k,1) = -1d100
+        Scon(i,j,k,2) = -1d100
+        Scon(i,j,k,3) = -1d100
+        tau(i,j,k) = -1d100
+        Bcons(i,j,k,1) = -1d100
+        Bcons(i,j,k,2) = -1d100
+        Bcons(i,j,k,3) = -1d100
+      end do
+    end do
+  end do
+  do i = cctk_lsh(1)-GRHydro_stencil*mol_substep+1, cctk_lsh(1)
+    do k = 1, cctk_lsh(3)
+      do j = 1, cctk_lsh(2)
+        dens(i,j,k) = -1d100
+        Scon(i,j,k,1) = -1d100
+        Scon(i,j,k,2) = -1d100
+        Scon(i,j,k,3) = -1d100
+        tau(i,j,k) = -1d100
+        Bcons(i,j,k,1) = -1d100
+        Bcons(i,j,k,2) = -1d100
+        Bcons(i,j,k,3) = -1d100
+      end do
+    end do
+  end do
+  do j = 1, GRHydro_stencil*mol_substep
+    do i = 1, cctk_lsh(1)
+      do k = 1, cctk_lsh(3)
+        dens(i,j,k) = -1d100
+        Scon(i,j,k,1) = -1d100
+        Scon(i,j,k,2) = -1d100
+        Scon(i,j,k,3) = -1d100
+        tau(i,j,k) = -1d100
+        Bcons(i,j,k,1) = -1d100
+        Bcons(i,j,k,2) = -1d100
+        Bcons(i,j,k,3) = -1d100
+      end do
+    end do
+  end do
+  do j = cctk_lsh(2)-GRHydro_stencil*mol_substep+1, cctk_lsh(2)
+    do i = 1, cctk_lsh(1)
+      do k = 1, cctk_lsh(3)
+        dens(i,j,k) = -1d100
+        Scon(i,j,k,1) = -1d100
+        Scon(i,j,k,2) = -1d100
+        Scon(i,j,k,3) = -1d100
+        tau(i,j,k) = -1d100
+        Bcons(i,j,k,1) = -1d100
+        Bcons(i,j,k,2) = -1d100
+        Bcons(i,j,k,3) = -1d100
+      end do
+    end do
+  end do
+#endif
+
 end subroutine SourceTermsM