fixed staggering with new global reduction, will add more tests to globalreduction

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

1 files changed, 45 insertions, 14 deletions

diff --git a/src/ReduceGA.c b/src/ReduceGA.c
index 4de66f1..045db2f 100644
--- a/src/ReduceGA.c
+++ b/src/ReduceGA.c
@@ -135,13 +135,16 @@ static int ReduceGridArrays (const cGH *GH,
   int total_num_points = 1;
   int perform_division = 1;
   int perform_all_reduce = 1;
+  pGA *GA;
+  const void *data;
 
   CCTK_INT * lower_array_bounds;
-  CCTK_INT * upper_array_bounds;
+  CCTK_INT * min_array_subscript;
+  CCTK_INT * array_lsh;
   CCTK_INT * input_array_dims;
   CCTK_POINTER * input_arrays;
   CCTK_INT * input_array_type_codes;
-
+  CCTK_INT * input_array_gz;
 
 #ifdef CCTK_MPI
   int nprocs = 0, myproc =0, global_operation = 0;
@@ -176,9 +179,11 @@ static int ReduceGridArrays (const cGH *GH,
   dim = CCTK_GroupDimFromVarI(input_array_variable_indices[0]);
  
   /* allocate memory for the array of the dimensions of the input arrays */
-  lower_array_bounds = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
-  upper_array_bounds = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
-  input_array_dims = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
+  lower_array_bounds  = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
+  min_array_subscript = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
+  array_lsh           = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
+  input_array_dims    = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
+  input_array_gz      = (CCTK_INT *)malloc (dim *sizeof(CCTK_INT));
 
   /* find out the types of the input arrays and put that */
   /* in an array */
@@ -188,37 +193,63 @@ static int ReduceGridArrays (const cGH *GH,
     input_arrays[i]           = CCTK_VarDataPtrI(GH, 0, input_array_variable_indices[i] );
   }
 
-  /* fill out the dimensions of the input arrays */
+  /* fill out the lsh of the input arrays, the lower array bounds and ghostzone numbers*/
+  ierr = CCTK_GrouplshVI(GH, dim, array_lsh, input_array_variable_indices[0] );
   ierr = CCTK_GrouplbndVI(GH, dim, lower_array_bounds, input_array_variable_indices[0] );
-  ierr = CCTK_GroupubndVI(GH, dim, upper_array_bounds, input_array_variable_indices[0] );
+  i    = CCTK_GroupIndexFromVarI(input_array_variable_indices[0]);
+  ierr = CCTK_GroupnghostzonesGI(GH, dim, input_array_gz,i );
 
-  /* incement max subscript to use < instead of <= */
+  /* set lower bounds to zero and increment max subscript to use < instead of <= */
   for ( i=0; i < dim; i++)
   {
-    upper_array_bounds[i]++;
+    min_array_subscript[i] = 0;
   }
+
+
+  /* Fix subscripts for staggering */
+  /* get the GA structure for the current timelevel */
+  GA   = ((pGA ***) PUGH_pGH (GH)->variables)[input_array_variable_indices[0]][0];
+  data = GA->data;
+
+  for ( i=0; i < dim; i++)
+  {
+    /* get the start- and endpoint to iterate over in this dimension */
+    min_array_subscript[i] = GA->extras->ownership[GA->stagger[i]][0][i];
+    array_lsh[i]           = GA->extras->ownership[GA->stagger[i]][1][i];
+
+    /* get the number of points in this dimension (excluding ghostzones) */
+    if (GA->connectivity->perme[i])
+    {
+      array_lsh[i] -= 2*GA->extras->nghostzones[i];
+    }
+    if (GA->stagger[i])
+    {
+      array_lsh[i]--;
+    }
+  }
+
   /* Create the parameter table if it is not there and add the bounds to it */
   if ( Util_TableQueryNKeys(param_table_handle) == 0)
   {
     param_table_handle = Util_TableCreate (UTIL_TABLE_FLAGS_DEFAULT);
-    ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,lower_array_bounds, "input_array_min_subscripts"); 
-    ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,upper_array_bounds, "input_array_max_subscripts"); 
+    ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,min_array_subscript, "input_array_min_subscripts"); 
+    ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,array_lsh, "input_array_max_subscripts"); 
   }
   else
   {
     if (!Util_TableQueryValueInfo(param_table_handle, NULL, NULL, "input_array_min_subscripts"))
     {
-      ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,lower_array_bounds, "input_array_min_subscripts"); 
+      ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,min_array_subscript, "input_array_min_subscripts"); 
     }
     if (!Util_TableQueryValueInfo(param_table_handle, NULL, NULL, "input_array_max_subscripts"))
     {
-      ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,lower_array_bounds, "input_array_max_subscripts"); 
+      ierr = Util_TableSetGenericArray (param_table_handle, CCTK_VARIABLE_INT, dim,array_lsh, "input_array_max_subscripts"); 
     }
   }
 
   for (i = 0; i<dim; i++)
   {
-    input_array_dims[i] = upper_array_bounds[i] - lower_array_bounds[i];
+    input_array_dims[i] = array_lsh[i] - min_array_subscript[i];
   }
 
   /* Set flag to tell local reduction not to divide by num_points */