Work with a physical metric

git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/Extract/trunk@54 5301f0c2-dbc4-4cee-b2f5-8d7afba4d129

6 files changed, 35 insertions, 16 deletions

diff --git a/src/D3_extract.F b/src/D3_extract.F
index 4fb6899..9f6fcad 100644
--- a/src/D3_extract.F
+++ b/src/D3_extract.F
@@ -3,7 +3,7 @@
 
 c     ==================================================================
 
-      SUBROUTINE D3_extract(cctkGH,do_ADMmass,do_momentum,do_spin,igrid,
+      SUBROUTINE D3_extract(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,igrid,
      &     origin,myproc,interpolation_order,Nt,Np,all_modes,
      &     l,m,x,y,z,Dx,Dy,Dz,Psi_power,Psi,
      &     g00,gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
@@ -71,7 +71,7 @@ c     Input variables
       CCTK_POINTER :: cctkGH
 
       INTEGER,INTENT(IN) :: 
-     &     igrid,l,m,Psi_power,myproc
+     &     conformal_state,igrid,l,m,Psi_power,myproc
       CCTK_INT,INTENT(IN) :: 
      &     Nt,Np,all_modes,do_momentum,do_spin,interpolation_order
       INTEGER,INTENT(IN) ::
@@ -161,7 +161,7 @@ c     2. Project quantities onto the 2-surface
 c
 c     ------------------------------------------------------------------
 
-      CALL D3_to_D2(cctkGH,do_ADMmass,do_momentum,do_spin,
+      CALL D3_to_D2(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,
      &     Psi_power,origin,myproc,interpolation_order,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,
@@ -189,9 +189,10 @@ c     4. Calculate physical quantities on 2-surface
 c
 c     ------------------------------------------------------------------
 
-         CALL unphysical_to_physical(grr,grt,grp,gtt,gtp,gpp,dgtt,dgtp,
+         if (conformal_state > 0) then
+           CALL unphysical_to_physical(grr,grt,grp,gtt,gtp,gpp,dgtt,dgtp,
      &        dgpp,Psis,dPsis,Psi_power)
-
+         end if
 
 c     ------------------------------------------------------------------
 c
diff --git a/src/D3_extract_int.F b/src/D3_extract_int.F
index 357f995..1d7244a 100644
--- a/src/D3_extract_int.F
+++ b/src/D3_extract_int.F
@@ -7,7 +7,7 @@ c     ------------------------------------------------------------------
 
       INTERFACE
 
-      SUBROUTINE D3_extract(cctkGH,do_ADMmass,do_momentum,do_spin,igrid,
+      SUBROUTINE D3_extract(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,igrid,
      &        origin,myproc,interpolation_order,
      &        Nt,Np,all_modes,l,m,x,y,z,Dx,Dy,Dz,Psi_power,Psi,
      &        g00,gxx,gxy,gxz,gyy,gyz,gzz,hxx,hxy,hxz,hyy,hyz,hzz,
@@ -19,7 +19,7 @@ c     ------------------------------------------------------------------
       CCTK_POINTER :: cctkGH
 
       INTEGER,INTENT(IN) :: 
-     &     igrid,l,m,Psi_power,myproc
+     &     conformal_state,igrid,l,m,Psi_power,myproc
       CCTK_INT,INTENT(IN) ::
      &     Nt,Np,all_modes,do_momentum,do_spin,interpolation_order
       INTEGER,INTENT(IN) ::
diff --git a/src/D3_to_D2.F b/src/D3_to_D2.F
index 3203a39..bc684dc 100644
--- a/src/D3_to_D2.F
+++ b/src/D3_to_D2.F
@@ -3,7 +3,7 @@
 
 c     ========================================================================
 
-      SUBROUTINE D3_to_D2(cctkGH,do_ADMmass,do_momentum,do_spin,
+      SUBROUTINE D3_to_D2(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,
      &     Psi_power,origin,myproc,interpolation_order,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,
@@ -31,7 +31,7 @@ c     Input variables
       CCTK_POINTER :: cctkGH
 
       INTEGER,INTENT(IN) :: 
-     &     myproc,Psi_power
+     &     conformal_state,myproc,Psi_power
       CCTK_INT, INTENT(IN) ::
      &     Nt,Np,do_momentum,do_spin,interpolation_order
       CCTK_REAL,INTENT(IN) :: 
@@ -142,7 +142,8 @@ c     ---------------------------------
 
 c     Project un-physical metric and conformal factor  onto sphere
 c     ------------------------------------------------------------
-      call CCTK_InterpGV (ierror, cctkGH, interp_handle, coord_system_handle,
+      if (conformal_state > 0) then
+        call CCTK_InterpGV (ierror, cctkGH, interp_handle, coord_system_handle,
      $                    npoints, 8, 8,
      $                    xs, ys, zs,
      $                    CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL,
@@ -156,6 +157,22 @@ c     ------------------------------------------------------------
      $                    CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL,
      $                    CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL,
      $                    CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL)
+      else
+        call CCTK_InterpGV (ierror, cctkGH, interp_handle, coord_system_handle,
+     $                    npoints, 7, 7,
+     $                    xs, ys, zs,
+     $                    CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL,
+     $                    CCTK_VARIABLE_REAL,
+     $                    in_array_indices(2),
+     $                    in_array_indices(3), in_array_indices(4),
+     $                    in_array_indices(5), in_array_indices(6),
+     $                    in_array_indices(7), in_array_indices(8),
+     $                    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)
+      end if
 
 
 c     Calculate radial derivatives and project onto sphere
diff --git a/src/D3_to_D2_int.F b/src/D3_to_D2_int.F
index 8f19412..8ca896c 100644
--- a/src/D3_to_D2_int.F
+++ b/src/D3_to_D2_int.F
@@ -7,7 +7,7 @@ c     ------------------------------------------------------------------
 
       INTERFACE
 
-      SUBROUTINE D3_to_D2(cctkGH,do_ADMmass,do_momentum,do_spin,
+      SUBROUTINE D3_to_D2(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,
      &        Psi_power,origin,myproc,interpolation_order,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,
@@ -19,7 +19,7 @@ c     ------------------------------------------------------------------
       IMPLICIT NONE
       CCTK_POINTER :: cctkGH
       INTEGER,INTENT(IN) :: 
-     &     myproc,Psi_power
+     &     conformal_state,myproc,Psi_power
       CCTK_INT, INTENT(IN) ::
      &     Nt,Np,do_momentum,do_spin,interpolation_order
       INTEGER ::
diff --git a/src/Extract.F b/src/Extract.F
index 7e95ac6..e7d9f56 100644
--- a/src/Extract.F
+++ b/src/Extract.F
@@ -377,7 +377,7 @@ c           Do extraction at each radius
 
             extract_at_each_radius: DO WHILE (radius < r2)
 
-               CALL D3_extract(cctkGH,do_ADMmass,do_momentum,do_spin,
+               CALL D3_extract(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,
      &           igrid,orig,myproc,interpolation_order,Nt,Np,all_modes,lmode,
      &           mmode,x_1d,y_1d,z_1d,Dx,Dy,Dz,Psi_power,Psi,g00,
      &           gxx,gxy,gxz,gyy,gyz,gzz,kxx,kxy,kxz,kyy,kyz,kzz,
@@ -548,7 +548,7 @@ c     Cannot use the conformal equation for ADM mass now
               WRITE(*,*) "Calling extract at radius ... ",radius
             END IF
 
-            CALL D3_extract(cctkGH,do_ADMmass,do_momentum,do_spin,
+            CALL D3_extract(cctkGH,conformal_state,do_ADMmass,do_momentum,do_spin,
      &           igrid,orig,myproc,interpolation_order,Nt,Np,all_modes,
      &           lmode,mmode,x_1d,y_1d,z_1d,
      &           Dx,Dy,Dz,Psi_power,Psi,g00,gxx,gxy,gxz,gyy,gyz,gzz,kxx,kxy,kxz,
diff --git a/src/ParamCheck.c b/src/ParamCheck.c
index c5a4684..ae68e81 100644
--- a/src/ParamCheck.c
+++ b/src/ParamCheck.c
@@ -37,8 +37,9 @@ void Extract_ParamCheck(CCTK_ARGUMENTS)
   DECLARE_CCTK_ARGUMENTS
   DECLARE_CCTK_PARAMETERS
 
-  if(! CCTK_EQUALS(metric_type, "static conformal"))
+  if(! (CCTK_EQUALS(metric_type, "physical") || 
+	CCTK_EQUALS(metric_type, "static conformal")))
   {
-    CCTK_PARAMWARN("Extract only works currently with metric_type \"static conformal\"");
+    CCTK_PARAMWARN("Extract only works currently with metric_type \"static conformal\" or \"physical\"");
   }
 }