aboutsummaryrefslogtreecommitdiff
path: root/src/D3_to_D2.F
diff options
context:
space:
mode:
Diffstat (limited to 'src/D3_to_D2.F')
-rw-r--r--src/D3_to_D2.F398
1 files changed, 398 insertions, 0 deletions
diff --git a/src/D3_to_D2.F b/src/D3_to_D2.F
new file mode 100644
index 0000000..884f8b8
--- /dev/null
+++ b/src/D3_to_D2.F
@@ -0,0 +1,398 @@
+
+#include "cctk.h"
+
+c ========================================================================
+
+ SUBROUTINE D3_to_D2(cctkGH,do_ADMmass,do_momentum,do_spin,
+ & Psi_power,origin,myproc,Dx,Dy,Dz,Psi,
+ & g00,gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & x,y,z,eta,Nt,Np,theta,phi,
+ & Psis,g00s,gxxs,gxys,gxzs,gyys,gyzs,gzzs,dPsis,dgxxs,dgxys,dgxzs,
+ & dgyys,dgyzs,dgzzs,ADMmass_int1,ADMmass_int2,
+ & momentum_int1,momentum_int2,momentum_int3,
+ & spin_int1,spin_int2,spin_int3,Extract_temp3d)
+
+c ------------------------------------------------------------------------
+c
+c Project the 3-metric and its 1st radial derivatives onto a
+c 2-sphere.
+c
+c ------------------------------------------------------------------------
+
+ USE met_rad_der_int
+ USE ADMmass_integrand3D_int
+ USE momentum_integrand3D_int
+ USE spin_integrand3D_int
+
+ IMPLICIT NONE
+
+c Input variables
+
+ CCTK_POINTER :: cctkGH
+
+ INTEGER,INTENT(IN) ::
+ & myproc,Psi_power
+ CCTK_INT, INTENT(IN) ::
+ & Nt,Np,do_momentum,do_spin
+ CCTK_REAL,INTENT(IN) ::
+ & origin(3),Dx,Dy,Dz,eta
+ CCTK_REAL,INTENT(IN),DIMENSION(:) ::
+ & theta,phi,x,y,z
+ CCTK_REAL,INTENT(IN),DIMENSION(:,:,:) ::
+ & Psi,g00,gxx,gxy,gxz,gyy,gyz,gzz,
+ & hxx,hxy,hxz,hyy,hyz,hzz
+ CCTK_REAL,INTENT(INOUT),DIMENSION(:,:,:) ::
+ & Extract_temp3d
+ LOGICAL ::
+ & do_ADMmass(2)
+
+c Output variables
+
+ CCTK_REAL,INTENT(OUT),DIMENSION(:,:) ::
+ & Psis,g00s,gxxs,gxys,gxzs,gyys,
+ & gyzs,gzzs,dPsis,dgxxs,dgxys,dgxzs,dgyys,dgyzs,dgzzs,
+ & ADMmass_int1,ADMmass_int2,
+ & momentum_int1,momentum_int2,momentum_int3,
+ & spin_int1,spin_int2,spin_int3
+
+
+c Local variables, passed on
+
+ LOGICAL ::
+ & err_flag
+ INTEGER ::
+ & iorder,npoints,Nx,Ny,Nz
+ INTEGER,DIMENSION(Nt,Np) ::
+ & ib,jb,kb
+ CCTK_REAL,DIMENSION(Nt,Np) ::
+ & xs,ys,zs,ux,uy,uz,xb,yb,zb
+
+c Local variables, here only
+
+ INTEGER ::
+ & i,j,handle,ierror
+ CCTK_REAL ::
+ & xmin,xmax,ymin,ymax,zmin,zmax
+
+c ------------------------------------------------------------------------
+
+c Initial Stuff
+c -------------
+ Nx = SIZE(x)
+ Ny = SIZE(y)
+ Nz = SIZE(z)
+
+
+c Compute Cartesian coordinates on the surface
+c --------------------------------------------
+ DO j = 1, Np
+ DO i = 1, Nt
+
+ xs(i,j) = origin(1)+eta*SIN(theta(i))*COS(phi(j))
+ ys(i,j) = origin(2)+eta*SIN(theta(i))*SIN(phi(j))
+ zs(i,j) = origin(3)+eta*COS(theta(i))
+
+ ENDDO
+ ENDDO
+
+
+c Only do interpolation on one processor
+c --------------------------------------
+ SELECT CASE (myproc)
+
+ CASE (0)
+ npoints = Np*Nt
+
+ CASE DEFAULT
+ npoints = 0
+
+ END SELECT
+
+
+c Choose the interpolator
+c -----------------------
+ call CCTK_GetInterpHandle (handle, "simple")
+
+c Find the origin of the spatial coordinates
+c ------------------------------------------
+ call CCTK_CoordRange(ierror,cctkGH,xmin,xmax,"x")
+ call CCTK_CoordRange(ierror,cctkGH,ymin,ymax,"y")
+ call CCTK_CoordRange(ierror,cctkGH,zmin,zmax,"z")
+
+
+c Project un-physical metric and conformal factor onto sphere
+c ------------------------------------------------------------
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 8, 8,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL,
+ $ CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Psi,g00,gxx,gxy,gxz,gyy,gyz,gzz,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ Psis,g00s,gxxs,gxys,gxzs,gyys,gyzs,gzzs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL
+ $ )
+
+
+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(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dPsis,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gxx,Extract_temp3d)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dgxxs,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gxy,Extract_temp3d)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dgxys,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gxz,Extract_temp3d)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dgxzs,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gyy,Extract_temp3d)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dgyys,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gyz,Extract_temp3d)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dgyzs,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+
+ CALL met_rad_der(origin,Dx,Dy,Dz,x,y,z,gzz,Extract_temp3d)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ dgzzs,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+
+c Calculate integrands for ADM masses
+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)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ ADMmass_int1,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ END IF
+
+c Conformal equation
+ IF (do_ADMmass(2)) THEN
+ ADMmass_int2 = -eta**2/2D0/3.1416D0*dPsis
+ ENDIF
+
+c Calculate integrands for momentum
+c ---------------------------------
+ IF (do_momentum==1) THEN
+
+ CALL momentum_integrand3D(origin,Dx,Dy,Dz,x,y,z,
+ & gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & Extract_temp3d,Psi,Psi_power)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ momentum_int1,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL momentum_integrand3D(origin,Dx,Dy,Dz,x,y,z,
+ & gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & Extract_temp3d,Psi,Psi_power)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ momentum_int2,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+
+ CALL momentum_integrand3D(origin,Dx,Dy,Dz,x,y,z,
+ & gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & Extract_temp3d,Psi,Psi_power)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ momentum_int3,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ END IF
+
+c Calculate integrands for spin
+c -----------------------------
+ IF (do_momentum==1) THEN
+
+ CALL spin_integrand3D(origin,Dx,Dy,Dz,x,y,z,
+ & gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & Extract_temp3d,Psi,Psi_power)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ spin_int1,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ CALL spin_integrand3D(origin,Dx,Dy,Dz,x,y,z,
+ & gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & Extract_temp3d,Psi,Psi_power)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ spin_int2,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+
+ CALL spin_integrand3D(origin,Dx,Dy,Dz,x,y,z,
+ & gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
+ & Extract_temp3d,Psi,Psi_power)
+ CALL CCTK_SyncGroup(cctkGH,"extract::temps")
+ call CCTK_Interp (ierror, cctkGH, handle, npoints, 3, 1, 1,
+ $ nx, ny, nz,
+ $ xs, ys, zs,
+ $ CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,CCTK_VARIABLE_REAL,
+ $ xmin, ymin, zmin,
+ $ dx, dy, dz,
+ $ Extract_temp3d,
+ $ CCTK_VARIABLE_REAL,
+ $ spin_int3,
+ $ CCTK_VARIABLE_REAL
+ $ )
+
+ END IF
+
+ END SUBROUTINE D3_to_D2
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
generated by cgit v1.2.3 (git 2.39.1) at 2024-06-01 17:14:51 +0000