Many changes that accumulated while Cactus and Carpet diverged.

Many changes that accumulated while Cactus and Carpet diverged.
Add processor splitting mechanism "along-dir" that splits along a
specified direction.

Rename group PostRestrict to bin POSTRESTRICT.

Prolongate initial data only when desired.  This saves much time.
Sorry, Ian.

Fix bug in time level cycling of grid arrays.  (Note: grid arrays
should not have time levels.)

Fix time_t bug on IRIX.

Make sure that there is no integer overflow when there are many
refinement levels.  Always put parentheses around
(maxreflevelfact/reflevelfact).

Fix typo in Carpet verbose output.

Add debug output in processor splitting.

Communicate in three stages: Irecv, (work), Isend, Wait.  This might
be more efficient.  Much more, potentially.

Fix bug in processor layout of grid arrays.  Sorry, Ian.

Make the interpolator interpolate between time levels.  Untested.

Fix bug in processor communication in interpolator.  Sorry, Ian.

Rewrite prolongation operators to make them twice as fast.  There you
are, Ian.

Move prolongation operator kind handling from data to gdata.

Add official hyperslabbing interfaces to CarpetSlab.

Adapt to new cGH * handling.

darcs-hash:20031105151837-07bb3-758a87ff0355dba053269df4b7d7d79bea018669.gz

3 files changed, 696 insertions, 32 deletions

diff --git a/Carpet/CarpetSlab/interface.ccl b/Carpet/CarpetSlab/interface.ccl
index 96d182492..29a7298be 100644
--- a/Carpet/CarpetSlab/interface.ccl
+++ b/Carpet/CarpetSlab/interface.ccl
@@ -1,5 +1,5 @@
 # Interface definition for thorn CarpetSlab
-# $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/interface.ccl,v 1.6 2003/09/04 16:23:22 tradke Exp $
+# $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/interface.ccl,v 1.7 2003/11/05 16:18:39 schnetter Exp $
 
 implements: Hyperslab
 
@@ -16,3 +16,82 @@ uses include header: gdata.hh
 uses include header: dh.hh
 uses include header: ggf.hh
 uses include header: gh.hh
+
+
+
+CCTK_INT FUNCTION \
+  Hyperslab_Get (CCTK_POINTER_TO_CONST IN cctkGH, \
+                 CCTK_INT IN mapping_handle, \
+                 CCTK_INT IN proc, \
+                 CCTK_INT IN vindex, \
+                 CCTK_INT IN timelevel, \
+                 CCTK_INT IN hdatatype, \
+                 CCTK_POINTER IN hdata)
+
+CCTK_INT FUNCTION \
+  Hyperslab_GetList (CCTK_POINTER_TO_CONST IN cctkGH, \
+                     CCTK_INT IN mapping_handle, \
+                     CCTK_INT IN num_arrays, \
+                     CCTK_INT ARRAY IN procs, \
+                     CCTK_INT ARRAY IN vindices, \
+                     CCTK_INT ARRAY IN timelevels, \
+                     CCTK_INT ARRAY IN hdatatypes, \
+                     CCTK_POINTER ARRAY IN hdata, \
+                     CCTK_INT ARRAY OUT retvals)
+
+#CCTK_INT FUNCTION \
+#  Hyperslab_DefineLocalMappingByIndex (CCTK_POINTER_TO_CONST IN cctkGH, \
+#                                       CCTK_INT IN vindex, \
+#                                       CCTK_INT IN hdim, \
+#                                       CCTK_INT ARRAY IN direction, \
+#                                       CCTK_INT ARRAY IN origin, \
+#                                       CCTK_INT ARRAY IN extent, \
+#                                       CCTK_INT ARRAY IN downsample, \
+#                                       CCTK_INT IN table_handle, \
+#                                       CCTK_INT CCTK_FPOINTER IN \
+#                                         conversion_fn (CCTK_INT IN nelems, \
+#                                                        CCTK_INT IN src_stride, \
+#                                                        CCTK_INT IN dst_stride, \
+#                                                        CCTK_INT IN src_type, \
+#                                                        CCTK_INT IN dst_type, \
+#                                                        CCTK_POINTER_TO_CONST IN from, \
+#                                                        CCTK_POINTER IN to), \
+#                                       CCTK_INT ARRAY OUT hsize_local, \
+#                                       CCTK_INT ARRAY OUT hsize_global, \
+#                                       CCTK_INT ARRAY OUT hoffset_global)
+
+CCTK_INT FUNCTION \
+  Hyperslab_DefineGlobalMappingByIndex (CCTK_POINTER_TO_CONST IN cctkGH, \
+                                        CCTK_INT IN vindex, \
+                                        CCTK_INT IN hdim, \
+                                        CCTK_INT ARRAY IN direction, \
+                                        CCTK_INT ARRAY IN origin, \
+                                        CCTK_INT ARRAY IN extent, \
+                                        CCTK_INT ARRAY IN downsample, \
+                                        CCTK_INT IN table_handle, \
+                                        CCTK_INT CCTK_FPOINTER IN \
+                                          conversion_fn (CCTK_INT IN nelems, \
+                                                         CCTK_INT IN src_stride, \
+                                                         CCTK_INT IN dst_stride, \
+                                                         CCTK_INT IN src_type, \
+                                                         CCTK_INT IN dst_type, \
+                                                         CCTK_POINTER_TO_CONST IN from, \
+                                                         CCTK_POINTER IN to), \
+                                        CCTK_INT ARRAY OUT hsize)
+
+CCTK_INT FUNCTION Hyperslab_FreeMapping (CCTK_INT IN mapping_handle)
+
+
+
+PROVIDES FUNCTION Hyperslab_Get \
+  WITH CarpetSlab_Get LANGUAGE C
+PROVIDES FUNCTION Hyperslab_GetList \
+  WITH CarpetSlab_GetList LANGUAGE C
+PROVIDES FUNCTION Hyperslab_Get \
+  WITH CarpetSlab_Fill LANGUAGE C
+PROVIDES FUNCTION Hyperslab_DefineGlobalMappingByIndex \
+  WITH CarpetSlab_DefineGlobalMappingByIndex LANGUAGE C
+#PROVIDES FUNCTION Hyperslab_DefineLocalMappingByIndex \
+#  WITH CarpetSlab_DefineLocalMappingByIndex LANGUAGE C
+PROVIDES FUNCTION Hyperslab_FreeMapping \
+  WITH CarpetSlab_FreeMapping LANGUAGE C
diff --git a/Carpet/CarpetSlab/src/slab.cc b/Carpet/CarpetSlab/src/slab.cc
index f6d613441..d72499b79 100644
--- a/Carpet/CarpetSlab/src/slab.cc
+++ b/Carpet/CarpetSlab/src/slab.cc
@@ -1,4 +1,4 @@
-// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/src/slab.cc,v 1.12 2003/10/14 16:39:17 schnetter Exp $
+// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/src/slab.cc,v 1.13 2003/11/05 16:18:39 schnetter Exp $
 
 #include <assert.h>
 #include <stdlib.h>
@@ -8,6 +8,8 @@
 
 #include "cctk.h"
 
+#include "util_Table.h"
+
 #include "bbox.hh"
 #include "bboxset.hh"
 #include "dh.hh"
@@ -21,7 +23,7 @@
 #include "slab.hh"
 
 extern "C" {
-  static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/src/slab.cc,v 1.12 2003/10/14 16:39:17 schnetter Exp $";
+  static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/src/slab.cc,v 1.13 2003/11/05 16:18:39 schnetter Exp $";
   CCTK_FILEVERSION(Carpet_CarpetSlab_slab_cc);
 }
 
@@ -33,7 +35,249 @@ namespace CarpetSlab {
   
   
   
-  void* GetSlab (cGH* const cgh,
+  // Mapping object
+  // (just store the mapping)
+  struct mapping {
+    int vindex;
+    int hdim;
+    int * origin;               // [vdim]
+    int * dirs;                 // [hdim]
+    int * stride;               // [hdim]
+    int * length;               // [hdim]
+  };
+  
+  
+  
+  int StoreMapping (mapping * const mp)
+  {
+    int const table = Util_TableCreate (UTIL_TABLE_FLAGS_DEFAULT);
+    assert (table>=0);
+    int const ierr = Util_TableSetPointer (table, mp, "mapping");
+    assert (ierr>=0);
+    return table;
+  }
+  
+  mapping * RetrieveMapping (int const table)
+  {
+    CCTK_POINTER mp;
+    int const ierr = Util_TableGetPointer (table, &mp, "mapping");
+    assert (ierr>=0);
+    return (mapping *)mp;
+  }
+  
+  void DeleteMapping (int const table)
+  {
+    int const ierr = Util_TableDestroy (table);
+    assert (ierr>=0);
+  }
+  
+  
+  
+  void FillSlab (const cGH* const cgh,
+		 const int dest_proc,
+		 const int n,
+		 const int ti,
+		 const int hdim,
+		 const int origin[/*vdim*/],
+		 const int dirs[/*hdim*/],
+		 const int stride[/*hdim*/],
+		 const int length[/*hdim*/],
+                 void* const hdata)
+  {
+    // Check Cactus grid hierarchy
+    assert (cgh);
+    
+    // Check destination processor
+    assert (dest_proc>=-1 && dest_proc<CCTK_nProcs(cgh));
+    
+    // Check variable index
+    assert (n>=0 && n<CCTK_NumVars());
+    
+    // Get info about variable
+    const int group = CCTK_GroupIndexFromVarI(n);
+    assert (group>=0);
+    const int n0 = CCTK_FirstVarIndexI(group);
+    assert (n0>=0);
+    const int var = n - n0;
+    assert (var>=0);
+    
+    // Get info about group
+    cGroup gp;
+    CCTK_GroupData (group, &gp);
+    assert (gp.dim<=dim);
+    assert (CCTK_QueryGroupStorageI(cgh, group));
+    const int typesize = CCTK_VarTypeSize(gp.vartype);
+    assert (typesize>0);
+    
+    if (gp.grouptype==CCTK_GF && reflevel==-1) {
+      CCTK_WARN (0, "It is not possible to use hyperslabbing for a grid function in global mode");
+    }
+    const int rl = gp.grouptype==CCTK_GF ? reflevel : 0;
+    
+    // Check dimension
+    assert (hdim>=0 && hdim<=gp.dim);
+    
+    // Check timelevel
+    const int num_tl = gp.numtimelevels;
+    assert (ti>=0 && ti<num_tl);
+    const int tl = -ti;
+    
+    // Check origin
+//     for (int d=0; d<dim; ++d) {
+//       assert (origin[d]>=0 && origin[d]<=sizes[d]);
+//     }
+    
+    // Check directions
+    for (int dd=0; dd<hdim; ++dd) {
+      assert (dirs[dd]>=1 && dirs[dd]<=dim);
+    }
+    
+    // Check stride
+    for (int dd=0; dd<hdim; ++dd) {
+      assert (stride[dd]>0);
+    }
+    
+    // Check length
+    for (int dd=0; dd<hdim; ++dd) {
+      assert (length[dd]>=0);
+    }
+    
+    // Check extent
+//     for (int dd=0; dd<hdim; ++dd) {
+//       assert (origin[dirs[dd]-1] + length[dd] <= sizes[dirs[dd]]);
+//     }
+    
+    // Get insider information about variable
+    const gh<dim>* myhh;
+    const dh<dim>* mydd;
+    const ggf<dim>* myff;
+    assert (group < (int)arrdata.size());
+    myhh = arrdata[group].hh;
+    assert (myhh);
+    mydd = arrdata[group].dd;
+    assert (mydd);
+    assert (var < (int)arrdata[group].data.size());
+    myff = arrdata[group].data[var];
+    assert (myff);
+    
+    // Detemine collecting processor
+    const int collect_proc = dest_proc<0 ? 0 : dest_proc;
+    
+    // Determine own rank
+    const int rank = CCTK_MyProc(cgh);
+    
+    // Calculate global size
+    int totalsize = 1;
+    for (int dd=0; dd<hdim; ++dd) {
+      totalsize *= length[dd];
+    }
+    
+    // Allocate memory
+    assert (hdata);
+    if (dest_proc==-1 || rank==dest_proc) {
+      memset (hdata, 0, totalsize * typesize);
+    }
+    
+    // Get sample data
+    const gdata<dim>* mydata;
+    mydata = (*myff)(tl, rl, 0, 0);
+    
+    // Stride of data in memory
+    const vect<int,dim> str = mydata->extent().stride();
+    
+    // Stride of collected data
+    vect<int,dim> hstr = str;
+    for (int dd=0; dd<hdim; ++dd) {
+      hstr[dirs[dd]-1] *= stride[dd];
+    }
+    
+    // Lower bound of collected data
+    vect<int,dim> hlb(0);
+    for (int d=0; d<gp.dim; ++d) {
+      hlb[d] = origin[d] * str[d];
+    }
+    
+    // Upper bound of collected data
+    vect<int,dim> hub = hlb;
+    for (int dd=0; dd<hdim; ++dd) {
+      hub[dirs[dd]-1] += (length[dd]-1) * hstr[dirs[dd]-1];
+    }
+    
+    // Calculate extent to collect
+    const bbox<int,dim> hextent (hlb, hub, hstr);
+    assert (hextent.size() == totalsize);
+    
+    // Create collector data object
+    void* myhdata = rank==collect_proc ? hdata : 0;
+    gdata<dim>* const alldata = mydata->make_typed(-1);
+    alldata->allocate (hextent, collect_proc, myhdata);
+    
+    // Done with the temporary stuff
+    mydata = 0;
+    
+    for (comm_state<dim> state; !state.done(); state.step()) {
+      
+      // Loop over all components, copying data from them
+      BEGIN_LOCAL_COMPONENT_LOOP (cgh, gp.grouptype) {
+        
+        // Get data object
+        mydata = (*myff)(tl, rl, component, mglevel);
+        
+        // Calculate overlapping extents
+        const bboxset<int,dim> myextents
+          = ((mydd->boxes[rl][component][mglevel].sync_not
+              | mydd->boxes[rl][component][mglevel].interior)
+             & hextent);
+        
+        // Loop over overlapping extents
+        for (bboxset<int,dim>::const_iterator ext_iter = myextents.begin();
+             ext_iter != myextents.end();
+             ++ext_iter) {
+          
+          // Copy data
+          alldata->copy_from (state, mydata, *ext_iter);
+          
+        }
+        
+      } END_LOCAL_COMPONENT_LOOP;
+      
+    } // for step
+    
+    // Copy result to all processors
+    if (dest_proc == -1) {
+      vector<gdata<dim>*> tmpdata(CCTK_nProcs(cgh));
+      vector<comm_state<dim> > state;
+      
+      for (int proc=0; proc<CCTK_nProcs(cgh); ++proc) {
+        if (proc != collect_proc) {
+          void* myhdata = rank==proc ? hdata : 0;
+          tmpdata[proc] = mydata->make_typed(-1);
+          tmpdata[proc]->allocate (alldata->extent(), proc, myhdata);
+          tmpdata[proc]->copy_from (state[proc], alldata, alldata->extent());
+        }
+      }
+      
+      for (int proc=0; proc<CCTK_nProcs(cgh); ++proc) {
+        if (proc != collect_proc) {
+          tmpdata[proc]->copy_from (state[proc], alldata, alldata->extent());
+        }
+      }
+      
+      for (int proc=0; proc<CCTK_nProcs(cgh); ++proc) {
+        if (proc != collect_proc) {
+          tmpdata[proc]->copy_from (state[proc], alldata, alldata->extent());
+          delete tmpdata[proc];
+        }
+      }
+      
+    } // Copy result
+    
+    delete alldata;
+  }
+  
+  
+  
+  void* GetSlab (const cGH* const cgh,
 		 const int dest_proc,
 		 const int n,
 		 const int ti,
@@ -134,6 +378,7 @@ namespace CarpetSlab {
     // Allocate memory
     void* hdata = 0;
     if (dest_proc==-1 || rank==dest_proc) {
+      assert (0);
       hdata = malloc(totalsize * typesize);
       assert (hdata);
       memset (hdata, 0, totalsize * typesize);
@@ -176,45 +421,61 @@ namespace CarpetSlab {
     // Done with the temporary stuff
     mydata = 0;
     
-    // Loop over all components, copying data from them
-    BEGIN_LOCAL_COMPONENT_LOOP (cgh, gp.grouptype) {
-      
-      const int c = gp.grouptype==CCTK_GF ? component : 0;
+    for (comm_state<dim> state; !state.done(); state.step()) {
       
-      // Get data object
-      mydata = (*myff)(tl, rl, c, mglevel);
-      
-      // Calculate overlapping extents
-      const bboxset<int,dim> myextents
-        = ((mydd->boxes[rl][c][mglevel].sync_not
-            | mydd->boxes[rl][c][mglevel].interior)
-           & hextent);
-      
-      // Loop over overlapping extents
-      for (bboxset<int,dim>::const_iterator ext_iter = myextents.begin();
-           ext_iter != myextents.end();
-           ++ext_iter) {
+      // Loop over all components, copying data from them
+      BEGIN_LOCAL_COMPONENT_LOOP (cgh, gp.grouptype) {
         
-        // Copy data
-        alldata->copy_from (mydata, *ext_iter);
+        // Get data object
+        mydata = (*myff)(tl, rl, component, mglevel);
         
-      }
+        // Calculate overlapping extents
+        const bboxset<int,dim> myextents
+          = ((mydd->boxes[rl][component][mglevel].sync_not
+              | mydd->boxes[rl][component][mglevel].interior)
+             & hextent);
+        
+        // Loop over overlapping extents
+        for (bboxset<int,dim>::const_iterator ext_iter = myextents.begin();
+             ext_iter != myextents.end();
+             ++ext_iter) {
+          
+          // Copy data
+          alldata->copy_from (state, mydata, *ext_iter);
+          
+        }
+        
+      } END_LOCAL_COMPONENT_LOOP;
       
-    } END_LOCAL_COMPONENT_LOOP;
+    } // for step
     
     // Copy result to all processors
     if (dest_proc == -1) {
+      vector<gdata<dim>*> tmpdata(CCTK_nProcs(cgh));
+      vector<comm_state<dim> > state;
+      
       for (int proc=0; proc<CCTK_nProcs(cgh); ++proc) {
         if (proc != collect_proc) {
-          
           void* myhdata = rank==proc ? hdata : 0;
-          gdata<dim>* const tmpdata = mydata->make_typed(-1);
-          tmpdata->allocate (alldata->extent(), proc, myhdata);
-          tmpdata->copy_from (alldata, alldata->extent());
-          delete tmpdata;
-          
+          tmpdata[proc] = mydata->make_typed(-1);
+          tmpdata[proc]->allocate (alldata->extent(), proc, myhdata);
+          tmpdata[proc]->copy_from (state[proc], alldata, alldata->extent());
         }
       }
+      
+      for (int proc=0; proc<CCTK_nProcs(cgh); ++proc) {
+        if (proc != collect_proc) {
+          tmpdata[proc]->copy_from (state[proc], alldata, alldata->extent());
+        }
+      }
+      
+      for (int proc=0; proc<CCTK_nProcs(cgh); ++proc) {
+        if (proc != collect_proc) {
+          tmpdata[proc]->copy_from (state[proc], alldata, alldata->extent());
+          delete tmpdata[proc];
+        }
+      }
+      
     } // Copy result
     
     delete alldata;
@@ -225,6 +486,271 @@ namespace CarpetSlab {
   
   
   
+  CCTK_INT CarpetSlab_Get (cGH const * const cctkGH,
+                           CCTK_INT const mapping_handle,
+                           CCTK_INT const proc,
+                           CCTK_INT const vindex,
+                           CCTK_INT const timelevel,
+                           CCTK_INT const hdatatype,
+                           void * const hdata)
+  {
+    // Check arguments
+    assert (cctkGH);
+    assert (mapping_handle>=0);
+    assert (proc==-1 || proc>=0 && proc<CCTK_nProcs(cctkGH));
+    assert (vindex>=0 && vindex<CCTK_NumVars());
+    assert (timelevel>=0);
+    assert (hdatatype>=0);
+    assert (hdata);
+    
+    // Get mapping
+    const mapping * const mp = RetrieveMapping (mapping_handle);
+    assert (mp);
+    
+    // Calculate total size
+    size_t size = 1;
+    for (size_t d=0; d<(size_t)mp->hdim; ++d) {
+      size *= mp->length[d];
+    }
+    
+    // Get type size
+    size_t const sz = CCTK_VarTypeSize (hdatatype);
+    assert (sz>0);
+    
+    // Forward call
+    FillSlab (cctkGH, proc, vindex, timelevel,
+              mp->hdim, mp->origin, mp->dirs, mp->stride, mp->length, hdata);
+    
+    return 0;
+  }
+  
+  
+  
+  CCTK_INT CarpetSlab_GetList (cGH const * const cctkGH,
+                               CCTK_INT const mapping_handle,
+                               CCTK_INT const num_arrays,
+                               CCTK_INT const * const procs,
+                               CCTK_INT const * const vindices,
+                               CCTK_INT const * const timelevels,
+                               CCTK_INT const * const hdatatypes,
+                               void * const * const hdata,
+                               CCTK_INT * const retvals)
+  {
+    // Check arguments
+    assert (cctkGH);
+    assert (mapping_handle>=0);
+    assert (num_arrays>=0);
+    assert (procs);
+    assert (vindices);
+    assert (timelevels);
+    assert (hdatatypes);
+    assert (hdata);
+    assert (retvals);
+    
+    // Remember whether there were errors
+    bool everyting_okay = true;
+    
+    // Loop over all slabs
+    for (int n=0; n<num_arrays; ++n) {
+      // Forward call
+      retvals[n] = CarpetSlab_Get (cctkGH, mapping_handle, procs[n],
+                                   vindices[n], timelevels[n], hdatatypes[n],
+                                   hdata[n]);
+      everyting_okay = everyting_okay && retvals[n];
+    }
+    
+    return everyting_okay ? 0 : -1;
+  }
+  
+  
+  
+  CCTK_INT CarpetSlab_DefineLocalMappingByIndex (cGH const * const cctkGH,
+                                                 CCTK_INT const vindex,
+                                                 CCTK_INT const hdim,
+                                                 CCTK_INT const * const direction,
+                                                 CCTK_INT const * const origin,
+                                                 CCTK_INT const * const extent,
+                                                 CCTK_INT const * const downsample,
+                                                 CCTK_INT const table_handle,
+                                                 CCTK_INT (* const conversion_fn) (CCTK_INT const nelems,
+                                                                                   CCTK_INT const src_stride,
+                                                                                   CCTK_INT const dst_stride,
+                                                                                   CCTK_INT const src_type,
+                                                                                   CCTK_INT const dst_type,
+                                                                                   void const * const from,
+                                                                                   void * const to),
+                                                 CCTK_INT * const hsize_local,
+                                                 CCTK_INT * const hsize_global,
+                                                 CCTK_INT * const hoffset_global)
+  {
+    CCTK_WARN (0, "not implemented");
+    return 0;
+  }
+  
+  
+  
+  CCTK_INT CarpetSlab_DefineGlobalMappingByIndex (cGH const * const cctkGH,
+                                                  CCTK_INT const vindex,
+                                                  CCTK_INT const hdim,
+                                                  CCTK_INT const * const direction,
+                                                  CCTK_INT const * const origin,
+                                                  CCTK_INT const * const extent,
+                                                  CCTK_INT const * const downsample,
+                                                  CCTK_INT const table_handle,
+                                                  CCTK_INT (* const conversion_fn) (CCTK_INT const nelems,
+                                                                                    CCTK_INT const src_stride,
+                                                                                    CCTK_INT const dst_stride,
+                                                                                    CCTK_INT const src_type,
+                                                                                    CCTK_INT const dst_type,
+                                                                                    void const * const from,
+                                                                                    void * const to),
+                                                  CCTK_INT * const hsize)
+  {
+    // Check arguments
+    assert (cctkGH);
+    assert (vindex>=0 && vindex<CCTK_NumVars());
+    assert (hdim>=0 && hdim<=dim);
+    assert (direction);
+    assert (origin);
+    assert (extent);
+    assert (downsample);
+    assert (table_handle>=0);
+    assert (hsize);
+    
+    // Get more information
+    int const vdim = CCTK_GroupDimFromVarI (vindex);
+    assert (vdim>=0 && vdim<dim);
+    assert (hdim<=vdim);
+    
+    // Not implemented
+    assert (! conversion_fn);
+    
+    // Allocate memory
+    mapping * mp = new mapping;
+    mp->origin = new int[vdim];
+    mp->dirs   = new int[hdim];
+    mp->stride = new int[hdim];
+    mp->length = new int[hdim];
+    
+    // Calculate more convenient representation of the direction
+    int dirs[dim];              // should really be dirs[hdim]
+    // The following if statement is written according to the
+    // definition of "dir".
+    if (hdim==1) {
+      // 1-dimensional hyperslab
+      int mydir = 0;
+      for (int d=0; d<vdim; ++d) {
+	if (direction[d]!=0) {
+	  mydir = d+1;
+	  break;
+	}
+      }
+      assert (mydir>0);
+      for (int d=0; d<vdim; ++d) {
+	if (d == mydir-1) {
+	  assert (direction[d]!=0);
+	} else {
+	  assert (direction[d]==0);
+	}
+      }
+      dirs[0] = mydir;
+    } else if (hdim==vdim) {
+      // vdim-dimensional hyperslab
+      for (int d=0; d<vdim; ++d) {
+	dirs[d] = d+1;
+      }
+    } else if (hdim==2) {
+      // 2-dimensional hyperslab with vdim==3
+      assert (vdim==3);
+      int mydir = 0;
+      for (int d=0; d<vdim; ++d) {
+	if (direction[d]==0) {
+	  mydir = d+1;
+	  break;
+	}
+      }
+      assert (mydir>0);
+      for (int d=0; d<vdim; ++d) {
+	if (d == mydir-1) {
+	  assert (direction[d]==0);
+	} else {
+	  assert (direction[d]!=0);
+	}
+      }
+      int dd=0;
+      for (int d=0; d<vdim; ++d) {
+	if (d != mydir-1) {
+	  dirs[dd] = d+1;
+	  ++dd;
+	}
+      }
+      assert (dd==hdim);
+    } else {
+      assert (0);
+    }
+    // Fill remaining length
+    for (int d=vdim; d<dim; ++d) {
+      dirs[d] = d+1;
+    }
+    
+    // Calculate lengths
+    for (int dd=0; dd<hdim; ++dd) {
+      if (extent[dd]<0) {
+	int gsh[dim];
+	int ierr = CCTK_GroupgshVI(cctkGH, dim, gsh, vindex);
+	assert (!ierr);
+	const int totlen = gsh[dirs[dd]-1];
+	assert (totlen>=0);
+	// Partial argument check
+	assert (origin[dirs[dd]-1]>=0);
+	assert (origin[dirs[dd]-1]<=totlen);
+	assert (downsample[dd]>0);
+	hsize[dd] = (totlen - origin[dirs[dd]-1]) / downsample[dd];
+      } else {
+	hsize[dd] = extent[dd];
+      }
+      assert (hsize[dd]>=0);
+    }
+    
+    // Store information
+    mp->vindex = vindex;
+    mp->hdim = hdim;
+    for (size_t d=0; d<(size_t)hdim; ++d) {
+      mp->origin[d] = origin[d];
+      mp->dirs[d]   = dirs[d];
+      mp->stride[d] = downsample[d];
+      mp->length[d] = hsize[d];
+    }
+    
+    return 0;
+  }
+  
+  
+  
+  CCTK_INT CarpetSlab_FreeMapping (CCTK_INT const mapping_handle)
+  {
+    // Check arguments
+    assert (mapping_handle>=0);
+    
+    // Get mapping
+    mapping * mp = RetrieveMapping (mapping_handle);
+    assert (mp);
+    
+    // Delete storage
+    DeleteMapping (mapping_handle);
+    
+    // Free memory
+    delete [] mp->origin;
+    delete [] mp->dirs;
+    delete [] mp->stride;
+    delete [] mp->length;
+    delete mp;
+    
+    return 0;
+  }
+  
+  
+  
   int Hyperslab_GetHyperslab (cGH* const GH,
 			      const int target_proc,
 			      const int vindex,
diff --git a/Carpet/CarpetSlab/src/slab.h b/Carpet/CarpetSlab/src/slab.h
index c72222700..e57205aa0 100644
--- a/Carpet/CarpetSlab/src/slab.h
+++ b/Carpet/CarpetSlab/src/slab.h
@@ -1,4 +1,4 @@
-/* $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/src/slab.h,v 1.1 2002/10/24 10:53:48 schnetter Exp $ */
+/* $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetSlab/src/slab.h,v 1.2 2003/11/05 16:18:39 schnetter Exp $ */
 
 #ifndef CARPETSLAB_H
 #define CARPETSLAB_H
@@ -10,6 +10,65 @@ namespace CarpetSlab {
   extern "C" {
 #endif
     
+    CCTK_INT CarpetSlab_Get (cGH const * const cctkGH,
+                             CCTK_INT const mapping_handle,
+                             CCTK_INT const proc,
+                             CCTK_INT const vindex,
+                             CCTK_INT const timelevel,
+                             CCTK_INT const hdatatype,
+                             void * const hdata);
+    
+    CCTK_INT CarpetSlab_GetList (cGH const * const cctkGH,
+                                 CCTK_INT const mapping_handle,
+                                 CCTK_INT const num_arrays,
+                                 CCTK_INT const * const procs,
+                                 CCTK_INT const * const vindices,
+                                 CCTK_INT const * const timelevels,
+                                 CCTK_INT const * const hdatatypes,
+                                 void * const * const hdata,
+                                 CCTK_INT * const retvals);
+    
+    CCTK_INT CarpetSlab_DefineLocalMappingByIndex (cGH const * const cctkGH,
+                                                   CCTK_INT const vindex,
+                                                   CCTK_INT const hdim,
+                                                   CCTK_INT const * const direction,
+                                                   CCTK_INT const * const origin,
+                                                   CCTK_INT const * const extent,
+                                                   CCTK_INT const * const downsample,
+                                                   CCTK_INT const table_handle,
+                                                   CCTK_INT (* const conversion_fn) (CCTK_INT const nelems,
+                                                                                     CCTK_INT const src_stride,
+                                                                                     CCTK_INT const dst_stride,
+                                                                                     CCTK_INT const src_type,
+                                                                                     CCTK_INT const dst_type,
+                                                                                     void const * const from,
+                                                                                     void * const to),
+                                                   CCTK_INT * const hsize_local,
+                                                   CCTK_INT * const hsize_global,
+                                                   CCTK_INT * const hoffset_global);
+
+    CCTK_INT CarpetSlab_DefineGlobalMappingByIndex (cGH const * const cctkGH,
+                                                    CCTK_INT const vindex,
+                                                    CCTK_INT const hdim,
+                                                    CCTK_INT const * const direction,
+                                                    CCTK_INT const * const origin,
+                                                    CCTK_INT const * const extent,
+                                                    CCTK_INT const * const downsample,
+                                                    CCTK_INT const table_handle,
+                                                    CCTK_INT (* const conversion_fn) (CCTK_INT const nelems,
+                                                                                      CCTK_INT const src_stride,
+                                                                                      CCTK_INT const dst_stride,
+                                                                                      CCTK_INT const src_type,
+                                                                                      CCTK_INT const dst_type,
+                                                                                      void const * const from,
+                                                                                      void * const to),
+                                                    CCTK_INT * const hsize);
+    
+    CCTK_INT CarpetSlab_FreeMapping (CCTK_INT const mapping_handle);
+    
+    
+    
+    /* Old interface -- don't use */
     int Hyperslab_GetHyperslab (cGH* const GH,
 				const int target_proc,
 				const int vindex,