[MOAB-dev] r1301 - in MOAB/trunk: . parallel tools/iMesh

tautges at mcs.anl.gov tautges at mcs.anl.gov
Mon Oct 8 14:15:01 CDT 2007


Author: tautges
Date: 2007-10-08 14:15:00 -0500 (Mon, 08 Oct 2007)
New Revision: 1301

Modified:
   MOAB/trunk/DenseTagCollections.hpp
   MOAB/trunk/MBBits.hpp
   MOAB/trunk/MBCore.cpp
   MOAB/trunk/MBCore.hpp
   MOAB/trunk/MBInterface.hpp
   MOAB/trunk/MBMeshSet.cpp
   MOAB/trunk/MBMeshSet.hpp
   MOAB/trunk/MBReadUtil.cpp
   MOAB/trunk/MBSkinner.cpp
   MOAB/trunk/ReadNCDF.cpp
   MOAB/trunk/ReadNCDF.hpp
   MOAB/trunk/SparseTagCollections.cpp
   MOAB/trunk/SparseTagCollections.hpp
   MOAB/trunk/TagServer.cpp
   MOAB/trunk/TagServer.hpp
   MOAB/trunk/mbparallelcomm_test.cpp
   MOAB/trunk/parallel/MBParallelComm.cpp
   MOAB/trunk/parallel/MBParallelComm.hpp
   MOAB/trunk/parallel/ReadParallel.cpp
   MOAB/trunk/tools/iMesh/Makefile.am
Log:
Fixing more bugs/issues with parallel communication, primarily with
loading files using broadcast and delete.

Passes make check.

MBCore:
- first cut at implementing replace_entities.  Not used here yet, but
will be needed when we start migrating entities.  This isn't a
complete implementation, but it's a start.

tools/iMesh/Makefile.am: adding flags to testc_cbind for extra C++
libs.

MBSkinner.cpp: treat structured mesh case (where you can't get
connectivity as a const*)

Sparse, Dense, Bit Tag Collections, TagServer: add functions to get
all entities with a given tag set, rather than by type.

ReadNCDF.cpp: 
- fix inclusion of read entities in a file set; now, it
just gets all entities at start and end and subtracts to get file
entities
- set global ids on elements as a block
- set vertex global ids too

mbparallelcomm_test:
- split into tests with overlapping linear vertex sequences, 3d
structured mesh sharing 2d faces, and file read partitioning by
material set

MBMeshSet, MBInterface, MBCore: adding a replace_entities function

MBCore: replacing #ifdef MPI with #ifdef USE_MPI, since mpi.h may
use MPI

parallel/MBParallelComm.cpp:
- better error reporting
- adding sets to entities returned from broadcast_entities
- changed communication of dense tags, since it is possible to find
out from sequences whether dense tags are really set
- adding a variant of resolve_shared_ents taking a dimension of
entities to look at for shared sub-entities
- changed resolve_shared_ents to handle dimension of shared entities
better
- added get_shared_ents function

parallel/ReadParallel.cpp:
- better error reporting
- fixed communication of file_set on non-reading processors
- added READ_DELETE enumeration and option for parallel reading (not
implemented yet)
- added ability to balance partition-designated sets around processors
in the case where a partition tag isn't used




Modified: MOAB/trunk/DenseTagCollections.hpp
===================================================================
--- MOAB/trunk/DenseTagCollections.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/DenseTagCollections.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -196,6 +196,9 @@
   //! get the entities
   MBErrorCode get_entities(MBEntityType type, MBRange& entities);
   
+  //! get the entities
+  MBErrorCode get_entities(MBRange& entities);
+  
   //! get the entities with a value
   MBErrorCode get_entities_with_tag_value(const MBEntityType type, const void* value, 
                                           MBRange &entities);
@@ -330,6 +333,27 @@
   return MB_SUCCESS;
 }
 
+//! get the entities
+inline MBErrorCode DensePageGroup::get_entities(MBRange& entities)
+{
+  std::vector<DensePage>::iterator iter;
+  int dum =0;
+  for (MBEntityType type = MBENTITYSET; type >= MBVERTEX; type--) {
+    const std::vector<DensePage>::iterator end = mDensePages[type].end();
+    MBEntityHandle handle = CREATE_HANDLE(type, 0, dum);
+    MBEntityID first_time = MB_START_ID; // Don't want zero-ID handle at start of range.
+    MBRange::iterator insert_pos = entities.begin();
+    for(iter = mDensePages[type].begin(); iter != end; ++iter, handle += DensePage::mPageSize)
+    {
+      if (iter->has_data())
+        insert_pos = entities.insert( insert_pos, handle + first_time, handle + DensePage::mPageSize - 1 );
+      first_time = 0;
+    }
+  }
+  
+  return MB_SUCCESS;
+}
+
 //! get number of entities of type
 inline MBErrorCode DensePageGroup::get_number_entities(MBEntityType type, int& entities)
 {
@@ -386,6 +410,9 @@
   MBErrorCode get_entities(const MBTagId tag_id, const MBEntityType type, MBRange& entities);
   
   //! get the entities with a tag
+  MBErrorCode get_entities(const MBTagId tag_id, MBRange& entities);
+  
+  //! get the entities with a tag
   MBErrorCode get_entities(const MBRange &range,
                             const MBTagId tag_id, const MBEntityType type, MBRange& entities);
   
@@ -472,7 +499,17 @@
   return MB_SUCCESS;
 }
 
+// get the entities with a tag
+inline MBErrorCode DenseTagSuperCollection::get_entities(const MBTagId tag_id, 
+                                                         MBRange& entities)
+{
+  std::vector<DensePageGroup*>::iterator group = mDensePageGroups.begin() + tag_id;
+  if (group >= mDensePageGroups.end() || !*group)
+    return MB_TAG_NOT_FOUND;
 
+  return (*group)->get_entities(entities);
+}
+
 #endif
 
 

Modified: MOAB/trunk/MBBits.hpp
===================================================================
--- MOAB/trunk/MBBits.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBBits.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -465,6 +465,8 @@
   //! set the bits associated with an entity handle, only if memory has been allocated
   MBErrorCode weak_set_bits(MBTagId tag_id, MBEntityHandle handle, unsigned char bits);
 
+  MBErrorCode get_entities(MBTagId tag_id, MBRange& entities);
+
   MBErrorCode get_entities(MBTagId tag_id, MBEntityType type, MBRange& entities);
 
   MBErrorCode get_entities_with_tag_value(MBTagId tag_id, 
@@ -539,6 +541,22 @@
   return mBitPageGroups[TYPE_FROM_HANDLE(handle)][tag_id]->weak_set_bits(handle, bits);
 }
 
+inline MBErrorCode MBBitServer::get_entities(MBTagId tag_id, MBRange& entities)
+{
+  --tag_id; // First ID is 1.
+  MBErrorCode result = MB_SUCCESS;
+  if(tag_id >= mBitPageGroupsSize || (*mBitPageGroups)[tag_id] == NULL)
+    result = MB_FAILURE;
+  else {
+    for (MBEntityType type = MBVERTEX; type != MBMAXTYPE; type++) {
+      MBErrorCode tmp_result = mBitPageGroups[type][tag_id]->get_entities(type, entities);
+      if (MB_SUCCESS != tmp_result) result = tmp_result;
+    }
+  }
+  
+  return result;
+}
+
 inline MBErrorCode MBBitServer::get_entities(MBTagId tag_id, MBEntityType type, MBRange& entities)
 {
   --tag_id; // First ID is 1.

Modified: MOAB/trunk/MBCore.cpp
===================================================================
--- MOAB/trunk/MBCore.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBCore.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -41,7 +41,7 @@
 #include "MBWriterIface.hpp"
 #include "MBHandleUtils.hpp"
 
-#ifdef MPI
+#ifdef USE_MPI
 #include "mpi.h"
 #include "MBParallelComm.hpp"
 #include "ReadParallel.hpp"
@@ -361,7 +361,7 @@
   std::string parallel_opt;
   rval = opts.get_option( "PARALLEL", parallel_opt);
   if (MB_SUCCESS == rval && !parallel_opt.empty()) {
-#ifdef MPI    
+#ifdef USE_MPI    
     return ReadParallel(this).load_file(file_name, file_set, opts,
                                         block_id_list, num_blocks);
 #else
@@ -2318,6 +2318,21 @@
     return MB_ENTITY_NOT_FOUND;
 }
 
+// replace entities in a meshset; entities appear in pairs,
+// old then new entity in each pair
+bool MBCore::replace_entities(MBEntityHandle meshset, 
+                              MBEntityHandle *entities,
+                              int num_entities) 
+{
+  if (0 == meshset) return false;
+  
+  MBMeshSet* set = get_mesh_set( sequence_manager(), meshset );
+  if (set)
+    return set->replace_entities( meshset, entities, num_entities, a_entity_factory() );
+  else
+    return false;
+}
+
 MBErrorCode MBCore::get_parent_meshsets(const MBEntityHandle meshset,
                                           std::vector<MBEntityHandle> &parents,
                                           const int num_hops) const

Modified: MOAB/trunk/MBCore.hpp
===================================================================
--- MOAB/trunk/MBCore.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBCore.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -767,6 +767,10 @@
                                        const MBEntityHandle *entities,
                                        const int num_entities);
 
+  virtual bool replace_entities(MBEntityHandle meshset, 
+                                MBEntityHandle *entities,
+                                int num_entities);
+
   //------MeshSet Parent/Child functions------
   
   //! get parent meshsets

Modified: MOAB/trunk/MBInterface.hpp
===================================================================
--- MOAB/trunk/MBInterface.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBInterface.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -1272,6 +1272,17 @@
                                       const MBEntityHandle *entities,
                                       const int num_entities) = 0;
 
+    //! Replace entities in a set with other entities
+    /** Replace entities in a set with other entities; entity list
+     * specified in pairs of (old, new)
+     * \param meshset Mesh set being modified
+     * \param entities Pairs of old/new entities
+     * \param num_entities Number of entities \em{total} in entities list
+     * \return bool If true, one or more entities were replaced
+    */
+  virtual bool replace_entities(MBEntityHandle meshset, 
+                                MBEntityHandle *entities,
+                                int num_entities) = 0;
     //@}
 
     //! \name MeshSet parent/child functions

Modified: MOAB/trunk/MBMeshSet.cpp
===================================================================
--- MOAB/trunk/MBMeshSet.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBMeshSet.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -279,6 +279,35 @@
 }
 
 
+bool MBMeshSet_MBRange::replace_entities(MBEntityHandle set_handle,
+                                         MBEntityHandle *entities,
+                                         int num_entities,
+                                         AEntityFactory* mAdjFact)
+{
+  bool was_contained = false;
+  MBRange::iterator rit;
+  for (int i = 0; i < num_entities; i += 2) {
+    if ((rit = mRange.find(entities[i])) != mRange.end()) {
+      mRange.erase(rit);
+      mRange.insert(entities[i+1]);
+      was_contained = true;
+
+      if (tracking() && mAdjFact)
+        mAdjFact->remove_adjacency(entities[i], set_handle),
+          mAdjFact->add_adjacency(entities[i+1], set_handle);
+    }
+    
+    if (remove_parent(entities[i]))
+      add_parent(entities[i+1]), was_contained = true;
+
+    if (remove_child(entities[i]))
+      add_child(entities[i+1]), was_contained = true;
+  }
+  
+  return was_contained;
+}
+
+
 MBErrorCode MBMeshSet_MBRange::add_entities( const MBEntityHandle *entities,
                                              const int num_entities,
                                              MBEntityHandle mEntityHandle,
@@ -430,6 +459,35 @@
   return MB_SUCCESS;
 }
 
+bool MBMeshSet_Vector::replace_entities(MBEntityHandle set_handle,
+                                        MBEntityHandle *entities,
+                                        int num_entities,
+                                        AEntityFactory* mAdjFact)
+{
+  bool was_contained = false;
+  std::vector<MBEntityHandle>::iterator vit;
+  for (int i = 0; i < num_entities; i += 2) {
+    if ((vit = std::find(mVector.begin(), mVector.end(), entities[i])) 
+        != mVector.end()) {
+      *vit = entities[i+1];
+      was_contained = true;
+
+      if(tracking() && mAdjFact)
+        mAdjFact->remove_adjacency(entities[i], set_handle),
+          mAdjFact->add_adjacency(entities[i+1], set_handle);
+    }
+    
+    if (remove_parent(entities[i]))
+      add_parent(entities[i+1]), was_contained = true;
+
+    if (remove_child(entities[i]))
+      add_child(entities[i+1]), was_contained = true;
+  }
+  
+  return was_contained;
+}
+
+
 void MBMeshSet_Vector::vector_to_range( std::vector<MBEntityHandle>& vect, MBRange& range )
 {
   std::sort( vect.begin(), vect.end() );

Modified: MOAB/trunk/MBMeshSet.hpp
===================================================================
--- MOAB/trunk/MBMeshSet.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBMeshSet.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -129,7 +129,7 @@
     //! remove a child from this meshset; returns true if child was removed, 0 if it was
     //! not a child of this meshset
   int remove_child(MBEntityHandle child);
-  
+
   unsigned flags() const { return mFlags; }
   //! returns whether entities of meshsets know this meshset 
   int tracking()     const { return mFlags & MESHSET_TRACK_OWNER; }
@@ -139,6 +139,12 @@
 
    //!  PURE VIRTUAL FUNCTIONS overwritten by derived classes  *******************
 
+    //! replace one entity with another in the set (contents and parent/child
+    //! lists); returns whether it was replaced or not
+  inline bool replace_entities(MBEntityHandle my_handle, 
+                               MBEntityHandle *entities, int num_entities,
+                               AEntityFactory* mAdjFact);
+  
   inline MBErrorCode clear( MBEntityHandle myhandle, AEntityFactory* adjacencies );
 
   inline MBErrorCode get_entities(std::vector<MBEntityHandle>& entities) const;
@@ -252,8 +258,11 @@
 };
 
 #define MESH_SET_VIRTUAL_FUNCTIONS      \
-  MBErrorCode clear( MBEntityHandle my_handle, AEntityFactory* adjacencies);                                                          \
+  MBErrorCode clear( MBEntityHandle my_handle, AEntityFactory* adjacencies); \
   \
+  bool replace_entities(MBEntityHandle my_handle, MBEntityHandle *entities, \
+                                       int num_entities, AEntityFactory* mAdjFact); \
+  \
   inline MBErrorCode get_entities(std::vector<MBEntityHandle>& entities) const;       \
   \
   inline MBErrorCode get_entities(MBRange& entities) const;                           \
@@ -374,6 +383,18 @@
     std::vector<MBEntityHandle> mVector;
 };
 
+inline bool
+MBMeshSet::replace_entities( MBEntityHandle my_handle, 
+                             MBEntityHandle *entities, int num_entities,
+                             AEntityFactory* mAdjFact)
+{ 
+  return vector_based() ?
+    reinterpret_cast<MBMeshSet_Vector* >(this)->replace_entities(my_handle, entities, num_entities,
+                                                                 mAdjFact) :
+    reinterpret_cast<MBMeshSet_MBRange*>(this)->replace_entities(my_handle, entities, num_entities,
+                                                                 mAdjFact);
+}
+
 inline MBErrorCode 
 MBMeshSet::clear( MBEntityHandle myhandle, AEntityFactory* adjacencies )
 { 

Modified: MOAB/trunk/MBReadUtil.cpp
===================================================================
--- MOAB/trunk/MBReadUtil.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBReadUtil.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -208,7 +208,7 @@
   std::pair<MBRange::const_iterator, MBRange::const_iterator> pair_it =
     partition.equal_range(MBENTITYSET);
 
-  MBErrorCode result;
+  MBErrorCode result = MB_SUCCESS;
   for (MBRange::const_iterator rit = pair_it.first; 
        rit != pair_it.second; rit++) {
     MBErrorCode tmp_result = 

Modified: MOAB/trunk/MBSkinner.cpp
===================================================================
--- MOAB/trunk/MBSkinner.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/MBSkinner.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -256,7 +256,7 @@
 
   MBRange::const_iterator iter, end_iter;
   end_iter = source_entities.end();
-  const MBEntityHandle *conn;
+  const MBEntityHandle *tmp_conn, *conn;
   MBEntityHandle match;
 
   direction direct;
@@ -264,14 +264,26 @@
     // assume we'll never have more than 32 vertices on a facet (checked
     // with assert later)
   static MBEntityHandle sub_conn[32];
+  static std::vector<MBEntityHandle> tmp_conn_vec;
   int num_nodes;
 
   // for each source entity
   for(iter = source_entities.begin(); iter != end_iter; ++iter)
   {
     // get the connectivity of this entity
-    result = thisMB->get_connectivity(*iter, conn, num_nodes);
-    assert(MB_SUCCESS == result);
+    result = thisMB->get_connectivity(*iter, tmp_conn, num_nodes, true);
+    if (MB_SUCCESS == result)
+      conn = tmp_conn;
+    else {
+        // that didn't work, possibly because it's a structured mesh
+        // which doesn't store connectivity explicitly; use a connect
+        // vector instead
+      tmp_conn_vec.clear();
+      result = thisMB->get_connectivity(&(*iter), 1, tmp_conn_vec, true);
+      if (MB_SUCCESS != result) return MB_FAILURE;
+      conn = &tmp_conn_vec[0];
+      num_nodes = tmp_conn_vec.size();
+    }
     
     type = thisMB->type_from_handle(*iter);
     MBRange::iterator seek_iter;

Modified: MOAB/trunk/ReadNCDF.cpp
===================================================================
--- MOAB/trunk/ReadNCDF.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/ReadNCDF.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -411,6 +411,8 @@
                                   const int *blocks_to_load,
                                   const int num_blocks) 
 {
+  MBErrorCode status;
+  
   file_set = 0;
     // this function directs the reading of an exoii file, but doesn't do any of
     // the actual work
@@ -419,7 +421,7 @@
   reset();
   bool previously_loaded = false;
   std::string filename( exodus_file_name );
-  MBErrorCode status = check_file_status(filename, previously_loaded);
+  status = check_file_status(filename, previously_loaded);
   if (MB_SUCCESS != status) 
     return status;
 
@@ -427,6 +429,9 @@
   status = read_exodus_header(exodus_file_name);
   if (MB_FAILURE == status) return status;
   
+  status = mdbImpl->get_entities_by_handle(0, initRange);
+  if (MB_FAILURE == status) return status;
+
     // 2. Read the nodes unless they've already been read before
   if (!previously_loaded)
   {
@@ -461,6 +466,14 @@
     if (MB_FAILURE == status) return status;
   }
 
+  
+  MBRange loaded_range;
+  status = mdbImpl->get_entities_by_handle(0, loaded_range);
+  if (MB_FAILURE == status) return status;
+  loaded_range = loaded_range.subtract(initRange);
+  status = mdbImpl->add_entities(mCurrentMeshHandle, loaded_range);
+  if (MB_FAILURE == status) return status;
+  
     // what about properties???
 
   file_set = mCurrentMeshHandle;
@@ -600,17 +613,7 @@
     }
   }
 
-    // loaded successfully; add these nodes to loading set
-  int ierr = 0;
-  MBEntityHandle end_handle = CREATE_HANDLE(MBVERTEX, 
-                                            ID_FROM_HANDLE(node_handle)+
-                                            numberNodes_loading-1,
-                                            ierr);
-  if (0 != ierr) return MB_FAILURE;
-  MBErrorCode result = mdbImpl->add_entities(mCurrentMeshHandle,
-                                             MBRange(node_handle, end_handle));
-  
-  return result;
+  return MB_SUCCESS;
 }
 
 MBErrorCode ReadNCDF::read_block_headers(const int *blocks_to_load,
@@ -736,9 +739,11 @@
 
   //get all the blocks of mCurrentMeshHandle
   MBRange child_meshsets; 
-  if(mdbImpl->get_entities_by_type(mCurrentMeshHandle, MBENTITYSET, child_meshsets ) != MB_SUCCESS )
+  if(mdbImpl->get_entities_by_type(0, MBENTITYSET, child_meshsets ) != MB_SUCCESS )
     return MB_FAILURE;
 
+  child_meshsets = child_meshsets.subtract(initRange);
+
   MBTag tag_handle;
 
   MBRange::iterator iter, end_iter;
@@ -881,10 +886,6 @@
       if( mdbImpl->add_entities( ms_handle, new_range) != MB_SUCCESS )
         return MB_FAILURE;
 
-      //add the meshset to the mCurrentMeshHandle
-      if( mdbImpl->add_entities( mCurrentMeshHandle, &ms_handle, 1 ) != MB_SUCCESS )
-        return MB_FAILURE;
-      
     // just a check because the following code won't work if this case fails
     assert(sizeof(MBEntityHandle) >= sizeof(int));
 
@@ -923,9 +924,6 @@
     if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS )
       return MB_FAILURE;
 
-      // success; add new elements to file set
-    if (mdbImpl->add_entities( mCurrentMeshHandle, new_range) != MB_SUCCESS)
-      return MB_FAILURE;
   }
 
   delete [] temp_string;
@@ -941,6 +939,7 @@
   NcVar *temp_var = ncFile->get_var("elem_map");
   if (NULL == temp_var || !temp_var->is_valid()) {
     readMeshIface->report_error("ReadNCDF:: Problem getting element number map variable.");
+    delete [] ptr;
     return MB_FAILURE;
   }
   NcBool status = temp_var->get(ptr, numberElements_loading);
@@ -951,28 +950,22 @@
   }
 
   std::vector<ReadBlockData>::iterator iter;
+  int ptr_pos = 0;
   for(iter = blocksLoading.begin(); iter != blocksLoading.end(); ++iter)
   {
     if (iter->reading_in)
     {
-      MBEntityHandle start_handle = iter->startMBId;
-      
-      if (start_handle != 0)
+      if (iter->startMBId != 0)
       {
-        int i;
-        for (i = 0; i < iter->numElements; i++)
+        MBRange range(iter->startMBId, iter->startMBId+iter->numElements-1);
+        MBErrorCode error = mdbImpl->tag_set_data(mGlobalIdTag, 
+                                                  range, &ptr[ptr_pos]);
+        if (error != MB_SUCCESS)
         {
-          // TODO:  The fact that this must be done one entity at a time is
-          // a tremendous bottle neck in the code.  We need a way of setting
-          // multiple handles and data.
-          MBEntityHandle handle = start_handle + i;
-          MBErrorCode error = mdbImpl->tag_set_data(mGlobalIdTag, &handle, 1, &(ptr[i]));
-          if (error != MB_SUCCESS)
-          {
-            delete [] ptr;
-            return error;
-          }
+          delete [] ptr;
+          return error;
         }
+        ptr_pos += iter->numElements;
       }
       else
       {
@@ -982,9 +975,34 @@
     }
   }
 
+    // read in node map next
+  if (numberNodes_loading > numberElements_loading) {
+    delete [] ptr;
+    ptr = new int [numberNodes_loading];
+  }
+
+  temp_var = ncFile->get_var("node_num_map");
+  if (NULL == temp_var || !temp_var->is_valid()) {
+    readMeshIface->report_error("ReadNCDF:: Problem getting node number map variable.");
+    delete [] ptr;
+    return MB_FAILURE;
+  }
+  status = temp_var->get(ptr, numberNodes_loading);
+  if (0 == status) {
+    readMeshIface->report_error("ReadNCDF:: Problem getting node number map data.");
+    delete [] ptr;
+    return MB_FAILURE;
+  }
+
+  MBRange range(1+vertexOffset, 1+vertexOffset+numberNodes_loading-1);
+  MBErrorCode error = mdbImpl->tag_set_data(mGlobalIdTag, 
+                                            range, &ptr[0]);
+  if (MB_SUCCESS != error)
+    readMeshIface->report_error("ReadNCDF:: Problem setting node global ids.");
+
   delete [] ptr;
   
-  return MB_SUCCESS;
+  return error;
 }
 
 MBErrorCode ReadNCDF::read_nodesets() 
@@ -1078,10 +1096,11 @@
 
       // Maybe there is already a nodesets meshset here we can append to 
     MBRange child_meshsets;
-    if( mdbImpl->get_entities_by_handle( mCurrentMeshHandle, 
-                                         child_meshsets ) != MB_SUCCESS ) 
+    if( mdbImpl->get_entities_by_handle(0, child_meshsets ) != MB_SUCCESS ) 
       return MB_FAILURE;
 
+    child_meshsets = child_meshsets.subtract(initRange);
+
     MBRange::iterator iter, end_iter;
     iter = child_meshsets.begin();
     end_iter = child_meshsets.end();
@@ -1140,8 +1159,6 @@
     {
       if( mdbImpl->create_meshset( MESHSET_ORDERED | MESHSET_TRACK_OWNER, ns_handle ) != MB_SUCCESS) 
         return MB_FAILURE;
-      if( mdbImpl->add_entities( mCurrentMeshHandle, &ns_handle, 1 ) != MB_SUCCESS )
-        return MB_FAILURE;
 
         // set a tag signifying dirichlet bc
         // TODO: create this tag another way
@@ -1228,10 +1245,12 @@
 
   // Maybe there is already a sidesets meshset here we can append to 
   MBRange child_meshsets;
-  if( mdbImpl->get_entities_by_type( mCurrentMeshHandle, MBENTITYSET, 
-                                     child_meshsets ) != MB_SUCCESS ) 
+  if( mdbImpl->get_entities_by_type(0, MBENTITYSET, 
+                                    child_meshsets ) != MB_SUCCESS ) 
     return MB_FAILURE;
 
+  child_meshsets = child_meshsets.subtract(initRange);
+
   MBRange::iterator iter, end_iter;
 
   int i;
@@ -1323,8 +1342,6 @@
           return MB_FAILURE;
         if( mdbImpl->tag_set_data(mGlobalIdTag, &ss_handle, 1, &sideset_id ) != MB_SUCCESS)
           return MB_FAILURE;
-        if( mdbImpl->add_entities( mCurrentMeshHandle, &ss_handle, 1 ) != MB_SUCCESS )
-          return MB_FAILURE;
 
         if (!reverse_entities.empty()) {
             // also make a reverse set to put in this set
@@ -1332,6 +1349,7 @@
           if (mdbImpl->create_meshset(MESHSET_SET | MESHSET_TRACK_OWNER, reverse_set) != MB_SUCCESS)
             return MB_FAILURE;
 
+
             // add the reverse set to the sideset set and the entities to the reverse set
           MBErrorCode result = mdbImpl->add_entities(ss_handle, &reverse_set, 1);
           if (MB_SUCCESS != result) return result;
@@ -1678,7 +1696,7 @@
 }
 
 MBErrorCode ReadNCDF::create_sideset_element( std::vector<MBEntityHandle> connectivity, 
-                                                     MBEntityType type, MBEntityHandle& handle)
+                                              MBEntityType type, MBEntityHandle& handle)
 {
   // get adjacent entities
   MBErrorCode error = MB_SUCCESS;
@@ -1746,9 +1764,7 @@
 
   // if we didn't find a match, create an element
   if(match_found == false)
-  {
     error = mdbImpl->create_element(type, &connectivity[0], connectivity.size(), handle);
-  }
 
   return error;
 }

Modified: MOAB/trunk/ReadNCDF.hpp
===================================================================
--- MOAB/trunk/ReadNCDF.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/ReadNCDF.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -40,6 +40,8 @@
 #include "MBForward.hpp"
 #include "MBReaderIface.hpp"
 #include "ExoIIInterface.hpp"
+#include "MBRange.hpp"
+
 class MBReadUtilIface;
 
 struct ReadBlockData
@@ -214,6 +216,9 @@
   MBTag mQaRecordTag;
 
   int max_line_length, max_str_length;
+
+    //! range of entities in initial mesh, before this read
+  MBRange initRange;
 };
 
 // inline functions

Modified: MOAB/trunk/SparseTagCollections.cpp
===================================================================
--- MOAB/trunk/SparseTagCollections.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/SparseTagCollections.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -36,8 +36,6 @@
   SparseTagSuperCollection functions -----------------------------
 */
 
-#define get_collection( A ) ((A) < mDataTags.size() ? mDataTags[(A)] : 0)
-
 SparseTagSuperCollection::~SparseTagSuperCollection()
 {
   std::vector<SparseTagCollection*>::iterator tag_iterator;

Modified: MOAB/trunk/SparseTagCollections.hpp
===================================================================
--- MOAB/trunk/SparseTagCollections.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/SparseTagCollections.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -44,8 +44,10 @@
 
 #include "MBTypes.h"
 #include "MBInternals.hpp"
-class MBRange;
+#include "MBRange.hpp"
 
+#define get_collection( A ) ((A) < mDataTags.size() ? mDataTags[(A)] : 0)
+
 //! allocator for tag data
 class SparseTagDataAllocator
 {
@@ -94,6 +96,9 @@
   //! gets all entity handles that match a type and tag
   MBErrorCode get_entities(MBEntityType type, MBRange &entities);
 
+  //! gets all entity handles that match a tag
+  MBErrorCode get_entities(MBRange &entities) const;
+
   //! gets all entity handles that match a type, tag, tag_value
   MBErrorCode get_entities_with_tag_value(MBEntityType type, 
                                            MBRange &entities, 
@@ -125,6 +130,15 @@
   return (mData.find(entity) == mData.end() ? false : true);
 }
 
+inline MBErrorCode SparseTagCollection::get_entities(MBRange &entities) const 
+{
+  for (std::map<MBEntityHandle,void*>::const_iterator mit = mData.begin();
+       mit != mData.end(); mit++) 
+    entities.insert((*mit).first);
+
+  return MB_SUCCESS;
+}
+
 //! a collection of SparseTagCollections
 class SparseTagSuperCollection
 {
@@ -158,6 +172,9 @@
   //! finds the entity handle with this data
   MBEntityHandle find_entity( const MBTagId tag_handle, const void* data );
 
+  //! gets all entity handles that match a type and tag
+  MBErrorCode get_entities(const MBTagId tag_handle, MBRange &entities);
+
   //! gets all entity handles that match a tag
   MBErrorCode get_entities(const MBTagId tag_handle, const MBEntityType type,
                             MBRange &entities);
@@ -199,7 +216,15 @@
   std::vector<SparseTagCollection*> mDataTags;
 };
 
+//! gets all entity handles that match a type and tag
+inline MBErrorCode SparseTagSuperCollection::get_entities(const MBTagId tag_handle, 
+                                                          MBRange &entities)
+{
+  SparseTagCollection* coll = get_collection(tag_handle);
+  return coll ? coll->get_entities(entities) : MB_TAG_NOT_FOUND;
+}
 
+
 #endif //SPARSE_TAG_COLLECTIONS_HPP
 
 

Modified: MOAB/trunk/TagServer.cpp
===================================================================
--- MOAB/trunk/TagServer.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/TagServer.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -726,6 +726,29 @@
   return result;
 }
 
+//! gets all entity handles that match a tag
+MBErrorCode TagServer::get_entities(const MBTag tag_handle, MBRange &entities)
+{
+  MBErrorCode result = MB_TAG_NOT_FOUND;
+  MBTagId id = ID_FROM_TAG_HANDLE(tag_handle);
+  switch (PROP_FROM_TAG_HANDLE(tag_handle)) {
+    case MB_TAG_SPARSE:
+      result = mSparseData->get_entities(id, entities);
+      break;
+    case MB_TAG_DENSE:
+      result = mDenseData->get_entities(id, entities);
+      break;
+    case MB_TAG_BIT:
+      result = mBitServer->get_entities(id, entities);
+      break;
+    case MB_TAG_MESH:
+      result = MB_TYPE_OUT_OF_RANGE;
+      break;
+  }
+  
+  return result;
+}
+
 //! gets all entity handles that match a type and tag
 MBErrorCode TagServer::get_entities(const MBRange &range,
                                      const MBTag tag_handle, const MBEntityType type,

Modified: MOAB/trunk/TagServer.hpp
===================================================================
--- MOAB/trunk/TagServer.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/TagServer.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -201,6 +201,10 @@
                              const MBEntityType type,
                              MBRange &entities);
 
+  //! gets all entity handles that match a tag
+  MBErrorCode get_entities( const MBTag tag_handle, 
+                             MBRange &entities);
+
   //! gets all entity handles that match a type and tag
   MBErrorCode get_entities( const MBRange &input_range,
                              const MBTag tag_handle, 

Modified: MOAB/trunk/mbparallelcomm_test.cpp
===================================================================
--- MOAB/trunk/mbparallelcomm_test.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/mbparallelcomm_test.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -10,11 +10,30 @@
 #include "MBParallelConventions.h"
 #include "MBTagConventions.hpp"
 #include "MBCore.hpp"
+#include "ScdVertexSeq.hpp"
+#include "ScdElementSeq.hpp"
+#include "EntitySequenceManager.hpp"
 #include "mpi.h"
 #include <iostream>
 
+#define REALTFI 1
+
 #define ERROR(a, b) {std::cerr << a << std::endl; return b;}
 
+#define PRINT_LAST_ERROR {\
+        std::string last_error;\
+        result = mbImpl->get_last_error(last_error);\
+        if (last_error.empty()) std::cerr << "(none)" << std::endl;\
+        else std::cerr << last_error << std::endl;\
+        }
+MBErrorCode create_linear_mesh(MBInterface *mbImpl,
+                               int N, int M, int &nshared);
+
+MBErrorCode create_scd_mesh(MBInterface *mbImpl,
+                            int IJK, int &nshared);
+
+MBErrorCode read_file(MBInterface *mbImpl, const char *filename);
+
 int main(int argc, char **argv) 
 {
     // need to init MPI first, to tell how many procs and rank
@@ -45,88 +64,453 @@
     N = atoi(argv[1]);
     M = atoi(argv[2]);
   }
+
+  int max_iter = 2;
+  if (argc > 3) max_iter = atoi(argv[3]);
+
+  for (int i = 0; i < max_iter; i++) {
+    int nshared;
+    MBErrorCode tmp_result = MB_SUCCESS;
+    if (0 == i) {
+      tmp_result = create_linear_mesh(mbImpl, N, M, nshared);
+      if (MB_SUCCESS != tmp_result) {
+        result = tmp_result;
+        std::cerr << "Couldn't create linear mesh; error message:." 
+                  << std::endl;
+        PRINT_LAST_ERROR
+        continue;
+      }
+    }
+    else if (1 == i) {
+      tmp_result = create_scd_mesh(mbImpl, N, nshared);
+      if (MB_SUCCESS != tmp_result) {
+        result = tmp_result;
+        std::cerr << "Couldn't create structured mesh; error message:" 
+                  << std::endl;
+        PRINT_LAST_ERROR
+        continue;
+      }
+    }
+    else if (2 == i && argc > 4) {
+        // read a file in parallel from the filename on the command line
+
+      tmp_result = read_file(mbImpl, argv[4]);
+      if (MB_SUCCESS != tmp_result) {
+        result = tmp_result;
+        std::cerr << "Couldn't read mesh; error message:" << std::endl;
+        PRINT_LAST_ERROR
+        continue;
+      }
+      nshared = -1;
+    }
+
+    if (MB_SUCCESS == tmp_result) {
+        // now figure out which vertices are shared
+      double wtime;
+      if (0 == rank) wtime = MPI_Wtime();
+      MBParallelComm *pcomm = new MBParallelComm(mbImpl);
+      tmp_result = pcomm->resolve_shared_ents();
+      if (MB_SUCCESS != tmp_result) {
+        std::cerr << "Couldn't resolve shared entities." << std::endl;
+        result = tmp_result;
+        continue;
+      }
+      
+      if (0 == rank) wtime = MPI_Wtime() - wtime;
+
+      MBRange shared_ents;
+      tmp_result = pcomm->get_shared_entities(0, shared_ents);
+      
+      if (MB_SUCCESS != tmp_result) {
+        std::cerr << "get_shared_entities returned error on proc " 
+                  << rank << "; message: " << std::endl;
+        PRINT_LAST_ERROR
+        result = tmp_result;
+      }
   
+        // check # shared entities
+      else if (0 <= nshared && nshared != (int) shared_ents.size()) {
+        std::cerr << "Didn't get correct number of shared vertices on "
+                  << "processor " << rank << std::endl;
+        result = MB_FAILURE;
+      }
+      else if (i < 2) {
+        std::cerr << "Proc " << rank;
+        if (0 == i) std::cerr << " linear mesh succeeded." << std::endl;
+        else std::cerr << " structured mesh succeeded." << std::endl;
+        if (0 == rank) std::cerr << "   Time = " << wtime << "." << std::endl;
+      }
+      else {
+        std::cerr << "Proc " << rank << " " << shared_ents.size()
+                  << " shared entities." << std::endl;
+      }
+  
+      delete pcomm;
+      tmp_result = mbImpl->delete_mesh();
+      if (MB_SUCCESS != tmp_result) {
+        result = tmp_result;
+        std::cerr << "Couldn't delete mesh on rank " << rank
+                  << "; error message: " << std::endl;
+        PRINT_LAST_ERROR
+      }
+    }
+  }
+  
+  
+  err = MPI_Finalize();
+
+  if (MB_SUCCESS == result)
+    std::cerr << "Proc " << rank << ": Success." << std::endl;
+    
+  return (MB_SUCCESS == result ? 0 : 1);
+}
+
+MBErrorCode read_file(MBInterface *mbImpl, const char *filename) 
+{
+  std::string options = "PARALLEL=BCAST_DELETE;PARTITION=MATERIAL_SET";
+  MBEntityHandle file_set;
+  MBErrorCode result = mbImpl->load_file(filename, file_set, 
+                                         options.c_str());
+  return result;
+}
+
+#define RR(a, b) {std::cerr << a; return b;}
+
+MBErrorCode create_linear_mesh(MBInterface *mbImpl,
+                               int N, int M, int &nshared) 
+{
+    /* create a mesh where each processor owns N vertices and shares
+     * M vertices with processors on either end; for non-root procs,
+     * that works out to N+M vertices.
+     *
+     * Number of vertices shared should be M on end procs, 2M on others
+     */
+  int my_rank = mbImpl->proc_rank();
+  int my_size = mbImpl->proc_size();
+  
   int nverts = N + M;
-  if (0 == rank) nverts = N;
+  if (0 == my_rank) nverts = N;
 
     // create my vertices and give them the right global ids
   MBRange my_verts;
   std::vector<double> coords(3*(nverts));
   std::fill(coords.begin(), coords.end(), 0.0);
-  result = mbImpl->create_vertices(&coords[0], nverts,
-                                   my_verts);
+  MBErrorCode result = mbImpl->create_vertices(&coords[0], nverts,
+                                               my_verts);
   if (MB_SUCCESS != 0)
-    ERROR("Failed to create vertices.", 1);
+    RR("Failed to create vertices.", MB_FAILURE);
   
   std::vector<int> global_ids(N+M);
   for (int i = 0; i < nverts; i++)
-    global_ids[i] = rank*N - (nverts-N) + i;
+    global_ids[i] = my_rank*N - (nverts-N) + i;
   
   int def_val = -1;
   MBTag gid_tag;
-  result = mbImpl->tag_create(GLOBAL_ID_TAG_NAME, 1, MB_TAG_DENSE,
+  result = mbImpl->tag_create(GLOBAL_ID_TAG_NAME, sizeof(int), MB_TAG_DENSE,
                               MB_TYPE_INTEGER, gid_tag,
                               &def_val, true);
   if (MB_SUCCESS != result && MB_ALREADY_ALLOCATED != result) 
-    ERROR("Failed to create tag.", 1);
+    RR("Failed to create tag.", MB_FAILURE);
   
   result = mbImpl->tag_set_data(gid_tag, my_verts, &global_ids[0]);
-  if (MB_SUCCESS != result) ERROR("Failed to set global_id tag.", 1);
+  if (MB_SUCCESS != result) RR("Failed to set global_id tag.", MB_FAILURE);
+
+  nshared = ((my_rank == 0 || my_rank == my_size-1) ? M : 2*M);
   
-    // now figure out what's shared
-  MBParallelComm pcomm(mbImpl);
-  result = pcomm.resolve_shared_ents(my_verts, 0);
-  if (MB_SUCCESS != result) ERROR("Couldn't resolve shared entities.", 1);
+  return MB_SUCCESS;
+}
+
+double LENGTH = 1.0;
+
+void build_coords(const int nelem, double *&coords);
+
+MBErrorCode create_scd_mesh(MBInterface *mbImpl,
+                            int IJK, int &nshared) 
+{
+    /* Create a 3d mesh of hexes, parameterized in i, j, k, where
+     * on processor of rank R the K parameterization starts at 
+     * i=0, j=0, k=RK, and there are I+1, J+1, K+1 vertices in the
+     * i, j, k directions, resp.; for now, make I, J, K equal
+     *
+     * Should share (I+1)(J+1) vertices on end procs, or twice that
+     * on interior procs.
+     */
+  int my_rank = mbImpl->proc_rank();
   
-    // check shared entities
-  MBTag sharedproc_tag, sharedprocs_tag;
-  result = mbImpl->tag_get_handle(PARALLEL_SHARED_PROC_TAG_NAME, 
-                                  sharedproc_tag);
-  if (MB_SUCCESS != result) ERROR("Shared processor tag not found.", 1);
+    // make a 3d block of vertices
+  MBEntitySequence *dum_seq = NULL;
+  ScdVertexSeq *vseq = NULL;
+  ScdElementSeq *eseq = NULL;
+  EntitySequenceManager *seq_mgr = 
+    dynamic_cast<MBCore*>(mbImpl)->sequence_manager();
+  HomCoord vseq_minmax[2] = {HomCoord(0,0,0), 
+                             HomCoord(IJK, IJK, IJK)};
+  MBEntityHandle vstart, estart;
+  
+  MBErrorCode result = 
+    seq_mgr->create_scd_sequence(vseq_minmax[0], vseq_minmax[1],
+                                 MBVERTEX, 1, vstart, dum_seq);
+  if (NULL != dum_seq) vseq = dynamic_cast<ScdVertexSeq*>(dum_seq);
+  assert (MB_FAILURE != result && vstart != 0 && dum_seq != NULL && 
+          vseq != NULL);
+    // now the element sequence
+  result = seq_mgr->create_scd_sequence(vseq_minmax[0], vseq_minmax[1], 
+                                        MBHEX, 1, estart, dum_seq);
+  if (NULL != dum_seq) eseq = dynamic_cast<ScdElementSeq*>(dum_seq);
+  assert (MB_FAILURE != result && estart != 0 && dum_seq != NULL && 
+          eseq != NULL);
+  
+    // only need to add one vseq to this, unity transform
+    // trick: if I know it's going to be unity, just input 3 sets of equivalent points
+  result = eseq->add_vsequence(vseq, vseq_minmax[0], vseq_minmax[0], 
+                               vseq_minmax[0], 
+                               vseq_minmax[0], vseq_minmax[0], vseq_minmax[0]);
+  assert(MB_SUCCESS == result);
 
-  result = mbImpl->tag_get_handle(PARALLEL_SHARED_PROCS_TAG_NAME, 
-                                  sharedprocs_tag);
-  if (MB_SUCCESS != result) 
-    ERROR("Shared processor*s* tag not found.", 1);
+    // set the coordinates of the vertices
+  double *coords = NULL;
+  build_coords(IJK, coords);
+
+    // offset coords by starting z distance
+  double deltaz = my_rank*LENGTH;
+  int num_verts = (IJK + 1)*(IJK + 1)*(IJK + 1);
+  for (int k = 0; k < num_verts; k++)
+    coords[3*k+2] += deltaz;
   
-    // get the tag values
-#define MAX_SHARING_PROCS 10
-  std::vector<int> shared_proc_tags(MAX_SHARING_PROCS*my_verts.size());
-  result = mbImpl->tag_get_data(sharedproc_tag, my_verts, 
-                                &shared_proc_tags[0]);
-  if (MB_SUCCESS != result) ERROR("Problem getting shared proc tag.", 1);
+  MBRange vrange(vstart, vstart+num_verts-1);
+  result = mbImpl->set_coords(vrange, coords);
+  if (MB_SUCCESS != result) RR("Couldn't build coords array.", MB_FAILURE);
+  
+    // set global ids for vertices; reuse coords space
+  int *gids = (int*) coords;
+  int start_gid = 1 + my_rank * (IJK)*(IJK+1)*(IJK+1);
+  for (int i = 0; i < num_verts; i++)
+    gids[i] = start_gid++;
+  
+  int def_val = -1;
+  MBTag gid_tag;
+  result = mbImpl->tag_create(GLOBAL_ID_TAG_NAME, sizeof(int), MB_TAG_DENSE,
+                              MB_TYPE_INTEGER, gid_tag,
+                              &def_val, true);
+  if (MB_SUCCESS != result && MB_ALREADY_ALLOCATED != result) 
+    RR("Failed to create tag.", MB_FAILURE);
+  
+  result = mbImpl->tag_set_data(gid_tag, vrange, gids);
+  if (MB_SUCCESS != result) RR("Failed to set global_id tag.", MB_FAILURE);
 
-    // interior procs should have 2*M shared, bdy procs should have M shared
-  int nshared = 0;
-  for (unsigned int nv = 0; nv < my_verts.size(); nv++)
-    if (shared_proc_tags[2*nv] > -1) nshared++;
+  nshared = (IJK+1) * (IJK+1);
   
-  if ((rank == 0 || rank == nprocs-1) && nshared != (unsigned int) M) {
-    std::cerr << "Didn't get correct number of shared vertices on "
-              << "processor " << rank << std::endl;
-    result = MB_FAILURE;
+  return result;
+}
+
+void compute_edge(double *start, const int nelem,  const double xint,
+                  const int stride) 
+{
+  for (int i = 1; i < nelem; i++) {
+    start[i*stride] = start[0]+i*xint;
+    start[nelem+1+i*stride] = start[nelem+1]+i*xint;
+    start[2*(nelem+1)+i*stride] = start[2*(nelem+1)]+i*xint;
   }
+}
+
+void compute_face(double *a, const int nelem,  const double xint,
+                  const int stride1, const int stride2) 
+{
+    // 2D TFI on a face starting at a, with strides stride1 in ada and stride2 in tse
+  for (int j = 1; j < nelem; j++) {
+    double tse = j * xint;
+    for (int i = 1; i < nelem; i++) {
+      double ada = i * xint;
+      
+      a[i*stride1+j*stride2] = (1.0 - ada)*a[i*stride1]
+        + ada*a[i*stride1+nelem*stride2]
+        + (1.0 - tse)*a[j*stride2]
+        + tse*a[j*stride2+nelem*stride1]
+        - (1.0 - tse)*(1.0 - ada)*a[0]
+        - (1.0 - tse)*ada*a[nelem*stride1]
+        - tse*(1.0 - ada)*a[nelem*stride2]
+        - tse*ada*a[nelem*(stride1+stride2)];
+      a[nelem+1+i*stride1+j*stride2] = (1.0 - ada)*a[nelem+1+i*stride1]
+        + ada*a[nelem+1+i*stride1+nelem*stride2]
+        + (1.0 - tse)*a[nelem+1+j*stride2]
+        + tse*a[nelem+1+j*stride2+nelem*stride1]
+        - (1.0 - tse)*(1.0 - ada)*a[nelem+1+0]
+        - (1.0 - tse)*ada*a[nelem+1+nelem*stride1]
+        - tse*(1.0 - ada)*a[nelem+1+nelem*stride2]
+        - tse*ada*a[nelem+1+nelem*(stride1+stride2)];
+      a[2*(nelem+1)+i*stride1+j*stride2] = (1.0 - ada)*a[2*(nelem+1)+i*stride1]
+        + ada*a[2*(nelem+1)+i*stride1+nelem*stride2]
+        + (1.0 - tse)*a[2*(nelem+1)+j*stride2]
+        + tse*a[2*(nelem+1)+j*stride2+nelem*stride1]
+        - (1.0 - tse)*(1.0 - ada)*a[2*(nelem+1)+0]
+        - (1.0 - tse)*ada*a[2*(nelem+1)+nelem*stride1]
+        - tse*(1.0 - ada)*a[2*(nelem+1)+nelem*stride2]
+        - tse*ada*a[2*(nelem+1)+nelem*(stride1+stride2)];
+    }
+  }
+}
+
+void build_coords(const int nelem, double *&coords) 
+{
+    // allocate the memory
+  int numv = nelem+1;
+  int numv_sq = numv*numv;
+  int tot_numv = numv*numv*numv;
+  coords = new double[3*tot_numv];
+
+// use FORTRAN-like indexing
+#define VINDEX(i,j,k) (i + (j*numv) + (k*numv_sq))
+  int idx;
+  double scale1, scale2, scale3;
+    // use these to prevent optimization on 1-scale, etc (real map wouldn't have
+    // all these equal)
+  scale1 = LENGTH/nelem;
+  scale2 = LENGTH/nelem;
+  scale3 = LENGTH/nelem;
+
+#ifdef REALTFI
+    // use a real TFI xform to compute coordinates
+    // compute edges
+    // i (stride=1)
+  compute_edge(&coords[VINDEX(0,0,0)], nelem, scale1, 1);
+  compute_edge(&coords[VINDEX(0,nelem,0)], nelem, scale1, 1);
+  compute_edge(&coords[VINDEX(0,0,nelem)], nelem, scale1, 1);
+  compute_edge(&coords[VINDEX(0,nelem,nelem)], nelem, scale1, 1);
+    // j (stride=numv)
+  compute_edge(&coords[VINDEX(0,0,0)], nelem, scale1, numv);
+  compute_edge(&coords[VINDEX(nelem,0,0)], nelem, scale1, numv);
+  compute_edge(&coords[VINDEX(0,0,nelem)], nelem, scale1, numv);
+  compute_edge(&coords[VINDEX(nelem,0,nelem)], nelem, scale1, numv);
+    // k (stride=numv^2)
+  compute_edge(&coords[VINDEX(0,0,0)], nelem, scale1, numv_sq);
+  compute_edge(&coords[VINDEX(nelem,0,0)], nelem, scale1, numv_sq);
+  compute_edge(&coords[VINDEX(0,nelem,0)], nelem, scale1, numv_sq);
+  compute_edge(&coords[VINDEX(nelem,nelem,0)], nelem, scale1, numv_sq);
+
+    // compute faces
+    // i=0, nelem
+  compute_face(&coords[VINDEX(0,0,0)], nelem, scale1, numv, numv_sq);
+  compute_face(&coords[VINDEX(nelem,0,0)], nelem, scale1, numv, numv_sq);
+    // j=0, nelem
+  compute_face(&coords[VINDEX(0,0,0)], nelem, scale1, 1, numv_sq);
+  compute_face(&coords[VINDEX(0,nelem,0)], nelem, scale1, 1, numv_sq);
+    // k=0, nelem
+  compute_face(&coords[VINDEX(0,0,0)], nelem, scale1, 1, numv);
+  compute_face(&coords[VINDEX(0,0,nelem)], nelem, scale1, 1, numv);
+
+    // initialize corner indices
+  int i000 = VINDEX(0,0,0);
+  int ia00 = VINDEX(nelem,0,0);
+  int i0t0 = VINDEX(0,nelem,0);
+  int iat0 = VINDEX(nelem,nelem,0);
+  int i00g = VINDEX(0,0,nelem);
+  int ia0g = VINDEX(nelem,0,nelem);
+  int i0tg = VINDEX(0,nelem,nelem);
+  int iatg = VINDEX(nelem,nelem,nelem);
+  double cX, cY, cZ;
+  int adaInts = nelem;
+  int tseInts = nelem;
+  int gammaInts = nelem;
   
-  else if ((rank != 0 && rank != nprocs-1) && nshared != (unsigned int) 2*M) 
-  {
-    std::cerr << "Didn't get correct number of shared vertices on "
-              << "processor " << rank << std::endl;
-    result = MB_FAILURE;
+  
+  for (int i=1; i < nelem; i++) {
+    for (int j=1; j < nelem; j++) {
+      for (int k=1; k < nelem; k++) {
+        idx = VINDEX(i,j,k);
+        double tse = i*scale1;
+        double ada = j*scale2;
+        double gamma = k*scale3;
+        double tm1 = 1.0 - tse;
+        double am1 = 1.0 - ada;
+        double gm1 = 1.0 - gamma;
+
+        cX = gm1 *   (am1*(tm1*coords[i000] + tse*coords[i0t0])  +
+                      ada*(tm1*coords[ia00] + tse*coords[iat0])) +
+          gamma * (am1*(tm1*coords[i00g] + tse*coords[i0tg])  +
+                   ada*(tm1*coords[ia0g] + tse*coords[iatg]));
+
+        cY = gm1 *   (am1*(tm1*coords[i000] + tse*coords[i0t0])  +
+                      ada*(tm1*coords[ia00] + tse*coords[iat0])) +
+          gamma * (am1*(tm1*coords[i00g] + tse*coords[i0tg])  +
+                   ada*(tm1*coords[ia0g] + tse*coords[iatg]));
+
+        cZ = gm1 *   (am1*(tm1*coords[i000] + tse*coords[i0t0])  +
+                      ada*(tm1*coords[ia00] + tse*coords[iat0])) +
+          gamma * (am1*(tm1*coords[i00g] + tse*coords[i0tg])  +
+                   ada*(tm1*coords[ia0g] + tse*coords[iatg]));
+
+        double *ai0k = &coords[VINDEX(k,0,i)];
+        double *aiak = &coords[VINDEX(k,adaInts,i)];
+        double *a0jk = &coords[VINDEX(k,j,0)];
+        double *atjk = &coords[VINDEX(k,j,tseInts)];
+        double *aij0 = &coords[VINDEX(0,j,i)];
+        double *aijg = &coords[VINDEX(gammaInts,j,i)];
+  
+        coords[VINDEX(i,j,k)] = (   am1*ai0k[0] 
+                                    + ada*aiak[0] 
+                                    + tm1*a0jk[0] 
+                                    + tse*atjk[0]
+                                    + gm1*aij0[0] 
+                                    + gamma*aijg[0] )/2.0 - cX/2.0;
+
+        coords[nelem+1+VINDEX(i,j,k)] = (   am1*ai0k[nelem+1] 
+                                            + ada*aiak[nelem+1] 
+                                            + tm1*a0jk[nelem+1] 
+                                            + tse*atjk[nelem+1]
+                                            + gm1*aij0[nelem+1] 
+                                            + gamma*aijg[nelem+1] )/2.0 - cY/2.0;
+
+        coords[2*(nelem+1)+VINDEX(i,j,k)] = (   am1*ai0k[2*(nelem+1)] 
+                                                + ada*aiak[2*(nelem+1)] 
+                                                + tm1*a0jk[2*(nelem+1)] 
+                                                + tse*atjk[2*(nelem+1)]
+                                                + gm1*aij0[2*(nelem+1)] 
+                                                + gamma*aijg[2*(nelem+1)] )/2.0 - cZ/2.0;
+      }
+    }
   }
+  
 
-    // now check sharedprocs; shouldn't be any 
-  MBErrorCode result2 = mbImpl->tag_get_data(sharedprocs_tag, my_verts, 
-                                             &shared_proc_tags[0]);
-  if (MB_SUCCESS == result2) {
-    std::cerr << "Shoudn't get shared proc*s* tag, but did on proc "
-              << rank << std::endl;
-    result = MB_FAILURE;
+#else
+  for (int i=0; i < numv; i++) {
+    for (int j=0; j < numv; j++) {
+      for (int k=0; k < numv; k++) {
+        idx = VINDEX(i,j,k);
+          // blocked coordinate ordering
+        coords[idx] = i*scale1;
+        coords[tot_numv+idx] = j*scale2;
+        coords[2*tot_numv+idx] = k*scale3;
+      }
+    }
   }
+#endif
+}
 
-  err = MPI_Finalize();
+void build_connect(const int nelem, const MBEntityHandle vstart, MBEntityHandle *&connect) 
+{
+    // allocate the memory
+  int nume_tot = nelem*nelem*nelem;
+  connect = new MBEntityHandle[8*nume_tot];
 
-  if (MB_SUCCESS == result)
-    std::cerr << "Proc " << rank << ": Success." << std::endl;
-    
-  return (MB_SUCCESS == result ? 0 : 1);
+  MBEntityHandle vijk;
+  int numv = nelem + 1;
+  int numv_sq = numv*numv;
+  int idx = 0;
+  for (int i=0; i < nelem; i++) {
+    for (int j=0; j < nelem; j++) {
+      for (int k=0; k < nelem; k++) {
+        vijk = vstart+VINDEX(i,j,k);
+        connect[idx++] = vijk;
+        connect[idx++] = vijk+1;
+        connect[idx++] = vijk+1+numv;
+        connect[idx++] = vijk+numv;
+        connect[idx++] = vijk+numv*numv;
+        connect[idx++] = vijk+1+numv*numv;
+        connect[idx++] = vijk+1+numv+numv*numv;
+        connect[idx++] = vijk+numv+numv*numv;
+        assert(i <= numv*numv*numv);
+      }
+    }
+  }
 }
+

Modified: MOAB/trunk/parallel/MBParallelComm.cpp
===================================================================
--- MOAB/trunk/parallel/MBParallelComm.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/parallel/MBParallelComm.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -9,7 +9,10 @@
 #include "MBSkinner.hpp"
 #include "MBParallelConventions.h"
 #include "MBCore.hpp"
+#include "MBError.hpp"
 
+#define MAX_SHARING_PROCS 10  
+
 extern "C" 
 {
 #include "gs.h"
@@ -56,7 +59,9 @@
 #define UNPACK_RANGE(buff, rng) {int num_subs; UNPACK_INTS(buff, &num_subs, 1); MBEntityHandle _eh[2]; \
           for (int i = 0; i < num_subs; i++) { UNPACK_EH(buff_ptr, _eh, 2); rng.insert(_eh[0], _eh[1]);}}
 
-#define RR if (MB_SUCCESS != result) return result
+#define RR(a) if (MB_SUCCESS != result) {\
+          dynamic_cast<MBCore*>(mbImpl)->get_error_handler()->set_last_error(a);\
+          return result;}
 
 MBParallelComm::MBParallelComm(MBInterface *impl, MPI_Comm comm) 
     : mbImpl(impl), procConfig(comm)
@@ -86,7 +91,8 @@
   MBErrorCode result;
   for (int dim = 0; dim <= dimension; dim++) {
     if (dim == 0 || !largest_dim_only || dim == dimension) {
-      result = mbImpl->get_entities_by_dimension(0, dim, entities[dim]); RR;
+      result = mbImpl->get_entities_by_dimension(0, dim, entities[dim]); 
+      RR("Failed to get vertices in assign_global_ids.");
     }
 
       // need to filter out non-locally-owned entities!!!
@@ -135,7 +141,8 @@
     for (MBRange::iterator rit = entities[dim].begin(); rit != entities[dim].end(); rit++)
       num_elements[i++] = total_elems[dim]++;
     
-    result = mbImpl->tag_set_data(gid_tag, entities[dim], &num_elements[0]); RR;
+    result = mbImpl->tag_set_data(gid_tag, entities[dim], &num_elements[0]); 
+    RR("Failed to set global id tag in assign_global_ids.");
   }
   
   return MB_SUCCESS;
@@ -167,7 +174,9 @@
 
     int buff_size;
     
-    result = pack_buffer(entities, adjacencies, tags, true, whole_range, buff_size); RR;
+    result = pack_buffer(entities, adjacencies, tags, true, 
+                         whole_range, buff_size); 
+    RR("Failed to compute buffer size in communicate_entities.");
 
       // if the message is large, send a first message to tell how large
     if (INITIAL_BUFF_SIZE < buff_size) {
@@ -183,7 +192,9 @@
 
       // pack the actual buffer
     int actual_buff_size;
-    result = pack_buffer(entities, adjacencies, tags, false, whole_range, actual_buff_size); RR;
+    result = pack_buffer(entities, adjacencies, tags, false, 
+                         whole_range, actual_buff_size); 
+    RR("Failed to pack buffer in communicate_entities.");
     
       // send it
     MPI_Request send_req;
@@ -212,7 +223,8 @@
     }
     
       // unpack the buffer
-    result = unpack_buffer(entities); RR;
+    result = unpack_buffer(entities); 
+    RR("Failed to unpack buffer in communicate_entities.");
   }
   
   return result;
@@ -244,7 +256,8 @@
   setPcs.clear();
 
   if ((int)procConfig.proc_rank() == from_proc) {
-    result = pack_buffer( entities, adjacencies, tags, true, whole_range, buff_size ); RR;
+    result = pack_buffer( entities, adjacencies, tags, true, whole_range, buff_size ); 
+    RR("Failed to compute buffer size in broadcast_entities.");
   }
 
   success = MPI_Bcast( &buff_size, 1, MPI_INT, from_proc, procConfig.proc_comm() );
@@ -258,7 +271,9 @@
   
   if ((int)procConfig.proc_rank() == from_proc) {
     int actual_buffer_size;
-    result = pack_buffer( entities, adjacencies, tags, false, whole_range, actual_buffer_size ); RR;
+    result = pack_buffer( entities, adjacencies, tags, false, 
+                          whole_range, actual_buffer_size );
+    RR("Failed to pack buffer in broadcast_entities.");
   }
 
   success = MPI_Bcast( &myBuffer[0], buff_size, MPI_UNSIGNED_CHAR, from_proc, procConfig.proc_comm() );
@@ -266,7 +281,8 @@
     return MB_FAILURE;
   
   if ((int)procConfig.proc_rank() != from_proc) {
-    result = unpack_buffer( entities ); RR;
+    result = unpack_buffer( entities );
+    RR("Failed to unpack buffer in broadcast_entities.");
   }
 
   return MB_SUCCESS;
@@ -289,22 +305,29 @@
   if (!just_count) buff_ptr = &myBuffer[0];
   
     // entities
-  result = pack_entities(entities, rit, whole_range, buff_ptr, buff_size, just_count); RR;
+  result = pack_entities(entities, rit, whole_range, buff_ptr, 
+                         buff_size, just_count); 
+  RR("Packing entities failed.");
   
     // sets
   int tmp_size;
-  result = pack_sets(entities, rit, whole_range, buff_ptr, tmp_size, just_count); RR;
+  result = pack_sets(entities, rit, whole_range, buff_ptr, tmp_size, just_count); 
+  RR("Packing sets failed.");
   buff_size += tmp_size;
   
     // adjacencies
   if (adjacencies) {
-    result = pack_adjacencies(entities, rit, whole_range, buff_ptr, tmp_size, just_count); RR;
+    result = pack_adjacencies(entities, rit, whole_range, buff_ptr, 
+                              tmp_size, just_count);
+    RR("Packing adjs failed.");
     buff_size += tmp_size;
   }
     
     // tags
   if (tags) {
-    result = pack_tags(entities, rit, whole_range, buff_ptr, tmp_size, just_count); RR;
+    result = pack_tags(entities, rit, whole_range, buff_ptr, 
+                       tmp_size, just_count);
+    RR("Packing tags failed.");
     buff_size += tmp_size;
   }
 
@@ -316,9 +339,12 @@
   if (myBuffer.capacity() == 0) return MB_FAILURE;
   
   unsigned char *buff_ptr = &myBuffer[0];
-  MBErrorCode result = unpack_entities(buff_ptr, entities); RR;
-  result = unpack_sets(buff_ptr, entities); RR;
-  result = unpack_tags(buff_ptr, entities); RR;
+  MBErrorCode result = unpack_entities(buff_ptr, entities);
+  RR("Unpacking entities failed.");
+  result = unpack_sets(buff_ptr, entities);
+  RR("Unpacking sets failed.");
+  result = unpack_tags(buff_ptr, entities);
+  RR("Unpacking tags failed.");
   
   return MB_SUCCESS;
 }
@@ -346,7 +372,10 @@
   MBErrorCode result;
   MBWriteUtilIface *wu = NULL;
   if (!just_count) {
-    result = mbImpl->query_interface(std::string("MBWriteUtilIface"), reinterpret_cast<void**>(&wu)); RR;
+    result = mbImpl->query_interface(std::string("MBWriteUtilIface"), 
+                                     reinterpret_cast<void**>(&wu));
+    RR("Couldn't get MBWriteUtilIface.");
+
   }
   
     // pack vertices
@@ -365,7 +394,8 @@
 
     assert(NULL != wu);
     
-    result = wu->get_node_arrays(3, num_verts, allRanges[0], 0, 0, coords); RR;
+    result = wu->get_node_arrays(3, num_verts, allRanges[0], 0, 0, coords);
+    RR("Couldn't allocate node space.");
 
     buff_ptr += 3 * num_verts * sizeof(double);
 
@@ -395,7 +425,9 @@
         // get the sequence holding this entity
       MBEntitySequence *seq;
       ElementEntitySequence *eseq;
-      result = sequenceManager->find(*start_rit, seq); RR;
+      result = sequenceManager->find(*start_rit, seq);
+      RR("Couldn't find entity sequence.");
+
       if (NULL == seq) return MB_FAILURE;
       eseq = dynamic_cast<ElementEntitySequence*>(seq);
 
@@ -419,7 +451,9 @@
     }
 
       // update vertex range and count those data, now that we know which entities get communicated
-    result = mbImpl->get_adjacencies(whole_range, 0, false, allRanges[0], MBInterface::UNION); RR;
+    result = mbImpl->get_adjacencies(whole_range, 0, false, allRanges[0], 
+                                     MBInterface::UNION);
+    RR("Failed get_adjacencies.");
     whole_range.merge(allRanges[0]);
     count += 3 * sizeof(double) * allRanges[0].size();
     
@@ -468,7 +502,8 @@
       if (*et_it == MBPOLYGON || *et_it == MBPOLYHEDRON) {
         std::vector<int> num_connects;
         for (MBRange::const_iterator rit = allr_it->begin(); rit != allr_it->end(); rit++) {
-          result = mbImpl->get_connectivity(*rit, connect, num_connect); RR;
+          result = mbImpl->get_connectivity(*rit, connect, num_connect);
+          RR("Failed to get connectivity.");
           num_connects.push_back(num_connect);
           PACK_EH(buff_ptr, &connect[0], num_connect);
         }
@@ -476,7 +511,8 @@
       }
       else {
         for (MBRange::const_iterator rit = allr_it->begin(); rit != allr_it->end(); rit++) {
-          result = mbImpl->get_connectivity(*rit, connect, num_connect); RR;
+          result = mbImpl->get_connectivity(*rit, connect, num_connect);
+          RR("Failed to get connectivity.");
           assert(num_connect == *nv_it);
           PACK_EH(buff_ptr, &connect[0], num_connect);
         }
@@ -500,7 +536,10 @@
   MBErrorCode result;
   bool done = false;
   MBReadUtilIface *ru = NULL;
-  result = mbImpl->query_interface(std::string("MBReadUtilIface"), reinterpret_cast<void**>(&ru)); RR;
+  result = mbImpl->query_interface(std::string("MBReadUtilIface"), 
+                                   reinterpret_cast<void**>(&ru));
+  RR("Failed to get MBReadUtilIface.");
+
   
   while (!done) {
     MBEntityType this_type;
@@ -527,7 +566,9 @@
         MBEntityHandle actual_start;
         int tmp_num_verts = (*pit).second - (*pit).first + 1;
         result = ru->get_node_arrays(3, tmp_num_verts, start_id, start_proc, actual_start,
-                                     coords); RR;
+                                     coords);
+        RR("Failed to allocate node arrays.");
+
         if (actual_start != (*pit).first)
           return MB_FAILURE;
 
@@ -562,14 +603,19 @@
         int num_elems = (*pit).second - (*pit).first + 1;
         MBEntityHandle *connect;
         int *connect_offsets;
-        if (this_type == MBPOLYGON || this_type == MBPOLYHEDRON)
+        if (this_type == MBPOLYGON || this_type == MBPOLYHEDRON) {
           result = ru->get_poly_element_array(num_elems, verts_per_entity, this_type,
                                               start_id, start_proc, actual_start,
-                                              connect_offsets, connect); RR;
-        else
+                                              connect_offsets, connect);
+          RR("Failed to allocate poly element arrays.");
+        }
+
+        else {
           result = ru->get_element_array(num_elems, verts_per_entity, this_type,
                                          start_id, start_proc, actual_start,
-                                         connect); RR;
+                                         connect);
+          RR("Failed to allocate element arrays.");
+        }
 
           // copy connect arrays
         if (this_type != MBPOLYGON && this_type != MBPOLYHEDRON) {
@@ -610,20 +656,24 @@
       count += sizeof(MBEntityHandle);
     
       unsigned int options;
-      result = mbImpl->get_meshset_options(*start_rit, options); RR;
+      result = mbImpl->get_meshset_options(*start_rit, options);
+      RR("Failed to get meshset options.");
       optionsVec.push_back(options);
       count += sizeof(unsigned int);
     
       if (options & MESHSET_SET) {
           // range-based set; count the subranges
         setRanges.push_back(MBRange());
-        result = mbImpl->get_entities_by_handle(*start_rit, *setRanges.rbegin()); RR;
+        result = mbImpl->get_entities_by_handle(*start_rit, *setRanges.rbegin());
+        RR("Failed to get set entities.");
         count += 2 * sizeof(MBEntityHandle) * num_subranges(*setRanges.rbegin()) + sizeof(int);
       }
       else if (options & MESHSET_ORDERED) {
           // just get the number of entities in the set
         int num_ents;
-        result = mbImpl->get_number_entities_by_handle(*start_rit, num_ents); RR;
+        result = mbImpl->get_number_entities_by_handle(*start_rit, num_ents);
+        RR("Failed to get number entities in ordered set.");
+        
         count += sizeof(int);
         
         setSizes.push_back(num_ents);
@@ -633,8 +683,11 @@
 
         // get numbers of parents/children
       int num_par, num_ch;
-      result = mbImpl->num_child_meshsets(*start_rit, &num_ch); RR;
-      result = mbImpl->num_parent_meshsets(*start_rit, &num_par); RR;
+      result = mbImpl->num_child_meshsets(*start_rit, &num_ch);
+      RR("Failed to get num children.");
+
+      result = mbImpl->num_parent_meshsets(*start_rit, &num_par);
+      RR("Failed to get num parents.");
       count += 2*sizeof(int) + (num_par + num_ch) * sizeof(MBEntityHandle);
     
     }
@@ -663,14 +716,16 @@
           // pack entities as vector, with length
         PACK_INT(buff_ptr, *mem_it);
         members.clear();
-        result = mbImpl->get_entities_by_handle(*set_it, members); RR;
+        result = mbImpl->get_entities_by_handle(*set_it, members);
+        RR("Failed to get set entities.");
         PACK_EH(buff_ptr, &members[0], *mem_it);
         mem_it++;
       }
       
         // pack parents
       members.clear();
-      result = mbImpl->get_parent_meshsets(*set_it, members); RR;
+      result = mbImpl->get_parent_meshsets(*set_it, members);
+      RR("Failed to pack parents.");
       PACK_INT(buff_ptr, members.size());
       if (!members.empty()) {
         PACK_EH(buff_ptr, &members[0], members.size());
@@ -678,7 +733,8 @@
       
         // pack children
       members.clear();
-      result = mbImpl->get_child_meshsets(*set_it, members); RR;
+      result = mbImpl->get_child_meshsets(*set_it, members);
+      RR("Failed to pack children.");
       PACK_INT(buff_ptr, members.size());
       if (!members.empty()) {
         PACK_EH(buff_ptr, &members[0], members.size());
@@ -703,37 +759,38 @@
   std::vector<MBRange>::const_iterator rit = setRanges.begin();
   std::vector<int>::const_iterator mem_it = setSizes.begin();
 
-  MBRange set_handles;
+  MBRange set_handles, new_sets;
   UNPACK_RANGE(buff_ptr, set_handles);
   std::vector<MBEntityHandle> members;
   
-  for (MBRange::const_iterator rit = set_handles.begin(); rit != set_handles.end(); rit++) {
+  for (MBRange::const_iterator rit = set_handles.begin(); 
+       rit != set_handles.end(); rit++) {
     
       // option value
     unsigned int opt;
     UNPACK_VOID(buff_ptr, &opt, sizeof(unsigned int));
       
       // create the set
-    MBEntityHandle set_handle;
-    result = mbImpl->create_meshset(opt, set_handle, 
-                                    mbImpl->handle_utils().id_from_handle(*rit), 
-                                    mbImpl->handle_utils().rank_from_handle(*rit)); RR;
-    if (set_handle != *rit)
-      return MB_FAILURE;
+    MBEntityHandle set_handle = *rit;
+    result = sequenceManager->allocate_mesh_set(set_handle, opt);
+    RR("Failed to create set in unpack.");
+    new_sets.insert(set_handle);
 
     int num_ents;
     if (opt & MESHSET_SET) {
         // unpack entities as a range
       MBRange set_range;
       UNPACK_RANGE(buff_ptr, set_range);
-      result = mbImpl->add_entities(*rit, set_range); RR;
+      result = mbImpl->add_entities(*rit, set_range);
+      RR("Failed to add ents to set in unpack.");
     }
     else if (opt & MESHSET_ORDERED) {
         // unpack entities as vector, with length
       UNPACK_INT(buff_ptr, num_ents);
       members.reserve(num_ents);
       UNPACK_EH(buff_ptr, &members[0], num_ents);
-      result = mbImpl->add_entities(*rit, &members[0], num_ents); RR;
+      result = mbImpl->add_entities(*rit, &members[0], num_ents);
+      RR("Failed to add ents to ordered set in unpack.");
     }
       
       // unpack parents/children
@@ -741,15 +798,19 @@
     members.reserve(num_ents);
     UNPACK_EH(buff_ptr, &members[0], num_ents);
     for (int i = 0; i < num_ents; i++) {
-      result = mbImpl->add_parent_meshset(*rit, members[i]); RR;
+      result = mbImpl->add_parent_meshset(*rit, members[i]);
+      RR("Failed to add parent to set in unpack.");
     }
     UNPACK_INT(buff_ptr, num_ents);
     members.reserve(num_ents);
     UNPACK_EH(buff_ptr, &members[0], num_ents);
     for (int i = 0; i < num_ents; i++) {
-      result = mbImpl->add_child_meshset(*rit, members[i]); RR;
+      result = mbImpl->add_child_meshset(*rit, members[i]);
+      RR("Failed to add child to set in unpack.");
     }
   }
+
+  entities.merge(new_sets);
   
   return MB_SUCCESS;
 }
@@ -785,25 +846,22 @@
   count = 0;
   unsigned char *orig_buff_ptr = buff_ptr;
   MBErrorCode result;
-  int whole_size = whole_range.size();
 
   if (just_count) {
 
     std::vector<MBTag> all_tags;
-    result = tagServer->get_tags(all_tags); RR;
+    result = tagServer->get_tags(all_tags);
+    RR("Failed to get tags in pack_tags.");
 
+
     for (std::vector<MBTag>::iterator tag_it = all_tags.begin(); tag_it != all_tags.end(); tag_it++) {
       const TagInfo *tinfo = tagServer->get_tag_info(*tag_it);
       int this_count = 0;
       MBRange tmp_range;
-      if (PROP_FROM_TAG_HANDLE(*tag_it) == MB_TAG_DENSE) {
-        this_count += whole_size * tinfo->get_size();
-      }
-      else {
-        result = tagServer->get_entities(*tag_it, MBMAXTYPE, tmp_range); RR;
-        tmp_range = tmp_range.intersect(whole_range);
-        if (!tmp_range.empty()) this_count = tmp_range.size() * tinfo->get_size();
-      }
+      result = tagServer->get_entities(*tag_it, tmp_range);
+      RR("Failed to get entities for tag in pack_tags.");
+      tmp_range = tmp_range.intersect(whole_range);
+      if (!tmp_range.empty()) this_count = tmp_range.size() * tinfo->get_size();
 
       if (0 == this_count) continue;
 
@@ -869,19 +927,12 @@
         // name
       PACK_CHAR_64(buff_ptr, tinfo->get_name().c_str());
       
-      if (PROP_FROM_TAG_HANDLE(*tag_it) == MB_TAG_DENSE) {
-        tag_data.reserve((whole_size+1) * tinfo->get_size() / sizeof(int));
-        result = mbImpl->tag_get_data(*tag_it, whole_range, &tag_data[0]);
-        PACK_VOID(buff_ptr, &tag_data[0], whole_size*tinfo->get_size());
-      }
-      else {
-        tag_data.reserve((tr_it->size()+1) * tinfo->get_size() / sizeof(int));
-        result = mbImpl->tag_get_data(*tag_it, *tr_it, &tag_data[0]); RR;
-        PACK_RANGE(buff_ptr, (*tr_it));
-        PACK_VOID(buff_ptr, &tag_data[0], tr_it->size()*tinfo->get_size());
-        tr_it++;
-      }
-      
+      tag_data.reserve((tr_it->size()+1) * tinfo->get_size() / sizeof(int));
+      result = mbImpl->tag_get_data(*tag_it, *tr_it, &tag_data[0]);
+      RR("Failed to get tag data in pack_tags.");
+      PACK_RANGE(buff_ptr, (*tr_it));
+      PACK_VOID(buff_ptr, &tag_data[0], tr_it->size()*tinfo->get_size());
+      tr_it++;
     }
 
     count = buff_ptr - orig_buff_ptr;
@@ -930,9 +981,10 @@
 
       // create the tag
     MBTagType tag_type;
-    result = mbImpl->tag_get_type(tag_handle, tag_type); RR;
-
-    result = mbImpl->tag_create(tag_name, tag_size, tag_type, (MBDataType) tag_data_type, tag_handle,
+    result = mbImpl->tag_get_type(tag_handle, tag_type);
+    RR("Failed to get tag type in unpack_tags.");
+    result = mbImpl->tag_create(tag_name, tag_size, tag_type, 
+                                (MBDataType) tag_data_type, tag_handle,
                                 def_val_ptr);
     if (MB_ALREADY_ALLOCATED == result) {
         // already allocated tag, check to make sure it's the same size, type, etc.
@@ -940,61 +992,24 @@
       if (tag_size != tag_info->get_size() ||
           tag_data_type != tag_info->get_data_type() ||
           (def_val_ptr && !tag_info->default_value() ||
-           !def_val_ptr && tag_info->default_value()))
-        return MB_FAILURE;
+           !def_val_ptr && tag_info->default_value())) {
+        RR("Didn't get correct tag info when unpacking tag.");
+      }
       MBTagType this_type;
       result = mbImpl->tag_get_type(tag_handle, this_type);
-      if (MB_SUCCESS != result || this_type != tag_type) return MB_FAILURE;
+      if (MB_SUCCESS != result || this_type != tag_type) {
+        RR("Didn't get correct tag type when unpacking tag.");
+      }
     }
     else if (MB_SUCCESS != result) return result;
     
-      // set the tag data
-    if (PROP_FROM_TAG_HANDLE(tag_handle) == MB_TAG_DENSE) {
-      if (NULL != def_val_ptr && tag_data_type != MB_TYPE_OPAQUE) {
-          // only set the tags whose values aren't the default value; only works
-          // if it's a known type
-        MBRange::iterator start_rit = entities.begin(), end_rit = start_rit;
-        MBRange set_ents;
-        while (end_rit != entities.end()) {
-          while (start_rit != entities.end() &&
-                 ((tag_data_type == MB_TYPE_INTEGER && *((int*)def_val_ptr) == *((int*)buff_ptr)) ||
-                  (tag_data_type == MB_TYPE_DOUBLE && *((double*)def_val_ptr) == *((double*)buff_ptr)) ||
-                  (tag_data_type == MB_TYPE_HANDLE && *((MBEntityHandle*)def_val_ptr) == *((MBEntityHandle*)buff_ptr)))) {
-            start_rit++;
-            buff_ptr += tag_size;
-          }
-          end_rit = start_rit;
-          void *end_ptr = buff_ptr;
-          while (start_rit != entities.end() && end_rit != entities.end() &&
-                 ((tag_data_type == MB_TYPE_INTEGER && *((int*)def_val_ptr) == *((int*)end_ptr)) ||
-                  (tag_data_type == MB_TYPE_DOUBLE && *((double*)def_val_ptr) == *((double*)end_ptr)) ||
-                  (tag_data_type == MB_TYPE_HANDLE && *((MBEntityHandle*)def_val_ptr) == *((MBEntityHandle*)end_ptr)))) {
-            set_ents.insert(*end_rit);
-            end_rit++;
-            buff_ptr += tag_size;
-          }
-          
-          if (!set_ents.empty()) {
-            result = mbImpl->tag_set_data(tag_handle, set_ents, buff_ptr); RR;
-          }
-          if (start_rit != entities.end()) {
-            end_rit++;
-            start_rit = end_rit;
-            buff_ptr += tag_size;
-          }
-        }
-      }
-      else {
-        result = mbImpl->tag_set_data(tag_handle, entities, buff_ptr); RR;
-        buff_ptr += entities.size() * tag_size;
-      }
-    }
-    else {
-      MBRange tag_range;
-      UNPACK_RANGE(buff_ptr, tag_range);
-      result = mbImpl->tag_set_data(tag_handle, tag_range, buff_ptr); RR;
-      buff_ptr += tag_range.size() * tag_size;
-    }
+      // set the tag data; don't have to worry about dense tags with default
+      // values, as those aren't sent
+    MBRange tag_range;
+    UNPACK_RANGE(buff_ptr, tag_range);
+    result = mbImpl->tag_set_data(tag_handle, tag_range, buff_ptr);
+    RR("Trouble setting range-based tag data when unpacking tag.");
+    buff_ptr += tag_range.size() * tag_size;
   }
   
   return MB_SUCCESS;
@@ -1012,30 +1027,66 @@
   new_buffer.swap(myBuffer);
 }
 
+MBErrorCode MBParallelComm::resolve_shared_ents(int dim,
+                                                int shared_dim) 
+{
+  MBErrorCode result;
+  MBRange proc_ents;
+  if (-1 == dim) {
+    int this_dim = 3;
+    while (proc_ents.empty() && this_dim >= 0) {
+      result = mbImpl->get_entities_by_dimension(0, this_dim, proc_ents);
+      if (MB_SUCCESS != result) return result;
+      this_dim--;
+    }
+  }
+  else {
+    result = mbImpl->get_entities_by_dimension(0, dim, proc_ents);
+    if (MB_SUCCESS != result) return result;
+  }
+
+  if (proc_ents.empty()) return MB_SUCCESS;
+  
+  return resolve_shared_ents(proc_ents, shared_dim);
+}
+  
 MBErrorCode MBParallelComm::resolve_shared_ents(MBRange &proc_ents,
-                                                const int dim) 
+                                                int shared_dim) 
 {
   MBRange::iterator rit;
   MBSkinner skinner(mbImpl);
   
     // get the skin entities by dimension
-  MBRange skin_ents[3];
+  MBRange skin_ents[4];
   MBErrorCode result;
   int upper_dim = MBCN::Dimension(TYPE_FROM_HANDLE(*proc_ents.begin()));
 
-  if (upper_dim > 0) {
-      // first get d-1 skin ents
-    result = skinner.find_skin(proc_ents, skin_ents[upper_dim-1],
-                               skin_ents[upper_dim-1]);
+  int skin_dim;
+  if (shared_dim < upper_dim) {
+      // if shared entity dimension is less than maximal dimension,
+      // start with skin entities
+    skin_dim = upper_dim-1;
+    result = skinner.find_skin(proc_ents, skin_ents[skin_dim],
+                               skin_ents[skin_dim]);
     if (MB_SUCCESS != result) return result;
-      // then get d-2, d-3, etc. entities adjacent to skin ents 
-    for (int this_dim = upper_dim-1; this_dim >= 0; this_dim--) {
-      result = mbImpl->get_adjacencies(skin_ents[upper_dim-1], this_dim,
-                                       true, skin_ents[this_dim]);
-      if (MB_SUCCESS != result) return result;
-    }
   }
-  else skin_ents[0] = proc_ents;
+  else {
+      // otherwise start with original entities
+    skin_ents[upper_dim] = proc_ents;
+    skin_dim = upper_dim;
+  }
+
+    // get entities adjacent to skin ents from shared_dim down to
+    // zero; don't create them if they don't exist already
+  for (int this_dim = shared_dim; this_dim >= 0; this_dim--) {
+
+    if (this_dim == skin_dim) continue;
+      
+    result = mbImpl->get_adjacencies(skin_ents[skin_dim], this_dim,
+                                     false, skin_ents[this_dim],
+                                     MBInterface::UNION);
+    if (MB_SUCCESS != result) return result;
+  }
   
     // global id tag
   MBTag gid_tag; int def_val = -1;
@@ -1043,9 +1094,10 @@
                               MB_TAG_DENSE, MB_TYPE_INTEGER, gid_tag,
                               &def_val, true);
   if (MB_FAILURE == result) return result;
+
   else if (MB_ALREADY_ALLOCATED != result) {
       // just created it, so we need global ids
-    result = assign_global_ids(dim);
+    result = assign_global_ids(upper_dim);
     if (MB_SUCCESS != result) return result;
   }
 
@@ -1083,7 +1135,6 @@
   if (NULL == gsd) return MB_FAILURE;
   
     // get shared proc tags
-#define MAX_SHARING_PROCS 10  
   int def_vals[2] = {-10*procConfig.proc_size(), -10*procConfig.proc_size()};
   MBTag sharedp_tag, sharedps_tag;
   result = mbImpl->tag_create(PARALLEL_SHARED_PROC_TAG_NAME, 2*sizeof(int), 
@@ -1135,7 +1186,7 @@
   
     // set sharing procs tags on other skin ents
   const MBEntityHandle *connect; int num_connect;
-  for (int d = dim-1; d > 0; d--) {
+  for (int d = shared_dim; d > 0; d--) {
     for (MBRange::iterator rit = skin_ents[d].begin();
          rit != skin_ents[d].end(); rit++) {
         // get connectivity
@@ -1185,9 +1236,58 @@
   return result;
 }
 
+MBErrorCode MBParallelComm::get_shared_entities(int dim,
+                                                MBRange &shared_ents) 
+{
+    // check shared entities
+  MBTag sharedproc_tag = 0, sharedprocs_tag = 0;
+  MBErrorCode result = mbImpl->tag_get_handle(PARALLEL_SHARED_PROC_TAG_NAME, 
+                                              sharedproc_tag);
 
+  result = mbImpl->tag_get_handle(PARALLEL_SHARED_PROCS_TAG_NAME, 
+                                  sharedprocs_tag);
+
+  if (0 == sharedproc_tag && 0 == sharedprocs_tag) 
+    return MB_SUCCESS;
+
+    // get the tag values
+  MBEntityType start_type = MBCN::TypeDimensionMap[dim].first,
+    end_type = MBCN::TypeDimensionMap[dim].second;
+  std::vector<int> proc_tags;
+  for (MBEntityType this_type = start_type; this_type <= end_type;
+       this_type++) {
+    MBRange tmp_ents;
+
+      // PARALLEL_SHARED_PROC is a dense tag, so all ents will have a
+      // value (the default value)
+    if (0 != sharedproc_tag) {
+      result = mbImpl->get_entities_by_type(0, this_type, tmp_ents);
+      RR("Trouble getting entities for shared entities.");
+      proc_tags.resize(2*tmp_ents.size());
+      if (!tmp_ents.empty()) {
+        result = mbImpl->tag_get_data(sharedproc_tag, 
+                                      tmp_ents, &proc_tags[0]);
+        RR("Trouble getting tag data for shared entities.");
+      }
+      int i;
+      MBRange::iterator rit;
+      for (i = 0, rit = tmp_ents.begin(); rit != tmp_ents.end(); i+=2, rit++) 
+        if (proc_tags[i] > -1) shared_ents.insert(*rit);
+    }
+    if (0 != sharedprocs_tag) {
+      // PARALLEL_SHARED_PROCS is a sparse tag, so only entities with this
+      // tag set will have one
+      result = mbImpl->get_entities_by_type_and_tag(0, this_type, 
+                                                    &sharedprocs_tag,
+                                                    NULL, 1, tmp_ents,
+                                                    MBInterface::UNION);
+      RR("Trouble getting sharedprocs_tag for shared entities.");
+    }
+  }
+
+  return MB_SUCCESS;
+}
   
-  
 #ifdef TEST_PARALLELCOMM
 
 #include <iostream>
@@ -1196,6 +1296,12 @@
 #include "MBParallelComm.hpp"
 #include "MBRange.hpp"
 
+#define PM {std::cerr << "Test failed; error message:" << std::endl;\
+          std::string errmsg; \
+          dynamic_cast<MBCore*>(my_impl)->get_last_error(errmsg); \
+          std::cerr << errmsg << std::endl;\
+          return 1;}
+
 int main(int argc, char* argv[])
 {
 
@@ -1224,14 +1330,18 @@
   if (MB_SUCCESS != result) return result;
     
   int buff_size;
-  result = pcomm.pack_buffer(all_mesh, false, true, true, whole_range, buff_size); RR;
+  result = pcomm.pack_buffer(all_mesh, false, true, true, whole_range, buff_size);
+  PM;
 
+
     // allocate space in the buffer
   pcomm.buffer_size(buff_size);
 
     // pack the actual buffer
   int actual_buff_size;
-  result = pcomm.pack_buffer(whole_range, false, true, false, all_mesh, actual_buff_size); RR;
+  result = pcomm.pack_buffer(whole_range, false, true, false, all_mesh, 
+                             actual_buff_size);
+  PM;
 
     // list the entities that got packed
   std::cout << "ENTITIES PACKED:" << std::endl;
@@ -1252,7 +1362,9 @@
 
     // unpack the results
   all_mesh.clear();
-  result = pcomm2.unpack_buffer(all_mesh); RR;
+  result = pcomm2.unpack_buffer(all_mesh);
+  PM;
+  
   std::cout << "ENTITIES UNPACKED:" << std::endl;
   mbImpl->list_entities(all_mesh);
   

Modified: MOAB/trunk/parallel/MBParallelComm.hpp
===================================================================
--- MOAB/trunk/parallel/MBParallelComm.hpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/parallel/MBParallelComm.hpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -71,11 +71,40 @@
      * number of sharing processors.  Values in these tags denote the ranks
      * of sharing processors, and the list ends with the value -10*#procs.
      *
+     * If shared_dim is input as -1 or not input, a value one less than the
+     * maximum dimension of entities in proc_ents is used.
+     *
      * \param proc_ents Entities for which to resolve shared entities
-     * \param dim Dimension of entities in proc_ents
+     * \param shared_dim Maximum dimension of shared entities to look for
      */
-  MBErrorCode resolve_shared_ents(MBRange &proc_ents, const int dim);
+  MBErrorCode resolve_shared_ents(MBRange &proc_ents, 
+                                  int shared_dim = 0);
   
+    /** Resolve shared entities between processors
+     * Same as resolve_shared_ents(MBRange&), except works for
+     * all entities in instance of dimension dim.  
+     *
+     * If dim = -1 or no dim input, uses entities of maximal
+     * dimension (3) in the instance.
+     *
+     * If shared_dim is input as -1 or not input, a value one less than the
+     * maximum dimension of entities is used.
+
+     * \param dim Maximal dimension of entities to be resolved
+     * \param shared_dim Maximum dimension of shared entities to look for
+     */
+  MBErrorCode resolve_shared_ents(int dim = -1, 
+                                  int shared_dim = 0);
+  
+    /** Get entities shared with other processors, based on 
+     * PARALLEL_SHARED_PROC_TAG_NAME and PARALLEL_SHARED_PROCS_TAG_NAME.
+     *
+     * \param dim Dimension of entities shared with other processors
+     * \param shared_ents Shared entities returned in this range
+     */
+  MBErrorCode get_shared_entities(int dim,
+                                  MBRange &shared_ents);
+
     //! pack a buffer (stored in this class instance) with ALL data for these entities
   MBErrorCode pack_buffer(MBRange &entities, 
                           const bool adjacencies,

Modified: MOAB/trunk/parallel/ReadParallel.cpp
===================================================================
--- MOAB/trunk/parallel/ReadParallel.cpp	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/parallel/ReadParallel.cpp	2007-10-08 19:15:00 UTC (rev 1301)
@@ -6,9 +6,12 @@
 #include "MBReaderWriterSet.hpp"
 #include "MBReadUtilIface.hpp"
 #include "MBParallelComm.hpp"
+#include "MBParallelConventions.h"
 #include "MBCN.hpp"
 
-#define RR if (MB_SUCCESS != result) return result
+#define RR(a) if (MB_SUCCESS != result) {\
+          dynamic_cast<MBCore*>(mbImpl)->get_error_handler()->set_last_error(a);\
+          return result;}
 
 MBErrorCode ReadParallel::load_file(const char *file_name,
                                     MBEntityHandle& file_set,
@@ -22,32 +25,34 @@
   
     // Get parallel settings
   int parallel_mode;
-  const char* parallel_opts[] = { "NONE", "BCAST", "BCAST_DELETE", "SCATTER", 
+  const char* parallel_opts[] = { "NONE", "BCAST", "BCAST_DELETE", 
+                                  "READ_DELETE", "SCATTER", 
                                   "FORMAT", 0 };
-  enum ParallelOpts {POPT_NONE=0, POPT_BCAST, POPT_BCAST_DELETE, POPT_SCATTER,
+  enum ParallelOpts {POPT_NONE=0, POPT_BCAST, POPT_BCAST_DELETE, 
+                     POPT_READ_DELETE, POPT_SCATTER,
                      POPT_FORMAT, POPT_LAST};
       
-  MBErrorCode rval = opts.match_option( "PARALLEL", parallel_opts, 
+  MBErrorCode result = opts.match_option( "PARALLEL", parallel_opts, 
                                         parallel_mode );
-  if (MB_FAILURE == rval) {
+  if (MB_FAILURE == result) {
     merror->set_last_error( "Unexpected value for 'PARALLEL' option\n" );
     return MB_FAILURE;
   }
-  else if (MB_ENTITY_NOT_FOUND == rval) {
+  else if (MB_ENTITY_NOT_FOUND == result) {
     parallel_mode = 0;
   }
     // Get partition setting
   std::string partition_tag_name;
-  rval = opts.get_option("PARTITION", partition_tag_name);
-  if (MB_ENTITY_NOT_FOUND == rval || partition_tag_name.empty())
+  result = opts.get_option("PARTITION", partition_tag_name);
+  if (MB_ENTITY_NOT_FOUND == result || partition_tag_name.empty())
     partition_tag_name += "PARTITION";
 
     // get MPI IO processor rank
   int reader_rank;
-  rval = opts.get_int_option( "MPI_IO_RANK", reader_rank );
-  if (MB_ENTITY_NOT_FOUND == rval)
+  result = opts.get_int_option( "MPI_IO_RANK", reader_rank );
+  if (MB_ENTITY_NOT_FOUND == result)
     reader_rank = 0;
-  else if (MB_SUCCESS != rval) {
+  else if (MB_SUCCESS != result) {
     merror->set_last_error( "Unexpected value for 'MPI_IO_RANK' option\n" );
     return MB_FAILURE;
   }
@@ -65,14 +70,15 @@
     return MB_NOT_IMPLEMENTED;
   }
 
-  if (parallel_mode != POPT_SCATTER || 
+  if ((parallel_mode != POPT_SCATTER && 
+       parallel_mode != POPT_BCAST_DELETE) || 
       reader_rank == (int)(mbImpl->proc_rank())) {
       // Try using the file extension to select a reader
     const MBReaderWriterSet* set = impl->reader_writer_set();
     MBReaderIface* reader = set->get_file_extension_reader( file_name );
     if (reader)
     { 
-      rval = reader->load_file( file_name, file_set, opts, 
+      result = reader->load_file( file_name, file_set, opts, 
                                 material_set_list, num_material_sets );
       delete reader;
     }
@@ -85,69 +91,134 @@
         MBReaderIface* reader = iter->make_reader( mbImpl );
         if (NULL != reader)
         {
-          rval = reader->load_file( file_name, file_set, opts, 
+          result = reader->load_file( file_name, file_set, opts, 
                                     material_set_list, num_material_sets );
           delete reader;
-          if (MB_SUCCESS == rval)
+          if (MB_SUCCESS == result)
             break;
         }
       }
     }
   }
   else {
-    rval = MB_SUCCESS;
+    result = MB_SUCCESS;
   }
   
+  if (MB_SUCCESS != result)
+    RR("Failed initial file load.");
+    
   if (parallel_mode == POPT_BCAST ||
       parallel_mode == POPT_BCAST_DELETE) {
     MBRange entities; 
-    if (MB_SUCCESS == rval && 
-        reader_rank == (int)(mbImpl->proc_rank())) {
-      rval = mbImpl->get_entities_by_handle( file_set, entities );
-      if (MB_SUCCESS != rval)
+
+      // get which entities need to be broadcast, only if I'm the reader
+    if (reader_rank == (int)(mbImpl->proc_rank())) {
+      result = mbImpl->get_entities_by_handle( file_set, entities );
+      if (MB_SUCCESS != result)
         entities.clear();
+
+        // add actual file set to entities too
+      entities.insert(file_set);
+
+        // mark the file set so the receivers know which one it is
+      MBTag file_set_tag;
+      int other_sets = 0;
+      result = mbImpl->tag_create("FILE_SET", sizeof(int), MB_TAG_SPARSE, 
+                                  MB_TYPE_INTEGER, file_set_tag, 0, true);
+      if (MB_ALREADY_ALLOCATED == result) {
+        MBRange other_file_sets;
+        result = mbImpl->get_entities_by_type_and_tag(0, MBENTITYSET,
+                                                      &file_set_tag, NULL, 1,
+                                                      other_file_sets);
+        if (MB_SUCCESS == result) other_sets = other_file_sets.size();
+      }
+      if (MB_SUCCESS == result)
+        result = mbImpl->tag_set_data(file_set_tag, &file_set, 1, &other_sets);
     }
-    
+
+      // do the actual broadcast; if single-processor, ignore error
     MBParallelComm tool( mbImpl);
-    MBErrorCode tmp_rval = tool.broadcast_entities( reader_rank, entities );
-    if (MB_SUCCESS != rval && mbImpl->proc_size() != 1)
-      tmp_rval = rval;
-    else if (MB_SUCCESS != rval) rval = MB_SUCCESS;
+    result = tool.broadcast_entities( reader_rank, entities );
+    if (mbImpl->proc_size() == 1 && MB_SUCCESS != result) 
+      result = MB_SUCCESS;
       
-    if (MB_SUCCESS == rval && 
+      // go get the file set if I'm not the reader
+    if (MB_SUCCESS == result && 
         reader_rank != (int)(mbImpl->proc_rank())) {
-      rval = mbImpl->create_meshset( MESHSET_SET, file_set );
-      if (MB_SUCCESS == rval) {
-        rval = mbImpl->add_entities( file_set, entities );
-        if (MB_SUCCESS != rval) {
-          mbImpl->delete_entities( &file_set, 1 );
-          file_set = 0;
+      MBTag file_set_tag;
+      result = mbImpl->tag_get_handle("FILE_SET", file_set_tag);
+      if (MB_SUCCESS == result) {
+        MBRange other_file_sets;
+        result = mbImpl->get_entities_by_type_and_tag(0, MBENTITYSET,
+                                                      &file_set_tag, NULL, 1,
+                                                      other_file_sets);
+        if (MB_SUCCESS == result && other_file_sets.size() > 1) {
+          int tag_val = other_file_sets.size(), *tag_val_ptr = &tag_val;
+          MBRange file_sets;
+          result = mbImpl->get_entities_by_type_and_tag(0, MBENTITYSET,
+                                                        &file_set_tag, 
+                                                        (void*const*) &tag_val_ptr, 1,
+                                                        file_sets);
+          if (!file_sets.empty()) other_file_sets = file_sets;
         }
+        if (!other_file_sets.empty()) file_set = *other_file_sets.rbegin();
       }
     }
-
-    if (parallel_mode == POPT_BCAST_DELETE)
-      rval = delete_nonlocal_entities(partition_tag_name, file_set);
     
+    RR("Failed to broadcast mesh.");
   }
+
+  if (parallel_mode == POPT_BCAST_DELETE ||
+      parallel_mode == POPT_READ_DELETE) {
+    result = delete_nonlocal_entities(partition_tag_name, file_set);
+    if (MB_SUCCESS != result) return result;
+  }
   
-  return rval;
+  return result;
 }
 
 MBErrorCode ReadParallel::delete_nonlocal_entities(std::string &ptag_name,
                                                    MBEntityHandle file_set) 
 {
-  MBRange partition_sets;
+  MBRange partition_sets, my_sets;
   MBErrorCode result;
 
   MBTag ptag;
-  result = mbImpl->tag_get_handle(ptag_name.c_str(), ptag); RR;
+  result = mbImpl->tag_get_handle(ptag_name.c_str(), ptag); 
+  RR("Failed getting tag handle in delete_nonlocal_entities.");
+
+  if (ptag_name == PARALLEL_PARTITION_TAG_NAME) {
+    int my_rank = mbImpl->proc_rank(), *my_rank_ptr = &my_rank;
+    
+    result = mbImpl->get_entities_by_type_and_tag(file_set, MBENTITYSET,
+                                                  &ptag, 
+                                                  (void* const*) &my_rank_ptr, 1,
+                                                  my_sets); 
+    RR("Couldn't get entities with partition tag.");
+  }
+  else {
+      // no explicit stategy to assign sets, so just balance # sets
+    result = mbImpl->get_entities_by_type_and_tag(file_set, MBENTITYSET,
+                                                  &ptag, NULL, 1,
+                                                  partition_sets);
+    RR("Failed to get sets with partition-type tag.");
+
+    int proc_sz = mbImpl->proc_size();
+    int proc_rk = mbImpl->proc_rank();
+    
+    int num_sets = partition_sets.size() / proc_sz, orig_numsets = num_sets;
+    if (partition_sets.size() % proc_sz != 0) {
+      num_sets++;
+      if (proc_rk == proc_sz-1) 
+        num_sets = partition_sets.size() % num_sets;
+    }
+
+    int istart = orig_numsets * proc_rk;
+    for (int i = 0; i < num_sets; i++) 
+      my_sets.insert(partition_sets[istart+i]);
+  }
   
-  result = mbImpl->get_entities_by_type_and_tag(file_set, MBENTITYSET,
-                                                &ptag, NULL, 1,
-                                                partition_sets); RR;
-
-  return delete_nonlocal_entities(partition_sets, file_set);
+  return delete_nonlocal_entities(my_sets, file_set);
 }
 
 MBErrorCode ReadParallel::delete_nonlocal_entities(MBRange &partition_sets,
@@ -164,11 +235,12 @@
   MBRange partition_ents, all_sets;
   result = read_iface->gather_related_ents(partition_sets, partition_ents,
                                            &all_sets);
-  RR;
+  RR("Failure gathering related entities.");
 
     // get pre-existing entities
   MBRange file_ents;
-  result = mbImpl->get_entities_by_handle(file_set, file_ents); RR;
+  result = mbImpl->get_entities_by_handle(file_set, file_ents); 
+  RR("Couldn't get pre-existing entities.");
 
     // get deletable entities by subtracting partition ents from file ents
   MBRange deletable_ents = file_ents.subtract(partition_ents);
@@ -180,14 +252,17 @@
     // remove deletable ents from all keepable sets
   for (MBRange::iterator rit = keepable_sets.begin();
        rit != keepable_sets.end(); rit++) {
-    result = mbImpl->remove_entities(*rit, deletable_ents); RR;
+    result = mbImpl->remove_entities(*rit, deletable_ents);
+    RR("Failure removing deletable entities.");
   }
 
     // delete sets, then ents
-  result = mbImpl->delete_entities(deletable_sets); RR;
+  result = mbImpl->delete_entities(deletable_sets);
+  RR("Failure deleting sets in delete_nonlocal_entities.");
 
   deletable_ents = deletable_ents.subtract(deletable_sets);
-  result = mbImpl->delete_entities(deletable_ents); RR;
+  result = mbImpl->delete_entities(deletable_ents);
+  RR("Failure deleting entities in delete_nonlocal_entities.");
   
   result = ((MBCore*)mbImpl)->check_adjacencies();
 

Modified: MOAB/trunk/tools/iMesh/Makefile.am
===================================================================
--- MOAB/trunk/tools/iMesh/Makefile.am	2007-10-05 23:10:57 UTC (rev 1300)
+++ MOAB/trunk/tools/iMesh/Makefile.am	2007-10-08 19:15:00 UTC (rev 1301)
@@ -16,7 +16,7 @@
 check_PROGRAMS = testc_cbind
 testc_cbind_SOURCES = testc_cbind.c
 
-LDADD = libiMesh.la $(top_builddir)/libMOAB.la
+LDADD = libiMesh.la $(top_builddir)/libMOAB.la ${MOAB_CXX_LDFLAGS} ${MOAB_CXX_LIBS}
 testc_cbind_DEPENDENCIES = libiMesh.la $(top_builddir)/libMOAB.la
 
 libiMesh_la_SOURCES = \




More information about the moab-dev mailing list