Fixed the case where an error code CCTK_ERROR_INTERP_POINT_OUTSIDE, as returned

by the local interpolator query call, should be translated by PUGHInterp into
CCTK_ERROR_INTERP_GHOST_SIZE_TOO_SMALL.
Also Jonathan Thornburg added some comments on how the query call interpolation
points are set up.


git-svn-id: http://svn.cactuscode.org/arrangements/CactusPUGH/PUGHInterp/trunk@56 1c20744c-e24a-42ec-9533-f5004cb800e5

1 files changed, 28 insertions, 9 deletions

diff --git a/src/InterpGridArrays.c b/src/InterpGridArrays.c
index 1fb3c1d..d3f09a7 100644
--- a/src/InterpGridArrays.c
+++ b/src/InterpGridArrays.c
@@ -565,21 +565,39 @@ int PUGHInterp_InterpGridArrays (const cGH *GH,
      interpolator to handle points on the processor's local bounding box.
      This is done as a query, simply by passing all arguments as for the
      real call to the local interpolator, except for the output_arrays[]
-     argument which contains NULL pointers. */
+     argument which contains NULL pointers.
+
+     To do this check, we want to set up (up to) 2*N_dims test points, one
+     in the center of each interprocessor face of this processor's local
+     bounding box; for non-interprocessor faces (i.e. the physical outer
+     boundaries of the grid) we don't want to make any test (since there's
+     no ghost zone there).
+
+     To simplify the bookkeeping, we actually always set up 2*N_dims points.
+     We initialize each point to be at the center of the bounding box, then
+     move the interprocessor-face points to the centers of their faces (i.e.
+     we move each interprocessor-face points to its bounding-box coordinate
+     in the direction perpendicular to its face.  The non-interprocessor-face
+     points are left at the center of the bounding box, where their interpolator
+     stencils will be "safe" and (we can be sure) never fall off the edge of
+     the bounding box.
+   */
   bbox_interp_coords = malloc (N_dims * sizeof (CCTK_REAL *));
   bbox_interp_coords[0] = malloc (2 * N_dims * N_dims * sizeof (CCTK_REAL));
   for (i = 1; i < N_dims; i++)
   {
     bbox_interp_coords[i] = bbox_interp_coords[i-1] + 2*N_dims;
   }
-  for (j = 0; j < 2 * N_dims; j++)
+  for (j = 0; j < 2 * N_dims; j++)	/* j indices faces of the bounding box */
   {
+    /* initialize this point to the middle of the bounding box */
     for (i = 0; i < N_dims; i++)
     {
       bbox_interp_coords[i][j] = 0.5 * (bbox_local[2*i+0] + bbox_local[2*i+1]);
     }
-    if (! GH->cctk_bbox[j])
+    if (! GH->cctk_bbox[j])		/* if this is an interprocessor face... */
     {
+      /* ... move the point to the center of this face of the bounding box */
       bbox_interp_coords[j/2][j] = bbox_local[j];
     }
   }
@@ -596,6 +614,13 @@ int PUGHInterp_InterpGridArrays (const cGH *GH,
                                     N_output_arrays, output_array_types,
                                     (void *const *) output_arrays_local);
 
+  /* if the query call returned 'point-outside' error
+     we turn this into a 'ghost-size-too-small' error here */
+  if (retval == CCTK_ERROR_INTERP_POINT_OUTSIDE)
+  {
+    retval = CCTK_ERROR_INTERP_GHOST_SIZE_TOO_SMALL;
+  }
+
   /* make sure that all points could be mapped onto a processor */
   if (retval == 0)
   {
@@ -640,12 +665,6 @@ int PUGHInterp_InterpGridArrays (const cGH *GH,
     free (input_array_dims);
     free (type_desc);
 
-    /* if the query call returned 'point-outside' error
-       we turn this into a 'ghost-size-too-small' error here */
-    if (retval == CCTK_ERROR_INTERP_POINT_OUTSIDE)
-    {
-      retval = CCTK_ERROR_INTERP_GHOST_SIZE_TOO_SMALL;
-    }
     return (retval);
   }