From 23c78444b4fdbea008ccca1637c4ab26b955237d Mon Sep 17 00:00:00 2001 From: allen Date: Tue, 13 Aug 2002 08:48:45 +0000 Subject: Fix extract for physical metrics Closes Cactus/1196 git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/Extract/trunk@60 5301f0c2-dbc4-4cee-b2f5-8d7afba4d129 --- src/ADMmass_integrand3D.F | 57 +++++++++++++++++++++++++++++++++---------- src/ADMmass_integrand3D_int.F | 4 +-- src/D3_to_D2.F | 25 ++++++++++--------- 3 files changed, 60 insertions(+), 26 deletions(-) diff --git a/src/ADMmass_integrand3D.F b/src/ADMmass_integrand3D.F index 06b1198..a501ef1 100644 --- a/src/ADMmass_integrand3D.F +++ b/src/ADMmass_integrand3D.F @@ -1,9 +1,10 @@ #include "cctk.h" +#include "cctk_Parameters.h" c ======================================================================== SUBROUTINE ADMmass_integrand3D(origin,Dx,Dy,Dz,x,y,z,gxx,gxy,gxz, - & gyy,gyz,gzz,ADMmass_int,Psi,Psi_power) + & gyy,gyz,gzz,ADMmass_int,Psi,Psi_power,conformal_state) c ------------------------------------------------------------------------ c @@ -16,7 +17,7 @@ c ------------------------------------------------------------------------ c Input variables INTEGER,INTENT(IN) :: - & Psi_power + & Psi_power,conformal_state CCTK_REAL,INTENT(IN) :: & Dx,Dy,Dz,origin(3) CCTK_REAL,DIMENSION(:),INTENT(IN) :: @@ -39,6 +40,9 @@ c Local variables, here only & dyy_x,dxz_x,dxz_y,dxz_z,dyz_x,dzz_x,dyy_z,dyz_y,dyz_z,dzz_y, & Pi,idet,p,pip,pim,pjp,pjm,pkp,pkm + DECLARE_CCTK_PARAMETERS + DECLARE_CCTK_FUNCTIONS + c ------------------------------------------------------------------------ Pi = ACOS(-1D0) @@ -74,11 +78,24 @@ c Because other codes evolve Psi**4 c Abbreviations for metric coefficients c ------------------------------------- - p = psi(i,j,k)**ip - dxx = p*gxx(i,j,k); dxy = p*gxy(i,j,k) - dxz = p*gxz(i,j,k); dyy = p*gyy(i,j,k) - dyz = p*gyz(i,j,k); dzz = p*gzz(i,j,k) + if (conformal_state > 0) then + + p = psi(i,j,k)**ip + + dxx = p*gxx(i,j,k); dxy = p*gxy(i,j,k) + dxz = p*gxz(i,j,k); dyy = p*gyy(i,j,k) + dyz = p*gyz(i,j,k); dzz = p*gzz(i,j,k) + + else + + p = 1 + + dxx = gxx(i,j,k); dxy = gxy(i,j,k) + dxz = gxz(i,j,k); dyy = gyy(i,j,k) + dyz = gyz(i,j,k); dzz = gzz(i,j,k) + + end if c Determinant of 3-metric c ----------------------- @@ -98,12 +115,26 @@ c ---------------- c Derivatives of the 3-metric c --------------------------- - pip = psi(i+1,j,k)**ip - pim = psi(i-1,j,k)**ip - pjp = psi(i,j+1,k)**ip - pjm = psi(i,j-1,k)**ip - pkp = psi(i,j,k+1)**ip - pkm = psi(i,j,k-1)**ip + + IF (conformal_state > 0) THEN + + pip = psi(i+1,j,k)**ip + pim = psi(i-1,j,k)**ip + pjp = psi(i,j+1,k)**ip + pjm = psi(i,j-1,k)**ip + pkp = psi(i,j,k+1)**ip + pkm = psi(i,j,k-1)**ip + + else + + pip = 1.0d0 + pim = 1.0d0 + pjp = 1.0d0 + pjm = 1.0d0 + pkp = 1.0d0 + pkm = 1.0d0 + + end if dxx_y = half/Dy*(pjp*gxx(i,j+1,k)-pjm*gxx(i,j-1,k)) dxx_z = half/Dz*(pkp*gxx(i,j,k+1)-pkm*gxx(i,j,k-1)) @@ -141,7 +172,7 @@ c --------------------------- ADMmass_int(i,j,k) = 0.0D0 ENDIF - + ENDDO ENDDO ENDDO diff --git a/src/ADMmass_integrand3D_int.F b/src/ADMmass_integrand3D_int.F index c103f55..b9f98d3 100644 --- a/src/ADMmass_integrand3D_int.F +++ b/src/ADMmass_integrand3D_int.F @@ -8,10 +8,10 @@ c ------------------------------------------------------------------ INTERFACE SUBROUTINE ADMmass_integrand3D(origin,Dx,Dy,Dz,x,y,z,gxx,gxy, - & gxz,gyy,gyz,gzz,ADMmass_int,Psi,Psi_power) + & gxz,gyy,gyz,gzz,ADMmass_int,Psi,Psi_power,conformal_state) IMPLICIT NONE INTEGER,INTENT(IN) :: - & Psi_power + & Psi_power,conformal_state CCTK_REAL,INTENT(IN) :: & Dx,Dy,Dz,origin(3) CCTK_REAL,DIMENSION(:),INTENT(IN) :: diff --git a/src/D3_to_D2.F b/src/D3_to_D2.F index bc684dc..2948517 100644 --- a/src/D3_to_D2.F +++ b/src/D3_to_D2.F @@ -177,17 +177,20 @@ c ------------------------------------------------------------ c Calculate radial derivatives and project onto sphere c ---------------------------------------------------- - CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,Psi,Extract_temp3d) - CALL CCTK_SyncGroup(ierror,cctkGH,"extract::temps") call CCTK_VarIndex(in_array_indices(1), "extract::temp3d") - call CCTK_InterpGV (ierror, cctkGH, interp_handle, coord_system_handle, - $ npoints, 1, 1, - $ xs, ys, zs, - $ CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL, - $ CCTK_VARIABLE_REAL, - $ in_array_indices(1), - $ dPsis, - $ CCTK_VARIABLE_REAL) + + if (conformal_state > 0) then + CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,Psi,Extract_temp3d) + CALL CCTK_SyncGroup(ierror,cctkGH,"extract::temps") + call CCTK_InterpGV (ierror, cctkGH, interp_handle, coord_system_handle, + $ npoints, 1, 1, + $ xs, ys, zs, + $ CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL, + $ CCTK_VARIABLE_REAL, + $ in_array_indices(1), + $ dPsis, + $ CCTK_VARIABLE_REAL) + end if CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gxx,Extract_temp3d) CALL CCTK_SyncGroup(ierror,cctkGH,"extract::temps") @@ -261,7 +264,7 @@ c ----------------------------------- c Standard equation IF (do_ADMmass(1)) THEN CALL ADMmass_integrand3D(origin,Dx,Dy,Dz,x,y,z,gxx,gxy, - & gxz,gyy,gyz,gzz,Extract_temp3d,Psi,Psi_power) + & gxz,gyy,gyz,gzz,Extract_temp3d,Psi,Psi_power,conformal_state) CALL CCTK_SyncGroup(ierror,cctkGH,"extract::temps") call CCTK_InterpGV (ierror, cctkGH, interp_handle, coord_system_handle, $ npoints, 1, 1, -- cgit v1.2.3