Fixed a bug for periodic BC.

This closes PR Cactus/803.


git-svn-id: http://svn.cactuscode.org/arrangements/CactusPUGH/PUGHReduce/trunk@18 d60812e6-3970-4df4-986e-c251b06effeb

1 files changed, 253 insertions, 133 deletions

diff --git a/src/Reduction.c b/src/Reduction.c
index 035b126..bce75e9 100644
--- a/src/Reduction.c
+++ b/src/Reduction.c
@@ -3,16 +3,16 @@
    @date      Thu Apr  3 11:54:53 1997
    @author    Thomas Radke, Paul Walker
    @desc 
-   Various MPI reduction operators. 
+              Various MPI reduction operators. 
    @enddesc 
-   @version $Header$
+   @version   $Id$
  @@*/
 
 #include <stdlib.h>
 
 #include "pugh_reductions.h"
 
-static const char *rcsid = "$Id$";
+static const char *rcsid = "$Header$";
 
 CCTK_FILEVERSION(CactusPUGH_PUGHReduce_Reduction_c)
 
@@ -21,12 +21,13 @@ CCTK_FILEVERSION(CactusPUGH_PUGHReduce_Reduction_c)
  ********************************************************************/
 
 static int PUGH_ReductionGA (cGH *GH, 
-			     int index, 
-			     int proc, 
-			     CCTK_REAL *outval,
+                             int index, 
+                             int proc, 
+                             CCTK_REAL *outval,
                              reduction_fn_t reduction_fn);
 
-static int PUGH_ReductionScalar (cGH *GH, int index, int proc, CCTK_REAL *outval,
+static int PUGH_ReductionScalar (cGH *GH, int index, int proc,
+                                 CCTK_REAL *outval,
                                  reduction_fn_t reduction_fn);
 
 static int copy_real_to_outtype (int num_elems,
@@ -41,26 +42,77 @@ static int copy_real_to_outtype (int num_elems,
 
  /*@@
    @routine    PUGH_ReductionArrays
-   @date       
-   @author     
+   @author     Thomas Radke
+   @date       19 Aug 1999
    @desc 
+               Wrapper to reduce a list of arrays.
+               Just calls the appropriate reduction operator and does
+               the type conversion of the results.
+   @enddesc
+   @calls      copy_real_to_outtype
+  
+   @var        GH
+   @vdesc      Pointer to CCTK grid hierarchy
+   @vtype      cGH *
+   @vio        in
+   @endvar
+   @var        proc
+   @vdesc      processor that should receive the result of operation
+               (negative value means all processors receive the result)
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        num_dims
+   @vdesc      number of dimensions of input arrays
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        dims
+   @vdesc      dimensions of input arrays
+   @vtype      int *
+   @vio        in
+   @endvar
+   @var        intype
+   @vdesc      (common) variable type of input arrays
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        num_inarrays
+   @vdesc      number of input arrays
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        inarrays
+   @vdesc      field of input arrays
+   @vtype      void **
+   @vio        in
+   @endvar
+   @var        outtype
+   @vdesc      (common) variable type of output arrays
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        num_outvals
+   @vdesc      number of values per output array
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        outvals
+   @vdesc      pointer to buffer holding the output values
+   @vtype      void *
+   @vio        in
+   @endvar
+   @var        reduction_fn
+   @vdesc      reduction operator callback
+   @vtype      reduction_fn_t
+   @vio        in
+   @endvar
 
-   @enddesc 
-   @calls      
-   @calledby   
-   @history 
- 
-   @endhistory 
-   @var     
-   @vdesc   
-   @vtype   
-   @vio     
-   @vcomment 
- 
-   @endvar 
-
+   @returntype int
+   @returndesc
+               the return code of the reduction operator
+   @endreturndesc
 @@*/
-
 int PUGH_ReductionArrays (cGH *GH,
                           int proc,
                           int num_dims,
@@ -73,16 +125,12 @@ int PUGH_ReductionArrays (cGH *GH,
                           void *outvals /* [num_outvals] */,
                           reduction_fn_t reduction_fn)
 {
-  int retval;
-  int i;
-  int num_points;
-  int from[1];
-  int to[1];
-  int iterator[1];
-  int points_per_dim[1];
+  int i, num_points, retval;
+  int from[1], to[1], iterator[1], points_per_dim[1];
   int *intypes;
   CCTK_REAL *buffer;
 
+
   points_per_dim[0] = 1;
   from[0] = 0;
   to[0] = dims[0];
@@ -119,26 +167,26 @@ int PUGH_ReductionArrays (cGH *GH,
 
   /* do the reduction on the input arrays */
   retval = reduction_fn (GH, 
-			 proc, 
-			 num_dims, 
-			 from, 
-			 to, 
-			 iterator, 
-			 points_per_dim,
+                         proc, 
+                         num_dims, 
+                         from, 
+                         to, 
+                         iterator, 
+                         points_per_dim,
                          num_points, 
-			 num_inarrays, 
-			 intypes, 
-			 inarrays,
+                         num_inarrays, 
+                         intypes, 
+                         inarrays,
                          num_outvals, 
-			 buffer);
+                         buffer);
 
   if (retval == 0 && (proc < 0 || proc == CCTK_MyProc (GH)))
   {
     /* type-cast the result to the requested datatype */
     retval = copy_real_to_outtype (num_inarrays * num_outvals, 
-				   buffer,
+                                   buffer,
                                    outtype, 
-				   outvals);
+                                   outvals);
   }
 
   free (intypes);
@@ -149,27 +197,67 @@ int PUGH_ReductionArrays (cGH *GH,
 
 
  /*@@
-   @routine    PUGH_ReductionArrays
-   @date       
-   @author     
+   @routine    PUGH_ReductionGVs
+   @author     Thomas Radke
+   @date       19 Aug 1999
    @desc 
+               Wrapper to reduce a list of grid variables.
+               Just calls the appropriate reduction operator and does
+               the type conversion of the results.
+   @enddesc
+   @calls      CCTK_VarTypeSize
+               CCTK_GroupTypeFromVarI
+               PUGH_ReductionGA
+               PUGH_ReductionScalar
+               copy_real_to_outtype
+  
+   @var        GH
+   @vdesc      Pointer to CCTK grid hierarchy
+   @vtype      cGH *
+   @vio        in
+   @endvar
+   @var        proc
+   @vdesc      processor that should receive the result of operation
+               (negative value means all processors receive the result)
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        num_invars
+   @vdesc      number of input arrays
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        inarrays
+   @vdesc      field of input arrays
+   @vtype      void **
+   @vio        in
+   @endvar
+   @var        outtype
+   @vdesc      (common) variable type of output arrays
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        num_outvals
+   @vdesc      number of values per output array
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        outvals
+   @vdesc      pointer to buffer holding the output values
+   @vtype      void *
+   @vio        in
+   @endvar
+   @var        reduction_fn
+   @vdesc      reduction operator callback
+   @vtype      reduction_fn_t
+   @vio        in
+   @endvar
 
-   @enddesc 
-   @calls      
-   @calledby   
-   @history 
- 
-   @endhistory 
-   @var     
-   @vdesc   
-   @vtype   
-   @vio     
-   @vcomment 
- 
-   @endvar 
-
+   @returntype int
+   @returndesc
+               the return code of the reduction operator
+   @endreturndesc
 @@*/
-
 int PUGH_ReductionGVs (cGH *GH,
                        int proc,
                        int num_invars,
@@ -240,30 +328,47 @@ int PUGH_ReductionGVs (cGH *GH,
 
  /*@@
    @routine    PUGH_ReductionGA
-   @date       
-   @author     
+   @date       Fri Jun 29 15:35:01 2001
+   @author     Tom Goodale
    @desc 
-
+               Reduction of a grid array variable.  
    @enddesc 
-   @calls      
-   @calledby   
-   @history 
- 
-   @endhistory 
-   @var     
-   @vdesc   
-   @vtype   
-   @vio     
-   @vcomment 
  
+   @var        GH
+   @vdesc      Pointer to CCTK grid hierarchy
+   @vtype      cGH *
+   @vio        in
+   @endvar 
+   @var        index
+   @vdesc      The GV index
+   @vtype      int
+   @vio        in
+   @endvar 
+   @var        proc
+   @vdesc      The processor at which we desire the result
+   @vtype      int
+   @vio        in
+   @endvar 
+   @var        outval
+   @vdesc      buffer to store the reduction result in
+   @vtype      CCTK_REAL *
+   @vio        out
+   @endvar 
+   @var        reduction_fn
+   @vdesc      The function which does the actual reduction
+   @vtype      reduction_fn_t
+   @vio        in
    @endvar 
 
+   @returntype int
+   @returndesc
+               The return code of the reduction operator.
+   @endreturndesc
 @@*/
-
 static int PUGH_ReductionGA (cGH *GH, int index, int proc, CCTK_REAL *outval,
                              reduction_fn_t reduction_fn)
 {
-  int i, retval, stagger_index, num_points;
+  int i, stagger_index, num_points, dir_points, retval;
   pGA *GA;
   int *from, *to, *iterator, *points_per_dim;
 
@@ -288,14 +393,16 @@ static int PUGH_ReductionGA (cGH *GH, int index, int proc, CCTK_REAL *outval,
     from[i] = GA->extras->ownership[stagger_index][0][i];
     to[i]   = GA->extras->ownership[stagger_index][1][i];
 
-    if (stagger_index == 0)
+    dir_points = GA->extras->nsize[i];
+    if (GA->connectivity->perme[i])
     {
-      num_points *= GA->extras->nsize[i];
+      dir_points -= 2*GA->extras->nghostzones[i];
     }
-    else
+    if (stagger_index)
     {
-      num_points *= GA->extras->nsize[i] - 1;
+      dir_points--;
     }
+    num_points *= dir_points;
 
     if (i > 0)
     {
@@ -314,51 +421,45 @@ static int PUGH_ReductionGA (cGH *GH, int index, int proc, CCTK_REAL *outval,
   return (retval);
 }
 
+
  /*@@
    @routine    PUGH_ReductionScalar
    @date       Fri Jun 29 15:35:01 2001
    @author     Tom Goodale
    @desc 
-   Reduction of a scalar GV.  
-   Basically a copy of PUGH_ReductionGA cut down to the scalar case.
+               Reduction of a scalar GV.  
+               Basically a copy of PUGH_ReductionGA cut down to the scalar case.
    @enddesc 
-   @calls     
-   @calledby   
-   @history 
- 
-   @endhistory 
-   @var     GH
-   @vdesc   A CCTK GH
-   @vtype   cGH *
-   @vio     in
-   @vcomment 
  
+   @var        GH
+   @vdesc      Pointer to CCTK grid hierarchy
+   @vtype      cGH *
+   @vio        in
    @endvar 
-   @var     index
-   @vdesc   The GV index
-   @vtype   int
-   @vio     in
-   @vcomment 
- 
+   @var        index
+   @vdesc      The GV index
+   @vtype      int
+   @vio        in
    @endvar 
-   @var     proc
-   @vdesc   The processor at which we desire the result
-   @vtype   int
-   @vio     in
-   @vcomment 
- 
+   @var        proc
+   @vdesc      The processor at which we desire the result
+   @vtype      int
+   @vio        in
    @endvar 
-   @var     reduction_fn
-   @vdesc   The function which does the actual reduction
-   @vtype   reduction_fn_t
-   @vio     in
-   @vcomment 
- 
+   @var        outval
+   @vdesc      buffer to store the reduction result in
+   @vtype      CCTK_REAL *
+   @vio        out
+   @endvar 
+   @var        reduction_fn
+   @vdesc      The function which does the actual reduction
+   @vtype      reduction_fn_t
+   @vio        in
    @endvar 
 
    @returntype int
    @returndesc
-      The return code of the reduction operator.
+               The return code of the reduction operator.
    @endreturndesc
 @@*/
 static int PUGH_ReductionScalar (cGH *GH, 
@@ -395,26 +496,45 @@ static int PUGH_ReductionScalar (cGH *GH,
 
  /*@@
    @routine    copy_real_to_outtype
-   @date       
-   @author     
+   @author     Thomas Radke
+   @date       19 Aug 1999
    @desc 
+               Does the type conversion from CCTK_REAL into the requested
+               datatype.
+   @enddesc
+   @calls      CCTK_VarTypeSize
+               CCTK_GroupTypeFromVarI
+               PUGH_ReductionGA
+               PUGH_ReductionScalar
+               copy_real_to_outtype
+  
+   @var        num_elems
+   @vdesc      number of elements to convert
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        inarray
+   @vdesc      input array with results of reductions
+   @vtype      CCTK_REAL *
+   @vio        in
+   @endvar
+   @var        outtype
+   @vdesc      requested output datatype
+   @vtype      int
+   @vio        in
+   @endvar
+   @var        outarray
+   @vdesc      pointer to output buffer to store the converted results
+   @vtype      void *
+   @vio        out
+   @endvar
 
-   @enddesc 
-   @calls      
-   @calledby   
-   @history 
- 
-   @endhistory 
-   @var     
-   @vdesc   
-   @vtype   
-   @vio     
-   @vcomment 
- 
-   @endvar 
-
+   @returntype int
+   @returndesc
+               0  for success, or<BR>
+               -1 if conversion into target datatype is not supported.
+   @endreturndesc
 @@*/
-
 static int copy_real_to_outtype (int num_elems,
                                  CCTK_REAL inarray[/* num_elems */],
                                  int outtype,
@@ -506,13 +626,13 @@ static int copy_real_to_outtype (int num_elems,
 #ifdef CCTK_REAL8
   else if (outtype == CCTK_VARIABLE_REAL8)
   {
-	CCTK_REAL8 *_outarray = (CCTK_REAL8 *) outarray;
+        CCTK_REAL8 *_outarray = (CCTK_REAL8 *) outarray;
 
 
-	for (i = 0; i < num_elems; i++)
-	{
-	  _outarray[i] = (CCTK_REAL8) inarray[i];
-	}
+        for (i = 0; i < num_elems; i++)
+        {
+          _outarray[i] = (CCTK_REAL8) inarray[i];
+        }
   }
 #endif
 #ifdef CCTK_REAL16