From 60fbb40f73b2b7393881cd557aa7deffc59120fe Mon Sep 17 00:00:00 2001 From: miguel Date: Thu, 29 Mar 2001 16:00:46 +0000 Subject: Removing obsolete routine. WE are quick, aren't we? git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinEvolve/LegoExcision/trunk@5 f75ba9e5-694f-0410-ac2c-87ea7ce7132b --- src/excision.F90 | 293 ------------------------------------------------------- 1 file changed, 293 deletions(-) delete mode 100644 src/excision.F90 diff --git a/src/excision.F90 b/src/excision.F90 deleted file mode 100644 index 2b7a670..0000000 --- a/src/excision.F90 +++ /dev/null @@ -1,293 +0,0 @@ -! $Header$ - -#include "cctk.h" -#include "cctk_Arguments.h" -#include "cctk_Parameters.h" - -subroutine excision_boundary (ierr, cctkgh, var, oldvar, mask) - - implicit none - - DECLARE_CCTK_FUNCTIONS - DECLARE_CCTK_PARAMETERS - - ! Arguments - integer, intent(out) :: ierr - CCTK_POINTER, intent(in) :: cctkgh - CCTK_REAL, intent(inout) :: var(:,:,:) - CCTK_REAL, intent(in) :: oldvar(:,:,:) - CCTK_REAL, intent(in) :: mask(:,:,:) - - ! distinguish between faces, edges, and corners. - - ! look at the excised points, and distinguish the following cases, - ! in this order: - - ! [trivial case] - ! 1: if (none) - ! case 0: interior - ! [box shape cases] - ! 2: if (only one corner) - ! case I: corner - ! 3: if (only points on one edge) - ! case II: edge - ! 4: if (only points on one face) - ! case III: face - ! [convex cases] - ! 5: if (only points on two adjacent faces) - ! case IIa: anti-edge - ! 6: if (one points on three adjacent faces) - ! case Ia: anti-corner - ! [illegal cases] - ! error - - logical mask_interior (-1:+1,-1:+1,-1:+1, 1) - logical mask_corner (-1:+1,-1:+1,-1:+1, 8) - logical mask_edge (-1:+1,-1:+1,-1:+1, 12) - logical mask_face (-1:+1,-1:+1,-1:+1, 6) - logical mask_antiedge (-1:+1,-1:+1,-1:+1, 12) - logical mask_anticorner(-1:+1,-1:+1,-1:+1, 8) - integer dir_interior (3, 1), num_interior - integer dir_corner (3, 8), num_corner - integer dir_edge (3, 12), num_edge - integer dir_face (3, 6), num_face - integer dir_antiedge (3, 12), num_antiedge - integer dir_anticorner (3, 8), num_anticorner - - logical msk(-1:+1,-1:+1,-1:+1) - - integer di,dj,dk - integer i,j,k - integer n - integer dir(3) - - character(1000) message - - - - ! - ! Check the mask - ! - do k=1,size(mask,3) - do j=1,size(mask,2) - do i=1,size(mask,1) - if (mask(i,j,k)/=0 .and. mask(i,j,k)/=1) then - write (message, '("Mask contains illegal value at (",i4,",",i4,",",i4,")")') i,j,k - call CCTK_WARN (0, trim(message)) - end if - end do - end do - end do - - - - ! - ! Create templates - ! - num_interior = 0 - num_corner = 0 - num_edge = 0 - num_face = 0 - num_antiedge = 0 - num_anticorner = 0 - - do dk=-1,+1 - do dj=-1,+1 - do di=-1,+1 - - select case (count( (/di,dj,dk/) /= 0 )) - - case (0) - ! interior (or error, which is not explicitely checked) - num_interior = num_interior + 1 - dir_interior(:,num_interior) = (/di,dj,dk/) - do k=-1,+1 - do j=-1,+1 - do i=-1,+1 - mask_interior(i,j,k,num_interior) = & - any( (/di,dj,dk/)/=0 & - .and. (/i,j,k/)==(/di,dj,dk/) ) - end do - end do - end do - - case (1) - ! face (or anti-face, which is the same thing) - num_face = num_face + 1 - dir_face(:,num_face) = (/di,dj,dk/) - do k=-1,+1 - do j=-1,+1 - do i=-1,+1 - mask_face(i,j,k,num_face) = & - any( (/di,dj,dk/)/=0 & - .and. (/i,j,k/)==(/di,dj,dk/) ) - end do - end do - end do - - case (2) - ! edge or anti-edge - num_edge = num_edge + 1 - dir_edge(:,num_edge) = (/di,dj,dk/) - do k=-1,+1 - do j=-1,+1 - do i=-1,+1 - mask_edge(i,j,k,num_edge) = & - all( (/di,dj,dk/)==0 & - .or. (/i,j,k/)==(/di,dj,dk/) ) - end do - end do - end do - - num_antiedge = num_antiedge + 1 - dir_antiedge(:,num_antiedge) = (/di,dj,dk/) - do k=-1,+1 - do j=-1,+1 - do i=-1,+1 - mask_antiedge(i,j,k,num_antiedge) = & - any( (/di,dj,dk/)/=0 & - .and. (/i,j,k/)==(/di,dj,dk/) ) - end do - end do - end do - - case (3) - ! corner or anti-corner - num_corner = num_corner + 1 - dir_corner(:,num_corner) = (/di,dj,dk/) - do k=-1,+1 - do j=-1,+1 - do i=-1,+1 - mask_corner(i,j,k,num_corner) = & - all( (/di,dj,dk/)==0 & - .or. (/i,j,k/)==(/di,dj,dk/) ) - end do - end do - end do - - num_anticorner = num_anticorner + 1 - dir_anticorner(:,num_anticorner) = (/di,dj,dk/) - do k=-1,+1 - do j=-1,+1 - do i=-1,+1 - mask_anticorner(i,j,k,num_anticorner) = & - any( (/di,dj,dk/)/=0 & - .and. (/i,j,k/)==(/di,dj,dk/) ) - end do - end do - end do - - case default - call CCTK_WARN (0, "internal error") - end select - - end do - end do - end do - - ! sanity check - if ( num_interior /= 1 & - .or. num_corner /= 8 & - .or. num_edge /= 12 & - .or. num_face /= 6 & - .or. num_antiedge /= 12 & - .or. num_anticorner /= 8 ) then - call CCTK_WARN (0, "internal error") - end if - - - - ! - ! Calculate normal directions - ! Extrapolate - ! - do k=2,size(mask,3)-1 - do j=2,size(mask,2)-1 - do i=2,size(mask,1)-1 - - if (mask(i,j,k)==0 .and. any(mask(i-1:i+1,j-1:j+1,k-1:k+1)/=0)) then - - ! - ! Calculate normal direction - ! - - ! make a snapshot - msk = mask(i-1:i+1,j-1:j+1,k-1:k+1) /= 0 - - ! for a template to match, all grid points have to be - ! either active, or set in the template. - search: do - - ! 1: case 0: interior - do n=1,num_interior - if (all(msk .or. mask_interior(:,:,:,n))) then - dir = dir_interior(:,n) - exit search - end if - end do - - ! 2: case I: corner - do n=1,num_corner - if (all(msk .or. mask_corner(:,:,:,n))) then - dir = dir_corner(:,n) - exit search - end if - end do - - ! 3: case II: edge - do n=1,num_edge - if (all(msk .or. mask_edge(:,:,:,n))) then - dir = dir_edge(:,n) - exit search - end if - end do - - ! 4: case III: face - do n=1,num_face - if (all(msk .or. mask_face(:,:,:,n))) then - dir = dir_face(:,n) - exit search - end if - end do - - ! 5: case IIa: anti-edge - do n=1,num_antiedge - if (all(msk .or. mask_antiedge(:,:,:,n))) then - dir = dir_antiedge(:,n) - exit search - end if - end do - - ! 6: case Ia: anticorner - do n=1,num_anticorner - if (all(msk .or. mask_anticorner(:,:,:,n))) then - dir = dir_anticorner(:,n) - exit search - end if - end do - - ! error: unclassified - call CCTK_WARN (0, "Mask is not convex") - - end do search - - dir(:) = -dir(:) - - - - ! - ! Extrapolate - ! - var(i,j,k) = oldvar(i,j,k) + var(i+dir(1), j+dir(2), k+dir(3)) & - - oldvar(i+dir(1), j+dir(2), k+dir(3)) - - end if - - end do - end do - end do - - ! Success - ierr = 0 - -end subroutine excision_boundary -- cgit v1.2.3