Improve periodic boundaries on small domains: allow the number of

ghost zones to be larger than the number of interior points if there
is only one processor in the corresponding direction.  This is
necessary to allow two ghost zones when there is only one interior
grid point, which is often used to set up lower-dimensional domains.

Second part of this patch.


git-svn-id: http://svn.cactuscode.org/arrangements/CactusPUGH/PUGH/trunk@470 b61c5cb5-eaca-4651-9a7a-d64986f99364

1 files changed, 91 insertions, 71 deletions

diff --git a/src/Comm.c b/src/Comm.c
index 3d04ae1..f739f6e 100644
--- a/src/Comm.c
+++ b/src/Comm.c
@@ -1005,95 +1005,115 @@ static int PUGH_SyncSingleProc(pGH *pughGH, pComm *comm)
     /* there's only something to do for periodic boundary conditions */
     if (GA->connectivity->perme[face / 2])
     {
-      /* get the index ranges for the nested loops */
-      istart_from   = GA->extras->overlap[GA->stagger[face>>1]][0][face];
-      iend_from     = GA->extras->overlap[GA->stagger[face>>1]][1][face];
-      if (face & 1)
+      /* if the domain is smaller than the number of ghost zones, then
+         we have to copy in several passes */
+      /* FIXME: these passes do not each need to copy the whole
+         regaion, they could copy a part of it only */
+      /* total number of points to copy */
+      pass_total = GA->extras->nghostzones[face>>1];
+      /* points copied in one pass */
+      pass_each = (GA->extras->lnsize[face>>1]
+                   - 2 * GA->extras->nghostzones[face>>1]);
+      /* necessary number of passes */
+      num_passes = (pass_total + pass_each - 1) / pass_each;
+      for (pass = 0; pass < num_passes; pass++)
       {
-        istart_to   = GA->extras->ghosts[GA->stagger[face>>1]][0][face - 1];
-      }
-      else
-      {
-        istart_to   = GA->extras->ghosts[GA->stagger[face>>1]][0][face + 1];
-      }
+        /* calculate what needs to be done in this pass */
+        /* FIXME: these values are currently ignored */
+        pass_from = pass * pass_each;
+        pass_to = pass_from + pass_each - 1;
+        if (pass_to > pass_total) pass_to = pass_total;
+
+        /* get the index ranges for the nested loops */
+        istart_from   = GA->extras->overlap[GA->stagger[face>>1]][0][face];
+        iend_from     = GA->extras->overlap[GA->stagger[face>>1]][1][face];
+        if (face & 1)
+        {
+          istart_to   = GA->extras->ghosts[GA->stagger[face>>1]][0][face - 1];
+        }
+        else
+        {
+          istart_to   = GA->extras->ghosts[GA->stagger[face>>1]][0][face + 1];
+        }
 
-      /* set iterators to the start vector */
-      for(dim = 0; dim < GA->extras->dim; dim++)
-      {
-        iterator_from[dim] = istart_from[dim];
-        iterator_to[dim] = istart_to[dim];
-      }
+        /* set iterators to the start vector */
+        for(dim = 0; dim < GA->extras->dim; dim++)
+        {
+          iterator_from[dim] = istart_from[dim];
+          iterator_to[dim] = istart_to[dim];
+        }
 
-      /* get the number of bytes to copy in the lowest dimension */
-      copy_bytes = (iend_from[0] - istart_from[0]) * GA->varsize;
+        /* get the number of bytes to copy in the lowest dimension */
+        copy_bytes = (iend_from[0] - istart_from[0]) * GA->varsize;
 
-      /* now do the nested loops starting with the innermost */
-      dim = 1;
-      while (1)
-      {
-        /* check for end of current loop */
-        if (GA->extras->dim > 1 && iterator_from[dim] >= iend_from[dim])
+        /* now do the nested loops starting with the innermost */
+        dim = 1;
+        while (1)
         {
-          /* increment outermost loopers */
-          for (dim++; dim < GA->extras->dim; dim++)
+          /* check for end of current loop */
+          if (GA->extras->dim > 1 && iterator_from[dim] >= iend_from[dim])
           {
-            iterator_from[dim]++;
-            iterator_to[dim]++;
-            if (iterator_from[dim] < iend_from[dim])
+            /* increment outermost loopers */
+            for (dim++; dim < GA->extras->dim; dim++)
+            {
+              iterator_from[dim]++;
+              iterator_to[dim]++;
+              if (iterator_from[dim] < iend_from[dim])
+              {
+                break;
+              }
+            }
+
+            /* done if beyond outermost loop */
+            if (dim >= GA->extras->dim)
             {
               break;
             }
+
+            /* reset innermost loopers */
+            for (dim--; dim > 0; dim--)
+            {
+              iterator_from[dim] = istart_from[dim];
+              iterator_to[dim] = istart_to[dim];
+            }
+            dim = 1;
           }
 
-          /* done if beyond outermost loop */
-          if (dim >= GA->extras->dim)
+          /* get the linear offsets */
+          offset_from = istart_from[0];
+          offset_to   = istart_to[0];
+          for (i = 1; i < GA->extras->dim; i++)
           {
-            break;
+            offset_from += iterator_from[i] * GA->extras->hyper_volume[i];
+            offset_to   += iterator_to[i]   * GA->extras->hyper_volume[i];
           }
+          offset_from *= GA->varsize;
+          offset_to   *= GA->varsize;
 
-          /* reset innermost loopers */
-          for (dim--; dim > 0; dim--)
+          /* copy the data for all the variables in the comm structure */
+          for (i = comm->first_var; i < comm->first_var + comm->n_vars; i++)
           {
-            iterator_from[dim] = istart_from[dim];
-            iterator_to[dim] = istart_to[dim];
+            data = ((pGA *) pughGH->variables[i][comm->sync_timelevel])->data;
+            memcpy (data + offset_to, data + offset_from, copy_bytes);
           }
-          dim = 1;
-        }
-
-        /* get the linear offsets */
-        offset_from = istart_from[0];
-        offset_to   = istart_to[0];
-        for (i = 1; i < GA->extras->dim; i++)
-        {
-          offset_from += iterator_from[i] * GA->extras->hyper_volume[i];
-          offset_to   += iterator_to[i]   * GA->extras->hyper_volume[i];
-        }
-        offset_from *= GA->varsize;
-        offset_to   *= GA->varsize;
 
-        /* copy the data for all the variables in the comm structure */
-        for (i = comm->first_var; i < comm->first_var + comm->n_vars; i++)
-        {
-          data = ((pGA *) pughGH->variables[i][comm->sync_timelevel])->data;
-          memcpy (data + offset_to, data + offset_from, copy_bytes);
-        }
-
-        /* increment current loopers */
-        if (GA->extras->dim > 1)
-        {
-          /* increment current looper */
-          iterator_from[dim]++;
-          iterator_to[dim]++;
-        }
-        else
-        {
-          /* exit loop if array dim is 1D */
-          break;
-        }
+          /* increment current loopers */
+          if (GA->extras->dim > 1)
+          {
+            /* increment current looper */
+            iterator_from[dim]++;
+            iterator_to[dim]++;
+          }
+          else
+          {
+            /* exit loop if array dim is 1D */
+            break;
+          }
 
-      } /* end of nested loops over all dimensions */
-    }
-  }
+        } /* end of nested loops over all dimensions */
+      } /* for pass */
+    } /* if periodic */
+  } /* for face */
 
   /* free the allocated iterator */
   free (iterator_to);