From 57e38b49b43b1501449900c68c1790ab0e57dc65 Mon Sep 17 00:00:00 2001 From: diener Date: Fri, 13 Jun 2003 09:14:36 +0000 Subject: Changed EHFinder_MarkSurfaces to check not only the nearest neighbours when marking points, but also the nearest neighbours of the nearest neighbours. The reason being that connected surfaces can contain points where none of the nearest neighbours are inside the surface, but at least one of the next nearest neighbours are. git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/EHFinder/trunk@118 2a26948c-0e4f-0410-aee8-f1d3e353619c --- src/EHFinder_FindSurface.F90 | 116 ++++++++++++++++++++++++++++++++++--------- 1 file changed, 92 insertions(+), 24 deletions(-) diff --git a/src/EHFinder_FindSurface.F90 b/src/EHFinder_FindSurface.F90 index 221d969..95e0e63 100644 --- a/src/EHFinder_FindSurface.F90 +++ b/src/EHFinder_FindSurface.F90 @@ -792,6 +792,7 @@ subroutine EHFinder_MarkSurfaces(CCTK_ARGUMENTS) DECLARE_CCTK_FUNCTIONS CCTK_INT :: n_loc, n_glob, i, j, k + CCTK_INT :: i1, i2, j1, j2, k1, k2 logical :: marked ! Find index ranges for points excluding ghost zones and symmetry points for @@ -812,40 +813,107 @@ subroutine EHFinder_MarkSurfaces(CCTK_ARGUMENTS) if ( ( eh_mask(i,j,k) .ge. 0 ) .and. & ( f(i,j,k) .lt. 0 ) .and. & ( sc(i,j,k) .ne. surface_counter ) ) then - if ( .not. btest(eh_mask(i,j,k),0) ) then - if ( sc(i-1,j,k) .eq. surface_counter ) then - marked = .true. - end if + + ! Find out if we are on the boundary + if ( i + cctk_lbnd(1) .eq. 1 ) then + i1 = 0 + else + i1 = -1 end if - if ( .not. btest(eh_mask(i,j,k),1) ) then - if ( sc(i+1,j,k) .eq. surface_counter ) then - marked = .true. - end if + if ( i + cctk_lbnd(1) .eq. cctk_gsh(1) ) then + i2 = 0 + else + i2 = 1 end if - if ( .not. btest(eh_mask(i,j,k),2) ) then - if ( sc(i,j-1,k) .eq. surface_counter ) then - marked = .true. - end if + if ( j + cctk_lbnd(2) .eq. 1 ) then + j1 = 0 + else + j1 = -1 end if - if ( .not. btest(eh_mask(i,j,k),3) ) then - if ( sc(i,j+1,k) .eq. surface_counter ) then - marked = .true. - end if + if ( j + cctk_lbnd(2) .eq. cctk_gsh(2) ) then + j2 = 0 + else + j2 = 1 end if - if ( .not. btest(eh_mask(i,j,k),4) ) then - if ( sc(i,j,k-1) .eq. surface_counter ) then - marked = .true. - end if + if ( k + cctk_lbnd(3) .eq. 1 ) then + k1 = 0 + else + k1 = -1 end if - if ( .not. btest(eh_mask(i,j,k),5) ) then - if ( sc(i,j,k+1) .eq. surface_counter ) then - marked = .true. - end if + if ( k + cctk_lbnd(3) .eq. cctk_gsh(3) ) then + k2 = 0 + else + k2 = 1 end if + +! Then check if any valid neighbours are marked + if ( any ( sc(i+i1:i+i2,j+j1:j+j2,k) .eq. surface_counter ) ) then + marked = .true. + end if + if ( any ( sc(i+i1:i+i2,j,k+k1:k+k2) .eq. surface_counter ) ) then + marked = .true. + end if + if ( any ( sc(i,j+j1:j+j2,k+k1:k+k2) .eq. surface_counter ) ) then + marked = .true. + end if + +! if ( any ( sc(i+i1:i+i2,j+j1:j+j2,k+k1:k+k2) & +! .eq. surface_counter ) ) then +! marked = .true. +! end if + +! ! First look in the -x direction. +! if ( .not. btest(eh_mask(i,j,k),0) ) then +! if ( sc(i-1,j,k) .eq. surface_counter ) then +! marked = .true. +! end if +! if ( .not. btest(eh_mask(i-1,j,k),2) ) then +! if ( sc(i-1,j-1,k) .eq. surface_counter ) then +! marked = .true. +! end if +! if ( .not. btest(eh_mask(i-1,j,k),3) ) then +! if ( sc(i-1,j+1,k) .eq. surface_counter ) then +! marked = .true. +! end if +! if ( .not. btest(eh_mask(i-1,j,k),4) ) then +! if ( sc(i-1,j,k-1) .eq. surface_counter ) then +! marked = .true. +! end if +! if ( .not. btest(eh_mask(i-1,j,k),5) ) then +! if ( sc(i-1,j,k+1) .eq. surface_counter ) then +! marked = .true. +! end if +! end if +! if ( .not. btest(eh_mask(i,j,k),1) ) then +! if ( sc(i+1,j,k) .eq. surface_counter ) then +! marked = .true. +! end if +! end if +! if ( .not. btest(eh_mask(i,j,k),2) ) then +! if ( sc(i,j-1,k) .eq. surface_counter ) then +! marked = .true. +! end if +! end if +! if ( .not. btest(eh_mask(i,j,k),3) ) then +! if ( sc(i,j+1,k) .eq. surface_counter ) then +! marked = .true. +! end if +! end if +! if ( .not. btest(eh_mask(i,j,k),4) ) then +! if ( sc(i,j,k-1) .eq. surface_counter ) then +! marked = .true. +! end if +! end if +! if ( .not. btest(eh_mask(i,j,k),5) ) then +! if ( sc(i,j,k+1) .eq. surface_counter ) then +! marked = .true. +! end if +! end if if ( marked ) then sc(i,j,k) = surface_counter n_loc = n_loc + 1 end if + end if end do end do -- cgit v1.2.3