A few more routines fleshed out.

Tom


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

1 files changed, 748 insertions, 13 deletions

diff --git a/src/SetupPGV.c b/src/SetupPGV.c
index 213855a..8b1c676 100644
--- a/src/SetupPGV.c
+++ b/src/SetupPGV.c
@@ -13,6 +13,18 @@ static char *rcsid = "$Header$";
 #include <stdlib.h>
 #include <math.h>
 
+/* Stuff to do standalone test */
+#define TEST_SETUPPGV
+
+#ifdef TEST_SETUPPGV
+#define PUGH_NSTAGGER 4
+#define PUGH_NO_STAGGER 0
+#define PUGH_STAGGER 1
+#define PUGH_VERTEXCTR 0
+#define PUGH_FACECTR(i) i+1
+#define CCTK_REAL double
+#endif /* TEST_SETUPPGV */
+
 typedef enum {pgv_none, pgv_scalar, pgv_array, pgv_gf} pgv_type;
 
 typedef struct PConnectivity
@@ -20,6 +32,7 @@ typedef struct PConnectivity
   int dim;
   int *nprocs;
   int **neighbours;
+  int periodic;		/* Is the system periodic? */
 } pConnectivity;
 
 typedef struct PGS
@@ -31,10 +44,12 @@ typedef struct PGS
 typedef struct PGExtras
 {
   int   dim;                    /* dimension of GF */
+
+  int *nsize;                  /* The global size of the array */
+
   /* Processor group layouts */
   CCTK_REAL maxskew;		/* Maximum point skew */
-  int       periodic;		/* Is the system periodic? */
-  int     **lb;			/* Lower bound (nprocs X 3) for each proc */
+  int     **lb;			/* Lower bound (nprocs X dim) for each proc */
   int     **ub;			/* Upper bound (same sizes) */
   int      *lnsize;		/* Size on this processor */
   int       npoints;		/* LOCAL number of points on this proc. */
@@ -42,6 +57,7 @@ typedef struct PGExtras
   int     **rnsize;	        /* [#points on a proc][in each dir] */
 
   /* Ghosts and overlaps. */
+  int      *nghostzones;	/* Width of ghost zone */
   int      *ownership[PUGH_NSTAGGER][2];	
                                 /* The box owned in each direction. */
                                 /* [stagger][min/max][dir] */
@@ -88,7 +104,7 @@ typedef struct PGA
   MPI_Status ms;
 #endif
 
-  pGExtra *extras;
+  pGExtras *extras;
 } pGA;
 
 typedef struct PGV
@@ -101,17 +117,34 @@ typedef struct PGV
 
 pGExtras *pugh_SetupPGExtras(int dim, 
                              int periodic,
+                             int staggertype,
                              int *sh,
                              int *nghosts,
-                             int *nprocs)
+                             int total_procs,
+                             int *nprocs,
+                             int this_proc)
 {
+  int error;
   pGExtras *this;
 
   this = (pGExtras *)malloc(sizeof(pGExtras));
 
+  /* Setup memory */
   if(this)
   {
-    pgExtras_SetupBasics();
+    error = pugh_SetupPGExtrasMemory(dim, 
+                                     total_procs,
+                                     nprocs,
+                                     this);
+
+    if(!error)
+    {
+      pugh_SetupPGExtrasSizes   (dim, periodic, staggertype, sh, nghosts, 
+                                 total_procs, nprocs, this_proc,this);
+      pugh_SetupPGExtrasOwnership(dim, periodic, staggertype, sh, nghosts, 
+                                 total_procs, nprocs, this_proc, this);
+    }      
+
   }
 
   return this;
@@ -135,14 +168,15 @@ pGExtras *pugh_SetupPGExtras(int dim,
 @@*/
 pConnectivity *pugh_SetupConnectivity(int dim, 
                                       int total_procs, 
-                                      int *nprocs)
+                                      int *nprocs,
+                                      int periodic)
 {
   pConnectivity *this;
 
   int i;
 
   /* Allocate memory */
-  this = (pConnectivity *)malloc(pConnectivity);
+  this = (pConnectivity *)malloc(sizeof(pConnectivity));
 
   if(this)
   {
@@ -151,15 +185,16 @@ pConnectivity *pugh_SetupConnectivity(int dim,
 
     if(this->nprocs && this->neighbours)
     {
-      this->neighbours[0] = (int *)malloc(nprocs*2*dim*sizeof(int));
+      this->neighbours[0] = (int *)malloc(total_procs*2*dim*sizeof(int));
 
       if(this->neighbours[0])
       {
-        for(i = 1; i < nprocs; i++)
+        for(i = 1; i < total_procs; i++)
         {
           this->neighbours[i] = this->neighbours[0]+(2*dim*i);
         }
 
+        this->periodic = periodic;
       }
       else
       {
@@ -188,7 +223,7 @@ pConnectivity *pugh_SetupConnectivity(int dim,
 
     pugh_GenerateTopology(dim, total_procs, this->nprocs);
 
-    pugh_GenerateNeighbours(dim, total_procs, this->nprocs, this->neighbours);
+    pugh_GenerateNeighbours(dim, total_procs, this->nprocs, this->neighbours, this->periodic);
 
   }
 
@@ -279,12 +314,28 @@ int pugh_GenerateTopology(int dim, int total_procs, int *nprocs)
   return retval;
 }
 
-int pugh_GenerateNeighbours(int dim, int total_procs, int *nprocs, int **neighbours)
+ /*@@
+   @routine    pugh_GenerateNeighbours
+   @date       Mon Nov  8 08:15:08 1999
+   @author     Tom Goodale
+   @desc 
+   Works out the array of neighbouring processors for
+   every processor.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+int pugh_GenerateNeighbours(int dim, int total_procs, int *nprocs, int **neighbours, int periodic)
 {
   int retval;
   int i;
   int idim;
   int *pos;
+  int temp;
   
   pos = (int *)malloc(dim*sizeof(int));
 
@@ -292,7 +343,7 @@ int pugh_GenerateNeighbours(int dim, int total_procs, int *nprocs, int **neighbo
   {
     for(i = 0; i < total_procs; i++)
     {
-      pugh_DecomposeIJK(i,nprocs, pos);
+      pugh_DecomposeIJK(dim,i,nprocs, pos);
       
       for(idim = 0; idim < dim ; idim++)
       {
@@ -301,7 +352,14 @@ int pugh_GenerateNeighbours(int dim, int total_procs, int *nprocs, int **neighbo
         
         if(pos[idim] > -1)
         {
-          neighbours[i][idim*2] = pugh_ComposeIJK(dim, i, nprocs, pos);
+          neighbours[i][idim*2] = pugh_ComposeIJK(dim, nprocs, pos);
+        }
+        else if(periodic)
+        {
+          temp = pos[idim];
+          pos[idim] = nprocs[idim]-1;
+          neighbours[i][idim*2] = pugh_ComposeIJK(dim, nprocs, pos);
+          pos[idim] = temp;
         }
 
         pos[idim]++;
@@ -313,6 +371,13 @@ int pugh_GenerateNeighbours(int dim, int total_procs, int *nprocs, int **neighbo
         {
           neighbours[i][idim*2+1] = pugh_ComposeIJK(dim, nprocs, pos);
         }
+        else if(periodic)
+        {
+          temp = pos[idim];
+          pos[idim] = 0;
+          neighbours[i][idim*2+1] = pugh_ComposeIJK(dim, nprocs, pos);
+          pos[idim] = temp;
+        }
         
         pos[idim]--;
       }
@@ -325,6 +390,8 @@ int pugh_GenerateNeighbours(int dim, int total_procs, int *nprocs, int **neighbo
     retval = 1;
   }
 
+  free(pos);
+
   return retval;
 }
   
@@ -401,3 +468,671 @@ int pugh_ComposeIJK(int dim,
 
   return ijk;
 }
+
+
+ /*@@
+   @routine    pugh_SetupPGExtrasMemory
+   @date       Mon Nov  8 08:16:02 1999
+   @author     Tom Goodale
+   @desc 
+   Allocate memory for the members of the pGExtras structure.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+int pugh_SetupPGExtrasMemory(int dim, 
+                             int total_procs,
+                             int *nprocs,
+                             pGExtras *this)
+{
+  int retcode;
+  int i,j,k;
+
+  retcode = 0;
+
+  if(this)
+  {
+    /* Do it in stages.
+     * First: things depending on the number of processors 
+     */
+    this->lb       = (int **)malloc(total_procs*sizeof(int *));
+    this->ub       = (int **)malloc(total_procs*sizeof(int *));
+    this->rnsize   = (int **)malloc(total_procs*sizeof(int *));
+    this->rnpoints = (int *) malloc(total_procs*sizeof(int));
+
+    /* Things just depending on dimension */
+    this->nghostzones = (int *)malloc(dim*sizeof(int));
+    this->nsize       = (int *)malloc(dim*sizeof(int));
+    this->lnsize      = (int *)malloc(dim*sizeof(int));
+
+    /* Check all the above succeeded and then get memory for 
+     * arrays hanging off the above. 
+     */
+    if(this->lb          &&
+       this->ub          &&
+       this->rnsize      &&
+       this->rnpoints    &&
+       this->nghostzones &&
+       this->nsize       &&
+       this->lnsize)
+    {
+      this->lb[0]     = (int *)malloc(total_procs  *dim*sizeof(int));
+      this->ub[0]     = (int *)malloc(total_procs  *dim*sizeof(int));
+      this->rnsize[0] = (int *)malloc(total_procs*2*dim*sizeof(int));
+
+      if(this->lb[0]     &&
+         this->ub[0]     &&
+         this->rnsize[0])
+      {
+        for (i = 1; i < total_procs; i++)
+        {
+          this->lb[i]     = this->lb[0]     + i*dim;
+          this->ub[i]     = this->ub[0]     + i*dim;
+          this->rnsize[i] = this->rnsize[0] + i*2*dim;
+        }
+      }
+      else
+      {
+        /* Free inner arrays */
+        free(this->lb[0]);
+        this->lb[0] = NULL;
+        free(this->ub[0]);
+        this->ub[0] = NULL;
+        free(this->rnsize[0]);
+        this->rnsize[0] = NULL;
+
+        /* Free toplevel arrays */
+        free(this->lb);
+        this->lb = NULL;
+        free(this->ub);
+        this->ub = NULL;
+        free(this->rnsize);
+        this->rnsize = NULL;
+        free(this->rnpoints);
+        this->rnpoints = NULL;
+        free(this->nghostzones);
+        this->nghostzones = NULL;               
+        free(this->nsize);
+        this->nsize = NULL;               
+        free(this->lnsize);
+        this->lnsize = NULL;               
+      }
+    }
+    else
+    {
+      /* Free toplevel arrays */
+      free(this->lb);
+      this->lb = NULL;
+      free(this->ub);
+      this->ub = NULL;
+      free(this->rnsize);
+      this->rnsize = NULL;
+      free(this->rnpoints);
+      this->rnpoints = NULL;
+      free(this->nghostzones);
+      this->nghostzones = NULL;               
+      free(this->nsize);
+      this->nsize = NULL;               
+      free(this->lnsize);
+      this->lnsize = NULL;               
+    }
+      
+
+    if(this->lb          &&
+       this->ub          &&
+       this->rnsize      &&
+       this->rnpoints    &&
+       this->nghostzones &&
+       this->nsize       &&
+       this->lnsize)
+    {
+      retcode = 0;
+      for (i = 0 ; i < PUGH_NSTAGGER; i++)
+      {
+        for (j = 0; j < 2; j++)
+        {
+          this->ownership[i][j]    = (int *) malloc(dim*sizeof(int));
+          this->ghosts[i][j]       = (int **)malloc(2*dim*sizeof(int *));
+          this->overlap[i][j]      = (int **)malloc(2*dim*sizeof(int *));
+          if(this->ghosts[i][j] &&
+             this->overlap[i][j])
+          {
+            this->ghosts[i][j][0]  = (int *)malloc(2*dim*dim*sizeof(int));
+            this->overlap[i][j][0] = (int *)malloc(2*dim*dim*sizeof(int));
+            for (k=1; k < 2*dim; k++)
+            {
+              this->ghosts[i][j][k]  = this->ghosts[i][j][0] + k*dim;
+              this->overlap[i][j][k] = this->overlap[i][j][0] + k*dim;
+            }
+          }
+          else
+          {
+            free(this->ownership[i][j]);
+            this->ownership[i][j] = NULL;
+
+            free(this->ghosts[i][j]);
+            this->ghosts[i][j] = NULL;
+
+            free(this->overlap[i][j]);
+            this->overlap[i][j] = NULL;
+            retcode = 1;
+            break;
+          }          
+        }
+        if(retcode)
+        {
+          for(j=1; j >=0 ; j--)
+          {
+              free(this->ownership[i][j]);
+              this->ownership[i][j] = NULL;
+              
+              free(this->ghosts[i][j]);
+              this->ghosts[i][j] = NULL;
+              
+              free(this->overlap[i][j]);
+              this->overlap[i][j] = NULL;
+              
+          }
+          break;
+        }
+      }
+      if(retcode)
+      {
+        /* Loop back through the arrays freeing things */
+        for(i--; i >=0; i--)
+        {
+          for(j=1; j >=0 ; j--)
+          {
+            free(this->ghosts[i][j][0]);
+            free(this->overlap[i][j][0]);
+
+            free(this->ownership[i][j]);
+            this->ownership[i][j] = NULL;
+            
+            free(this->ghosts[i][j]);
+            this->ghosts[i][j] = NULL;
+          
+            free(this->overlap[i][j]);
+            this->overlap[i][j] = NULL;
+          }
+
+          free(this->ownership[i]);
+          free(this->ghosts[i]);
+          free(this->overlap[i]);
+        }
+
+        /* Free the stuff originally allocated */
+
+        free(this->lb[0]);
+        this->lb[0] = NULL;
+        free(this->ub[0]);
+        this->ub[0] = NULL;
+        free(this->rnsize[0]);
+        this->rnsize[0] = NULL;
+
+        free(this->lb);
+        this->lb = NULL;
+        free(this->ub);
+        this->ub = NULL;
+        free(this->rnsize);
+        this->rnsize = NULL;
+        free(this->rnpoints);
+        this->rnpoints = NULL;
+        free(this->nghostzones);
+        this->nghostzones = NULL;               
+        free(this->nsize);
+        this->nsize = NULL;               
+        free(this->lnsize);
+        this->lnsize = NULL;               
+      }
+    }          
+  }
+  else
+  {
+    retcode = -1;
+  }
+
+  return retcode;
+}
+
+ /*@@
+   @routine    pugh_SetupPGExtrasSizes
+   @date       Mon Nov  8 08:59:33 1999
+   @author     Tom Goodale
+   @desc 
+   Sets up the size information in the pGExtras 
+   structure.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+int pugh_SetupPGExtrasSizes(int dim, 
+                            int periodic,
+                            int staggertype,
+                            int *sh,
+                            int *nghosts,
+                            int total_procs,
+                            int *nprocs,
+                            int this_proc,
+                            pGExtras *this)
+{
+  int dir;
+
+  /* Setup the global shape */
+  for (dir=0 ; dir < dim ; dir++)
+  {
+    /* A -ve size means constant load per proc */
+    if (sh[dir] < 0 && nprocs[dir] > 1) 
+    {
+      this->nsize[dir] = (nprocs[dir]-2) * 
+	(-sh[dir] - 2*nghosts[dir]) +
+	2 * (-sh[dir] - nghosts[dir]);
+
+      if (staggertype == PUGH_STAGGER) 
+      {
+	this->nsize[dir] -= nprocs[dir]-1;
+      }
+    } 
+    else 
+    {
+      this->nsize[dir] = abs(sh[dir]);
+    }
+  }
+
+  /* Set the number of ghost zones. */
+  for (dir = 0; dir < dim; dir++)
+  {
+    this->nghostzones[dir] = nghosts[dir];
+  }
+
+  /* Setup the bounding box stuff */
+
+  pugh_SetupBoundingBox(dim, periodic, staggertype, sh, nghosts, total_procs, nprocs, this); 
+
+  /* Set the remote sizes */
+
+  pugh_SetupRemoteSizes(dim, periodic, staggertype, sh, nghosts, total_procs, nprocs, this); 
+  
+  /* Set the local sizes */
+
+  for (dir = 0; dir < dim; dir++)
+  {
+    this->lnsize[dir] = this->rnsize[this_proc][dir];
+  }
+
+
+  return 0;
+}
+     
+
+ /*@@
+   @routine    pugh_SetupPGExtrasOwnership
+   @date       Mon Nov  8 09:00:10 1999
+   @author     Tom Goodale
+   @desc 
+   Sets up ownership, overlap, ghostzones, etc 
+   in a pGExtras structure.
+
+   Mostly taken from original SetupOwnership by Paul.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+int pugh_SetupPGExtrasOwnership(int dim, 
+                               int periodic,
+                               int staggertype,
+                               int *sh,
+                               int *nghosts,
+                               int total_procs,
+                               int *nprocs,
+                               int this_proc,
+                               pGExtras *this)
+{
+  int tmp;
+  int dir, idir;
+  int i,j,k;
+  int istart, iend;
+  int staglcv;
+
+  /* Ownership is pretty easy. Remember ownership is indexed as
+     [stagger][ijk][min/max]. See pGF_Reduction for a use of this,
+     among others.
+     Note: Ownership is same for staggered and non-staggered grids
+   */
+  for (dir = 0 ; dir < dim; dir++)
+  {
+    this->ownership[PUGH_VERTEXCTR][0][dir] = (this->lb[this_proc][dir] == 0 ? 
+                                               0 : this->nghostzones[dir]);
+    this->ownership[PUGH_VERTEXCTR][1][dir]=(this->ub[this_proc][dir] == this->nsize[dir]-1 ? 
+                                             this->lnsize[dir] : this->lnsize[dir] - 
+                                             this->nghostzones[dir]);
+  }
+
+  /* correct for periodic identification : Tue Jun 17 08:40:15 CDT 1997 */
+  if(periodic)
+  {
+    for (dir = 0; dir < dim; dir++)
+    {
+      this->ownership[PUGH_VERTEXCTR][0][dir] = this->nghostzones[dir];
+      this->ownership[PUGH_VERTEXCTR][1][dir] = this->lnsize[dir] - this->nghostzones[dir];
+    }
+  } 
+
+  for (dir = 0; dir < dim; dir++)
+  {
+    if (this->nsize[dir] == 1) 
+    {
+      this->ownership[PUGH_VERTEXCTR][0][dir] = 0;
+      this->ownership[PUGH_VERTEXCTR][1][dir] = 1;
+    }
+  }
+
+
+  /* Great now set up the overlap region. This is the region we own,
+     but which is covered by our neighbors ghost. So this is what we
+     want to SEND
+   */
+  for (dir = 0; dir < 2*dim; dir++) 
+  {
+    for (idir = 0; idir < dim; idir++)
+    {
+      istart = 0; 
+      iend   = this->lnsize[idir];
+
+      /* We want to send at the dir - 1 */
+
+      /* Correct to direction specific sweeps ... */
+      if (dir == 2*idir) 
+      {
+	/* Want to go from sw to sw + sw (eg, 1 to < 2, 2 to < 4) */
+	istart = this->nghostzones[idir];
+	iend   = istart + this->nghostzones[idir];
+      }
+      if (dir == 2*idir+1) {
+	/* Want to go from nx - 2 sw to nx - sw
+	   (eg nx-2 to < nx-1, nx-4 to < nx-2)
+	   Luckily iend is aready set to nx ... */
+	tmp = iend;
+	istart = tmp - 2 * this->nghostzones[idir];
+	iend   = tmp -     this->nghostzones[idir];
+      }
+
+      this->overlap[PUGH_VERTEXCTR][0][dir][idir] = istart;
+      this->overlap[PUGH_VERTEXCTR][1][dir][idir] = iend;
+    }
+  }
+
+  /* And finally the ghost indices. This used to be in
+     pGF_FinishRecv.c, but now it is here.
+   */
+  for (dir = 0; dir < 2*dim; dir++) 
+  {
+    /* We want to send at the dir - 1 */
+    for (idir = 0; idir < dim; idir++)
+    {
+      istart = 0; 
+      iend = this->lnsize[idir];
+
+      /* Correct to direction specific sweeps ... */
+      /* Remember this is SW less than the other way... */
+      if (dir == 2*idir) 
+      {
+      /* Want to go from sw to sw + sw (eg, 1 to < 2, 2 to < 4) */
+	istart = 0;
+	iend   = this->nghostzones[idir];
+      }
+      if (dir == 2*idir+1) 
+      {
+      /* Want to go from nx - 2 sw to nx - sw
+         (eg nx-2 to < nx-1, nx-4 to < nx-2)
+         Luckily iend is aready set to nx ... */
+	tmp = iend;
+	istart = tmp - this->nghostzones[idir];
+	iend   = tmp;
+      }
+
+      this->ghosts[PUGH_VERTEXCTR][0][dir][idir] = istart;
+      this->ghosts[PUGH_VERTEXCTR][1][dir][idir] = iend;
+
+    }
+
+  }
+
+  /* if staggering, set up ownership, overlaps, ghosts */
+  if(staggertype == PUGH_STAGGER) 
+  {
+    /* copy ownership from PUGH_VERTEXCTR */
+    for(i = 1; i < PUGH_NSTAGGER; i++)
+    {
+      for(j = 0; j < dim; j++)
+      {
+	for(k = 0; k < 2; k++)
+        {
+	  this->ownership[i][k][j] = this->ownership[PUGH_VERTEXCTR][k][j];
+        }
+      }
+    }
+
+    /* correct ownership for staggered grids */
+    for(i = 0; i < PUGH_NSTAGGER; i++)
+    {
+      for (idir = 0; idir < dim; idir++)
+      {	
+	if(this->ub[this_proc][idir] != this->nsize[idir]-1)
+        {
+	  this->ownership[i][1][idir] --;
+        }
+	else
+        {
+	  if(i == PUGH_FACECTR(idir))
+          {
+	    this->ownership[i][1][idir] --;
+          }
+        }
+      }
+    }
+
+    /* set overlaps for staggered grids */
+    for(staglcv = 0; staglcv < PUGH_NSTAGGER; staglcv++)
+    {
+      for(dir = 0; dir < 2*dim; dir++)
+      {
+        for (idir = 0; idir < dim; idir++)
+        {
+          istart = 0;  
+          iend   = this->lnsize[idir];
+          if(staglcv == PUGH_FACECTR(idir))
+          {
+            iend --;
+          }
+          if(dir == 2*idir)
+          {
+            istart =     this->nghostzones[idir];
+            iend   = 2 * this->nghostzones[idir];
+            if(staglcv != PUGH_FACECTR(idir))
+            {
+              istart ++;
+              iend ++;
+            }
+          }
+          if(dir == 2*idir+1)
+          {
+            istart = this->lnsize[idir] - 2 * this->nghostzones[idir] -1;
+            iend   = istart                 + this->nghostzones[idir];
+          }
+          
+          this->overlap[staglcv][0][dir][idir] = istart;
+          this->overlap[staglcv][1][dir][idir] = iend;
+        }
+      }
+    }
+       
+    /* set ghosts for staggered grids */
+    for(staglcv = 0; staglcv < PUGH_NSTAGGER; staglcv++)
+    {
+      for(dir = 0; dir < 2*dim; dir++)
+      {
+	for(idir=0;idir<dim;idir++)
+	{
+	  istart = 0;  
+	  iend   = this->lnsize[idir];
+	  if(staglcv == PUGH_FACECTR(idir))
+          {
+            iend --;
+          }
+	  if(dir == 2*idir)
+	  {
+	    istart = 0;
+	    iend   = this->nghostzones[idir];
+	  }
+	  if(dir == 2*idir+1)
+	  {
+	    istart = this->lnsize[idir] - this->nghostzones[idir];
+	    iend   = this->lnsize[idir];
+	    if(staglcv == PUGH_FACECTR(idir))
+	    {
+	      istart --;
+	      iend --;
+	    }
+	  }
+	  this->ghosts[staglcv][0][dir][idir] = istart;
+	  this->ghosts[staglcv][1][dir][idir] = iend;
+	}
+      }
+    }
+  }
+
+  return 0;
+}
+
+ /*@@
+   @routine    pugh_SetupBoundingBox
+   @date       Mon Nov  8 09:03:40 1999
+   @author     Tom Goodale
+   @desc 
+   Sets up the bounding box info for a pgExtras structure.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+int pugh_SetupBoundingBox(int dim, 
+                          int periodic,
+                          int staggertype,
+                          int *sh,
+                          int *nghosts,
+                          int total_procs,
+                          int *nprocs,
+                          pGExtras *this)
+{
+  int pnum,dir;
+
+  int *step;
+  int *pos;
+
+  step = (int *)malloc(dim*sizeof(int));
+  pos  = (int *)malloc(dim*sizeof(int));
+
+  if(step && pos)
+  {
+    /* Work out the step in each direction */
+    for (dir = 0 ; dir < dim; dir++)
+    {
+      step[dir] = (this->nsize[dir]-1) / nprocs[dir];
+    }
+
+    for(pnum = 0; pnum < total_procs; pnum++)
+    {
+      pugh_DecomposeIJK(dim, pnum,nprocs, pos);
+      
+      for(dir = 0 ; dir < dim; dir++)
+      {
+        if (pos[dir] == 0)
+	{
+	  this->lb[pnum][dir] = 0;
+	}
+	else            
+	{
+	  this->lb[pnum][dir] = pos[dir]*step[dir] +1 - nghosts[dir];
+	  if(staggertype == PUGH_STAGGER)
+          {
+            this->lb[pnum][dir] --;
+          }
+	}
+
+        if (pos[dir] == nprocs[dir]-1) 
+	{
+	  this->ub[pnum][dir] = this->nsize[dir]-1;
+	}
+	else
+	{
+	  this->ub[pnum][dir] = (pos[dir]+1)*step[dir] + this->nghostzones[dir];
+	}
+      }
+    }
+
+  }
+
+  free(step);
+  free(pos);
+
+  return 0;
+}
+
+
+ /*@@
+   @routine    pugh_SetupRemoteSizes
+   @date       Mon Nov  8 09:07:27 1999
+   @author     Tom Goodale
+   @desc 
+   Determines info about the sizes on each processor.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+int pugh_SetupRemoteSizes(int dim, 
+                          int periodic,
+                          int staggertype,
+                          int *sh,
+                          int *nghosts,
+                          int total_procs,
+                          int *nprocs,
+                          pGExtras *this)
+{
+  int pnum;
+  int dir;
+
+  /* Save the number of points on each processor */
+
+  for(pnum = 0; pnum < total_procs; pnum++)
+  {
+
+    this->rnpoints[pnum] = 1;
+    for (dir=0;dir<dim;dir++)
+    {
+      this->rnsize[pnum][dir] = (this->ub[pnum][dir]-this->lb[pnum][dir]+1);
+      this->rnpoints[pnum]  *=  this->rnsize[pnum][dir];
+    }
+  }
+
+  return 0;
+}