From a18a0d30a72672b037362064cad51c5ff6f503fc Mon Sep 17 00:00:00 2001 From: diener Date: Tue, 30 Sep 2003 08:59:17 +0000 Subject: Removed EHFinder_Integrate2.F90 (now EHFinder_Integrate.F90). git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/EHFinder/trunk@142 2a26948c-0e4f-0410-aee8-f1d3e353619c --- src/EHFinder_Integrate2.F90 | 596 -------------------------------------------- 1 file changed, 596 deletions(-) delete mode 100644 src/EHFinder_Integrate2.F90 diff --git a/src/EHFinder_Integrate2.F90 b/src/EHFinder_Integrate2.F90 deleted file mode 100644 index d822050..0000000 --- a/src/EHFinder_Integrate2.F90 +++ /dev/null @@ -1,596 +0,0 @@ -! Integrating over the surface(s). Right now only in full mode. -! $Header$ - -#include "cctk.h" -#include "cctk_Parameters.h" -#include "cctk_Arguments.h" - -subroutine EHFinder_FindSurfaceElement(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - - CCTK_INT :: i, j, k, l, im, jm - CCTK_INT :: interp_handle, table_handle, coord_system_handle - - character(len=200) :: area_interp - CCTK_INT :: area_interp_len - character(len=7) :: area_order - - CCTK_INT, dimension(4) :: bbox - CCTK_INT, dimension(2) :: gsh, lsh, lbnd, ubnd, nghost - - CCTK_REAL :: dtheta, dphi, dthetainv, dphiinv - CCTK_REAL :: dxdth, dxdph, dydth, dydph, dzdth, dzdph -! CCTK_REAL :: gtt, gtp, gpp - CCTK_POINTER, dimension(3) :: interp_coords - CCTK_POINTER, dimension(7) :: out_array - CCTK_INT, dimension(7) :: in_array - CCTK_INT, dimension(7), parameter :: out_types = (/ CCTK_VARIABLE_REAL, & - CCTK_VARIABLE_REAL, & - CCTK_VARIABLE_REAL, & - CCTK_VARIABLE_REAL, & - CCTK_VARIABLE_REAL, & - CCTK_VARIABLE_REAL, & - CCTK_VARIABLE_REAL /) - -! If finding of surface failed do not try to integrate but exit. - if ( find_surface_status .lt. 0 ) then - return - endif - - l = levelset_counter - - call CCTK_INFO ( 'Finding surface element' ) - call CCTK_GroupbboxGN ( ierr, cctkGH, 4, bbox, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get bounding box for surface arrays" ) - end if - call CCTK_GroupgshGN ( ierr, cctkGH, 2, gsh, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get global size for surface arrays" ) - end if - call CCTK_GrouplbndGN ( ierr, cctkGH, 2, lbnd, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get lower bounds for surface arrays" ) - end if - call CCTK_GroupubndGN ( ierr, cctkGH, 2, ubnd, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get upper bounds for surface arrays" ) - end if - call CCTK_GrouplshGN ( ierr, cctkGH, 2, lsh, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for surface arrays" ) - end if - call CCTK_GroupnghostzonesGN ( ierr, cctkGH, 2, nghost, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for surface arrays" ) - end if - - dtheta = ctheta(2,1) - ctheta(1,1) - dphi = cphi(1,2) - cphi(1,1) - dthetainv = one / dtheta - dphiinv = one / dphi - - call CCTK_FortranString ( area_interp_len, area_interpolator, & - area_interp ) - - call CCTK_InterpHandle ( interp_handle, area_interp(1:area_interp_len) ) - - if ( interp_handle .lt. 0 ) then - call CCTK_WARN( 0, "Cannot get handle for interpolation. Forgot to activate an implementation providing interpolation operators??" ) - end if - - write(area_order,'(a6,i1)') 'order=',area_interpolation_order - - call Util_TableCreateFromString ( table_handle, area_order ) - if ( table_handle .lt. 0 ) then - call CCTK_WARN( 0, "Cannot create parameter table for interpolator" ) - end if - - call CCTK_CoordSystemHandle ( coord_system_handle, "cart3d" ) - if ( coord_system_handle .lt. 0) then - call CCTK_WARN( 0, "Cannot get handle for cart3d coordinate system. Forgot to activate an implementation providing coordinates ??" ) - endif - - interp_coords(1) = CCTK_PointerTo(interp_x) - interp_coords(2) = CCTK_PointerTo(interp_y) - interp_coords(3) = CCTK_PointerTo(interp_z) - - out_array(1) = CCTK_PointerTo(gxxi) - out_array(2) = CCTK_PointerTo(gxyi) - out_array(3) = CCTK_PointerTo(gxzi) - out_array(4) = CCTK_PointerTo(gyyi) - out_array(5) = CCTK_PointerTo(gyzi) - out_array(6) = CCTK_PointerTo(gzzi) - out_array(7) = CCTK_PointerTo(psii) - -! find the cartesian coordinates for the interpolation points - - do j = 1, lsh(2) - do i = 1, lsh(1) - interp_x(i,j) = center(1) + rsurf(i,j) * sintheta(i,j) * cosphi(i,j) - interp_y(i,j) = center(2) + rsurf(i,j) * sintheta(i,j) * sinphi(i,j) - interp_z(i,j) = center(3) + rsurf(i,j) * costheta(i,j) - end do - end do - - call CCTK_VarIndex ( in_array(1), "admbase::gxx" ) - call CCTK_VarIndex ( in_array(2), "admbase::gxy" ) - call CCTK_VarIndex ( in_array(3), "admbase::gxz" ) - call CCTK_VarIndex ( in_array(4), "admbase::gyy" ) - call CCTK_VarIndex ( in_array(5), "admbase::gyz" ) - call CCTK_VarIndex ( in_array(6), "admbase::gzz" ) - - if ( CCTK_EQUALS ( metric_type, "static conformal" ) ) then - call CCTK_VarIndex ( in_array(7), "staticconformal::psi" ) - - call CCTK_InterpGridArrays ( ierr, cctkGH, 3, interp_handle, & - table_handle, coord_system_handle, & - lsh(1) * lsh(2), CCTK_VARIABLE_REAL, & - interp_coords, 7, in_array, & - 7, out_types, out_array ) - gxxi = psii**4 * gxxi - gxyi = psii**4 * gxyi - gxzi = psii**4 * gxzi - gyyi = psii**4 * gyyi - gyzi = psii**4 * gyzi - gzzi = psii**4 * gzzi - else - call CCTK_InterpGridArrays ( ierr, cctkGH, 3, interp_handle, & - table_handle, coord_system_handle, & - lsh(1) * lsh(2), CCTK_VARIABLE_REAL, & - interp_coords, 6, in_array(1:6), & - 6, out_types(1:6), out_array(1:6) ) - end if - - do j = 1, lsh(2) - do i = 1, lsh(1) - dxdth = ( drdtheta(i,j) * sintheta(i,j) + & - rsurf(i,j) * costheta(i,j) ) * cosphi(i,j) - dydth = ( drdtheta(i,j) * sintheta(i,j) + & - rsurf(i,j) * costheta(i,j) ) * sinphi(i,j) - dzdth = drdtheta(i,j) * costheta(i,j) - rsurf(i,j) * sintheta(i,j) - dxdph = ( drdphi(i,j) * cosphi(i,j) - & - rsurf(i,j) * sinphi(i,j) ) * sintheta(i,j) - dydph = ( drdphi(i,j) * sinphi(i,j) + & - rsurf(i,j) * cosphi(i,j) ) * sintheta(i,j) - dzdph = drdphi(i,j) * costheta(i,j) - - gtt(i,j) = dxdth**2 * gxxi(i,j) + dydth**2 * gyyi(i,j) + & - dzdth**2 * gzzi(i,j) + & - two * ( dxdth * dydth * gxyi(i,j) + & - dxdth * dzdth * gxzi(i,j) + & - dydth * dzdth * gyzi(i,j) ) - gtp(i,j) = dxdth * dxdph * gxxi(i,j) + & - ( dxdth * dydph + dydth * dxdph ) * gxyi(i,j) + & - ( dxdth * dzdph + dzdth * dxdph ) * gxzi(i,j) + & - dydth * dydph * gyyi(i,j) + & - ( dydth * dzdph + dzdth * dydph ) * gyzi(i,j) + & - dzdth * dzdph * gzzi(i,j) - gpp(i,j) = dxdph**2 * gxxi(i,j) + dydph**2 * gyyi(i,j) + & - dzdph**2 * gzzi(i,j) + & - two * ( dxdph * dydph * gxyi(i,j) + & - dxdph * dzdph * gxzi(i,j) + & - dydph * dzdph * gyzi(i,j) ) - da(i,j) = dtheta * dphi * sqrt ( gtt(i,j) * gpp(i,j) - gtp(i,j)**2 ) - dltheta(i,j) = dtheta * sqrt ( gtt(i,j) ) - dlphi(i,j) = dphi * sqrt ( gpp(i,j) ) - end do - end do - -end subroutine EHFinder_FindSurfaceElement - - -subroutine EHFinder_IntegrateArea(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - - CCTK_INT :: int_var, sn, l - CCTK_REAL :: eh_area_tmp - character(len=30) :: info_message - CCTK_INT, dimension(1) :: lbnd, ubnd, lsh - -! If finding of surface failed do not try to integrate but exit. - if ( find_surface_status .lt. 0 ) then - return - endif - - l = levelset_counter - - call CCTK_INFO ( 'Integrating area' ) - - sn = surface_counter - integrate_counter - int_tmp = sym_factor * weights * da - - call CCTK_VarIndex ( int_var, "ehfinder::int_tmp" ) - if ( int_var .lt. 0 ) then - call CCTK_WARN ( 0, 'Could not get index to grid array int_tmp' ) - end if - - call CCTK_Reduce ( ierr, cctkGH, -1, sum_handle, 1, CCTK_VARIABLE_REAL, & - eh_area_tmp, 1, int_var ) - write(info_message,'(a15,f10.5)') 'Horizon area = ',eh_area_tmp - - call CCTK_INFO(info_message) - - call CCTK_GrouplbndGN ( ierr, cctkGH, 1, lbnd, "ehfinder::eh_area2" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get lower bounds for area array" ) - end if - call CCTK_GroupubndGN ( ierr, cctkGH, 1, ubnd, "ehfinder::eh_area2" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get upper bounds for area array" ) - end if - call CCTK_GrouplshGN ( ierr, cctkGH, 1, lsh, "ehfinder::eh_area2" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for area array" ) - end if - - if ( ( sn .ge. lbnd(1) + 1 ) .and. ( sn .le. ubnd(1) + 1 ) ) then - eh_area2(sn-lbnd(1),l) = eh_area_tmp - end if - -end subroutine EHFinder_IntegrateArea - - -subroutine EHFinder_IntegrateCentroid(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - - CCTK_INT :: int_var, sn, k, l - character(len=50) :: info_message - CCTK_INT, dimension(1) :: lbnd, ubnd, lsh - CCTK_REAL :: centroid_x, centroid_y, centroid_z - -! If finding of surface failed do not try to integrate but exit. - if ( find_surface_status .lt. 0 ) then - return - endif - - l = levelset_counter - - call CCTK_INFO ( 'Integrating centroid' ) - - sn = surface_counter - integrate_counter - - interp_x = center(1) + rsurf * sintheta * cosphi - interp_y = center(2) + rsurf * sintheta * sinphi - interp_z = center(3) + rsurf * costheta - - call CCTK_VarIndex ( int_var, "ehfinder::int_tmp" ) - if ( int_var .lt. 0 ) then - call CCTK_WARN ( 0, 'Could not get index to grid array int_tmp' ) - end if - - int_tmp = sym_factor * interp_x * weights * da - call CCTK_Reduce ( ierr, cctkGH, -1, sum_handle, 1, CCTK_VARIABLE_REAL, & - centroid_x, 1, int_var ) - - int_tmp = sym_factor * interp_y * weights * da - call CCTK_Reduce ( ierr, cctkGH, -1, sum_handle, 1, CCTK_VARIABLE_REAL, & - centroid_y, 1, int_var ) - - int_tmp = sym_factor * interp_z * weights * da - call CCTK_Reduce ( ierr, cctkGH, -1, sum_handle, 1, CCTK_VARIABLE_REAL, & - centroid_z, 1, int_var ) - - if ( s_symx ) centroid_x = zero - if ( s_symy ) centroid_y = zero - if ( s_symz ) centroid_z = zero - - call CCTK_GrouplbndGN ( ierr, cctkGH, 1, lbnd, "ehfinder::eh_centroid2_x" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get lower bounds for centroid array" ) - end if - call CCTK_GroupubndGN ( ierr, cctkGH, 1, ubnd, "ehfinder::eh_centroid2_x" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get upper bounds for centroid array" ) - end if - call CCTK_GrouplshGN ( ierr, cctkGH, 1, lsh, "ehfinder::eh_centroid2_x" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for centroid array" ) - end if - - if ( ( sn .ge. lbnd(1) + 1 ) .and. ( sn .le. ubnd(1) + 1 ) ) then - k = sn - lbnd(1) - eh_centroid2_x(k,l) = centroid_x / eh_area2(k,l) - eh_centroid2_y(k,l) = centroid_y / eh_area2(k,l) - eh_centroid2_z(k,l) = centroid_z / eh_area2(k,l) - end if - - write(info_message,'(a19,3(f10.5))') 'Horizon centroid = ', & - eh_centroid2_x(sn,l), & - eh_centroid2_y(sn,l), & - eh_centroid2_z(sn,l) - call CCTK_INFO(info_message) - -end subroutine EHFinder_IntegrateCentroid - - -subroutine EHFinder_IntegrateEquatorial(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - - CCTK_INT :: int_var, sn, k, l - character(len=50) :: info_message - CCTK_INT, dimension(1) :: lbnd1, ubnd1, lsh1 - CCTK_INT, dimension(2) :: lbnd, lsh, nghost - CCTK_INT, dimension(4) :: bbox - CCTK_INT :: ithl, ithr, jphl, jphr, itha, jpha - CCTK_REAL :: eq_circ_loc, eq_circ, pol_circ_loc, pol_circ - -! If finding of surface failed do not try to integrate but exit. - if ( find_surface_status .lt. 0 ) then - return - endif - - l = levelset_counter - - sn = surface_counter - integrate_counter - - ! Find out the lower bounds of the distributed integration grid arrays. - call CCTK_GrouplbndGN ( ierr, cctkGH, 2, lbnd, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get lower bounds for surface arrays" ) - end if - - ! Find out the local size of the distributed integration grid arrays - call CCTK_GrouplshGN ( ierr, cctkGH, 2, lsh, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for surface arrays" ) - end if - - ! Find out the bounding box of the distributed integration grid arrays - call CCTK_GroupbboxGN ( ierr, cctkGH, 4, bbox, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get bounding box for surface arrays" ) - end if - - ! Find out the ghost size of the distributed integration grid arrays - call CCTK_GroupnghostzonesGN ( ierr, cctkGH, 2, nghost, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get ghost zones for surface arrays" ) - end if - - ithl = 1; ithr = lsh(1) - jphl = 1; jphr = lsh(2) - - if ( bbox(1) .eq. 0 ) ithl = ithl + nghost(1) - if ( bbox(2) .eq. 0 ) ithr = ithr - nghost(1) - if ( bbox(3) .eq. 0 ) jphl = jphl + nghost(2) - if ( bbox(4) .eq. 0 ) jphr = jphr - nghost(2) - - call CCTK_INFO ( 'Integrating equatorial circumference' ) - - int_tmp = phi_sym_factor * phiweights * dlphi - - if ( ( ithl+lbnd(1) .le. githeta ) .and. ( githeta .le. ithr+lbnd(1) ) ) then - itha = githeta - lbnd(1) - eq_circ_loc = sum ( int_tmp(itha,jphl:jphr) ) - else - eq_circ_loc = zero - end if - - call CCTK_ReduceLocScalar ( ierr, cctkGH, -1, sum_handle, & - eq_circ_loc, eq_circ, CCTK_VARIABLE_REAL ) - - call CCTK_GrouplbndGN ( ierr, cctkGH, 1, lbnd1, "ehfinder::eh_circ_eq2" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get lower bounds for area array" ) - end if - call CCTK_GroupubndGN ( ierr, cctkGH, 1, ubnd1, "ehfinder::eh_circ_eq2" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get upper bounds for area array" ) - end if - call CCTK_GrouplshGN ( ierr, cctkGH, 1, lsh1, "ehfinder::eh_circ_eq2" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for area array" ) - end if - - if ( ( sn .ge. lbnd1(1) + 1 ) .and. ( sn .le. ubnd1(1) + 1 ) ) then - k = sn - lbnd1(1) - eh_circ_eq2(k,l) = eq_circ - end if - - write(info_message,'(a35,f10.5)') 'Horizon equatorial circumference = ', & - eq_circ - call CCTK_INFO(info_message) - - call CCTK_INFO ( 'Integrating polar circumference' ) - - int_tmp = theta_sym_factor * thetaweights * dltheta - - if ( ( jphl+lbnd(2) .le. gjphi ) .and. ( gjphi .le. jphr+lbnd(2) ) ) then - jpha = gjphi - lbnd(2) - pol_circ_loc = sum ( int_tmp(ithl:ithr,jpha) ) - else - pol_circ_loc = zero - end if - - call CCTK_ReduceLocScalar ( ierr, cctkGH, -1, sum_handle, & - pol_circ_loc, pol_circ, CCTK_VARIABLE_REAL ) - - - if ( ( sn .ge. lbnd1(1) + 1 ) .and. ( sn .le. ubnd1(1) + 1 ) ) then - k = sn - lbnd1(1) - eh_circ_pol2(k,l) = pol_circ - end if - - write(info_message,'(a30,f10.5)') 'Horizon polar circumference = ',pol_circ - call CCTK_INFO(info_message) - -end subroutine EHFinder_IntegrateEquatorial - - -! This routine sets up the weights for the Simpsons rule integration -! over the surface. -subroutine EHFinder_InitWeights(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - - CCTK_INT :: i, j, k, im, jm - - CCTK_INT, dimension(4) :: bbox - CCTK_INT, dimension(2) :: gsh, lsh, lbnd, ubnd, nghost - - ! Find out the lower bounds of the distributed integration grid arrays. - call CCTK_GrouplbndGN ( ierr, cctkGH, 2, lbnd, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get lower bounds for surface arrays" ) - end if - - ! Find out the local size of the distributed integration grid arrays - call CCTK_GrouplshGN ( ierr, cctkGH, 2, lsh, "ehfinder::surface_arrays" ) - if ( ierr .lt. 0 ) then - call CCTK_WARN ( 0, "cannot get local size for surface arrays" ) - end if - - weights = one - - ! Initialise the weight grid array for the 2D Simpsons rule integration. - ! To do this I need to figure out the global location of the given point. - ! There are 3 cases in the one dimensional case. If the point is on the - ! boundary the weight is 1/3. If it is at an even position the weight - ! is 4/3 and if it is at an odd position the weight is 2/3. - - do j = 1, lsh(2) - - ! This is first done in the phi direction. Meaning that all points with - ! the same theta coordinate are set to the same weight. - if ( ( lbnd(2)+j .eq. 1 ) .or. ( lbnd(2)+j .eq. nphi ) ) then - weights(:,j) = onethird - phiweights(:,j) = onethird - else if ( mod(lbnd(2)+j,2) .eq. 0 ) then - weights(:,j) = fourthirds - phiweights(:,j) = fourthirds - else - weights(:,j) = twothirds - phiweights(:,j) = twothirds - end if - - ! Then it is done in the theta direction with the one-dimensional - ! weights beeing multiplied. - do i = 1, lsh(1) - if ( ( lbnd(1)+i .eq. 1 ) .or. ( lbnd(1)+i .eq. ntheta ) ) then - weights(i,j) = onethird * weights(i,j) - thetaweights(i,j) = onethird - else if ( mod(lbnd(1)+i,2) .eq. 0 ) then - weights(i,j) = fourthirds * weights(i,j) - thetaweights(i,j) = fourthirds - else - weights(i,j) = twothirds * weights(i,j) - thetaweights(i,j) = twothirds - end if - end do - end do - - ! The end result is a 2D array with the weights in the following pattern. - - ! 1/9 4/9 2/9 4/9 2/9 4/9 1/9 - ! 4/9 16/9 8/9 16/9 8/9 16/9 4/9 - ! 2/9 8/9 4/9 8/9 4/9 8/9 2/9 - ! 4/9 16/9 8/9 16/9 8/9 16/9 4/9 - ! 2/9 8/9 4/9 8/9 4/9 8/9 2/9 - ! 4/9 16/9 8/9 16/9 8/9 16/9 4/9 - ! 1/9 4/9 2/9 4/9 2/9 4/9 1/9 -end subroutine EHFinder_InitWeights - - -subroutine EHFinder_CopyArea(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - -! If finding of surface failed set area to zero. - if ( find_surface_status .lt. 0 ) then - eh_area = zero - return - else - call CCTK_INFO('Copying area data') - eh_area = eh_area2 - endif -end subroutine EHFinder_CopyArea - - -subroutine EHFinder_CopyCentroid(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - -! If finding of surface failed set centroids to zero. - if ( find_surface_status .lt. 0 ) then - eh_centroid_x = zero - eh_centroid_y = zero - eh_centroid_z = zero - return - else - call CCTK_INFO('Copying centroid data') - eh_centroid_x = eh_centroid2_x - eh_centroid_y = eh_centroid2_y - eh_centroid_z = eh_centroid2_z - end if - -end subroutine EHFinder_CopyCentroid - - -subroutine EHFinder_CopyCircumference(CCTK_ARGUMENTS) - - use EHFinder_mod - - implicit none - - DECLARE_CCTK_PARAMETERS - DECLARE_CCTK_ARGUMENTS - DECLARE_CCTK_FUNCTIONS - - if ( find_surface_status .lt. 0 ) then - eh_circ_eq = zero - eh_circ_pol = zero - return - else - call CCTK_INFO('Copying circumference data') - eh_circ_eq = eh_circ_eq2 - eh_circ_pol = eh_circ_pol2 - end if - -end subroutine EHFinder_CopyCircumference -- cgit v1.2.3